Merge pull request 'feat: skin sub budget neutrality notebook' (#291) from feat/skin-sub-budget-neutrality into main

notebooks/skin_sub_budget_neutrality.py

@@ -0,0 +1,564 @@
+import marimo
+
+__generated_with = "0.21.1"
+app = marimo.App(width="medium")
+
+
+@app.cell(hide_code=True)
+def _():
+    import marimo as mo
+
+    return (mo,)
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    # PFS Budget Neutrality Impact of Skin Substitute Codes
+
+    The Physician Fee Schedule is **budget neutral** — when CMS raises
+    RVUs for any service, it must lower them elsewhere (or reduce the
+    conversion factor) so aggregate spending stays flat. This notebook
+    quantifies how skin substitute application code revaluations
+    redistribute payment away from other services.
+
+    **Mechanism:** CMS publishes ~11,000 HCPCS codes with Work, PE, and
+    MP RVUs. When the sum of (RVU × frequency) grows, the conversion
+    factor or a budget neutrality adjustor (BNA) scales down to
+    compensate. Every code shares the compression proportionally.
+
+    **What this means:** A PE RVU increase on 15271–15278 is not free
+    money — it is a tax on every other service in the fee schedule.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _():
+    import altair as alt
+    import polars as pl
+    from conf import connect
+    from pfs.rules import RULES
+
+    con = connect.duckdb()
+
+    def q(sql):
+        return con.execute(sql).pl()
+
+    SKIN_CODES = "('15271','15272','15273','15274','15275','15276','15277','15278')"
+
+    return RULES, SKIN_CODES, alt, con, pl, q
+
+
+# ── 1. Skin sub RVU growth vs. total pool ────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 1. Skin Substitute Share of the RVU Pool
+
+    The unweighted RVU pool (sum of all base-mod RVUs across ~11k codes)
+    grows over time as CMS adds codes and revalues services. The skin
+    sub application codes (15271–15278) are a small but growing share.
+
+    > **Note:** Budget neutrality operates on *frequency-weighted* RVUs
+    > (RVU × utilization), not unweighted sums. Without CMS utilization
+    > data, this analysis uses unweighted RVUs as a structural proxy.
+    > The actual impact is amplified by the explosive volume growth in
+    > skin substitute claims documented by OIG.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(SKIN_CODES, alt, mo, q):
+    pool_share = q(f"""
+        WITH pool AS (
+            SELECT year,
+                   sum(work_rvu + non_fac_pe_rvu + mp_rvu) as total_rvu,
+                   sum(CASE WHEN hcpcs IN {SKIN_CODES}
+                       THEN work_rvu + non_fac_pe_rvu + mp_rvu ELSE 0 END) as skin_rvu,
+                   sum(CASE WHEN hcpcs NOT IN {SKIN_CODES}
+                       THEN work_rvu + non_fac_pe_rvu + mp_rvu ELSE 0 END) as other_rvu
+            FROM pfs.rvu
+            WHERE mod IS NULL OR mod = ''
+            GROUP BY year
+        )
+        SELECT year,
+               round(skin_rvu, 2) as skin_rvu,
+               round(total_rvu, 1) as total_rvu,
+               round(100.0 * skin_rvu / total_rvu, 4) as skin_pct,
+               round(skin_rvu - LAG(skin_rvu) OVER (ORDER BY year), 2) as skin_delta,
+               round(total_rvu - LAG(total_rvu) OVER (ORDER BY year), 1) as pool_delta
+        FROM pool ORDER BY year
+    """)
+
+    share_chart = (
+        alt.Chart(pool_share.to_pandas())
+        .mark_bar()
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("skin_rvu:Q", title="Skin Sub Total NF RVUs (8 codes)"),
+            tooltip=["year", "skin_rvu", "total_rvu", "skin_pct", "skin_delta"],
+        )
+        .properties(title="Skin Sub Application Codes — Total NF RVUs by Year", width=700, height=300)
+    )
+
+    mo.vstack([share_chart, pool_share])
+    return (pool_share,)
+
+
+# ── 2. PE RVU revaluation trajectory ─────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 2. Practice Expense RVU Revaluation
+
+    PE is the largest RVU component for skin sub application codes and
+    is where the budget neutrality tax bites hardest. When CMS increases
+    PE RVUs for these codes (e.g., to reflect updated clinical labor
+    rates or supply costs), **all other codes' PE RVUs must absorb a
+    compensating reduction** via the BNA.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(alt, q):
+    pe_trajectory = q("""
+        SELECT year, hcpcs,
+               CASE hcpcs
+                   WHEN '15271' THEN '15271 trunk <100cm²'
+                   WHEN '15272' THEN '15272 trunk add-on'
+                   WHEN '15275' THEN '15275 face <100cm²'
+                   WHEN '15276' THEN '15276 face add-on'
+               END as label,
+               non_fac_pe_rvu,
+               work_rvu,
+               mp_rvu,
+               non_fac_pe_rvu + work_rvu + mp_rvu as total_nf_rvu
+        FROM pfs.rvu
+        WHERE hcpcs IN ('15271','15272','15275','15276')
+          AND (mod IS NULL OR mod = '')
+        ORDER BY year, hcpcs
+    """)
+
+    pe_chart = (
+        alt.Chart(pe_trajectory.to_pandas())
+        .mark_line(point=True)
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("non_fac_pe_rvu:Q", title="Non-Facility PE RVU"),
+            color=alt.Color("label:N", title="Code"),
+            tooltip=["year", "hcpcs", "label", "non_fac_pe_rvu", "work_rvu", "total_nf_rvu"],
+        )
+        .properties(title="Practice Expense RVU Trajectory", width=700, height=350)
+    )
+    pe_chart
+    return
+
+
+# ── 3. Implied budget neutrality tax ─────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 3. Implied Budget Neutrality Tax
+
+    When skin sub PE RVUs increase by Δ, **every other code's effective
+    payment decreases** proportionally. The "tax rate" is:
+
+    ```
+    tax_rate = skin_sub_PE_delta / total_pool_PE
+    ```
+
+    This table shows the year-over-year PE RVU increase for the 8 skin
+    sub codes, the total PE pool, and the implied compression on all
+    other codes.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(SKIN_CODES, q):
+    bn_tax = q(f"""
+        WITH yearly AS (
+            SELECT year,
+                   sum(CASE WHEN hcpcs IN {SKIN_CODES} THEN non_fac_pe_rvu ELSE 0 END) as skin_pe,
+                   sum(non_fac_pe_rvu) as total_pe,
+                   sum(CASE WHEN hcpcs NOT IN {SKIN_CODES} THEN non_fac_pe_rvu ELSE 0 END) as other_pe
+            FROM pfs.rvu
+            WHERE mod IS NULL OR mod = ''
+            GROUP BY year
+        )
+        SELECT year,
+               round(skin_pe, 2) as skin_pe_rvu,
+               round(total_pe, 1) as total_pe_pool,
+               round(skin_pe - LAG(skin_pe) OVER (ORDER BY year), 2) as skin_pe_delta,
+               round(total_pe - LAG(total_pe) OVER (ORDER BY year), 1) as pool_pe_delta,
+               -- If skin PE grew and pool grew less, the difference is absorbed by others
+               round(CASE WHEN LAG(skin_pe) OVER (ORDER BY year) IS NOT NULL
+                     THEN (skin_pe - LAG(skin_pe) OVER (ORDER BY year))
+                          / NULLIF(LAG(total_pe) OVER (ORDER BY year), 0) * 100
+                     END, 4) as implied_tax_pct,
+               -- Dollar impact: tax_pct × average CF
+               round(CASE WHEN LAG(skin_pe) OVER (ORDER BY year) IS NOT NULL
+                     THEN (skin_pe - LAG(skin_pe) OVER (ORDER BY year))
+                          / NULLIF(LAG(total_pe) OVER (ORDER BY year), 0) * 100
+                     END, 4) as pct_compression
+        FROM yearly ORDER BY year
+    """)
+    bn_tax
+    return (bn_tax,)
+
+
+# ── 4. Which services bear the burden? ────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 4. Which Services Bear the Burden?
+
+    Budget neutrality compression is proportional to each service
+    category's share of the PE pool. Categories with large PE shares
+    (surgery, cardiology) absorb more dollars even though the per-code
+    reduction is tiny.
+
+    The table shows: if the skin sub PE increase in the most recent year
+    were fully offset by compressing other categories, how much does
+    each category lose?
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(SKIN_CODES, alt, mo, q):
+    # Get the latest year's skin sub PE delta
+    skin_pe_delta = q(f"""
+        WITH yearly AS (
+            SELECT year, sum(non_fac_pe_rvu) as skin_pe
+            FROM pfs.rvu
+            WHERE hcpcs IN {SKIN_CODES} AND (mod IS NULL OR mod = '')
+            GROUP BY year
+        )
+        SELECT year, skin_pe,
+               skin_pe - LAG(skin_pe) OVER (ORDER BY year) as delta
+        FROM yearly ORDER BY year DESC LIMIT 1
+    """)
+
+    delta_val = skin_pe_delta.select("delta").item()
+    delta_year = skin_pe_delta.select("year").item()
+
+    category_impact = q(f"""
+        SELECT
+            CASE
+                WHEN hcpcs BETWEEN '99201' AND '99499' THEN 'E/M'
+                WHEN hcpcs BETWEEN '10000' AND '19999' THEN 'Integumentary'
+                WHEN hcpcs BETWEEN '20000' AND '29999' THEN 'Musculoskeletal'
+                WHEN hcpcs BETWEEN '30000' AND '39999' THEN 'Resp/Cardiovascular'
+                WHEN hcpcs BETWEEN '40000' AND '49999' THEN 'Digestive'
+                WHEN hcpcs BETWEEN '50000' AND '59999' THEN 'Urinary/Genital'
+                WHEN hcpcs BETWEEN '60000' AND '69999' THEN 'Endocrine/Nervous'
+                WHEN hcpcs BETWEEN '70000' AND '79999' THEN 'Radiology'
+                WHEN hcpcs BETWEEN '80000' AND '89999' THEN 'Path/Lab'
+                WHEN hcpcs BETWEEN '90000' AND '99199' THEN 'Medicine'
+                WHEN hcpcs LIKE 'G%' THEN 'G-codes'
+                ELSE 'Other'
+            END as category,
+            count(*) as codes,
+            round(sum(non_fac_pe_rvu), 1) as category_pe,
+            round(100.0 * sum(non_fac_pe_rvu) /
+                  NULLIF((SELECT sum(non_fac_pe_rvu) FROM pfs.rvu
+                   WHERE year={delta_year} AND (mod IS NULL OR mod = '')
+                     AND hcpcs NOT IN {SKIN_CODES}), 0), 2) as pe_share_pct,
+            -- Implied PE reduction absorbed by this category
+            round({delta_val} * sum(non_fac_pe_rvu) /
+                  NULLIF((SELECT sum(non_fac_pe_rvu) FROM pfs.rvu
+                   WHERE year={delta_year} AND (mod IS NULL OR mod = '')
+                     AND hcpcs NOT IN {SKIN_CODES}), 0), 4) as implied_pe_loss
+        FROM pfs.rvu
+        WHERE year = {delta_year}
+          AND (mod IS NULL OR mod = '')
+          AND hcpcs NOT IN {SKIN_CODES}
+        GROUP BY category
+        ORDER BY category_pe DESC
+    """)
+
+    impact_chart = (
+        alt.Chart(category_impact.to_pandas())
+        .mark_bar()
+        .encode(
+            x=alt.X("implied_pe_loss:Q", title=f"Implied PE RVU Loss (from {delta_val:+.2f} skin sub PE delta)"),
+            y=alt.Y("category:N", title="", sort="-x"),
+            color=alt.Color("pe_share_pct:Q", title="PE Pool Share %",
+                            scale=alt.Scale(scheme="reds")),
+            tooltip=["category", "codes", "category_pe", "pe_share_pct", "implied_pe_loss"],
+        )
+        .properties(
+            title=f"Budget Neutrality Burden by Service Category (CY{delta_year})",
+            width=700, height=350,
+        )
+    )
+
+    mo.vstack([
+        mo.md(f"""
+        **CY{delta_year}:** Skin sub application codes gained **{delta_val:+.2f} PE RVUs**.
+        Under budget neutrality, this is redistributed across ~{category_impact.select('codes').sum().item():,} other codes
+        proportional to their PE share.
+        """),
+        impact_chart,
+        category_impact,
+    ])
+    return
+
+
+# ── 5. Conversion factor erosion ─────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 5. Conversion Factor Erosion
+
+    The conversion factor has declined from $36.09 (2020) to $32.35
+    (2025–2026). While this is driven primarily by MACRA spending
+    targets, **RVU pool growth contributes to the pressure.** When total
+    unweighted RVUs grow faster than allowed spending, the CF must
+    decline to maintain budget neutrality.
+
+    This chart overlays the CF trajectory with the total RVU pool
+    growth to show the inverse relationship.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(RULES, alt, pl, q):
+    cf_data = pl.DataFrame({
+        "year": list(RULES.keys()),
+        "conversion_factor": [r.conversion_factor for r in RULES.values()],
+    })
+
+    pool_growth = q("""
+        SELECT year, round(sum(work_rvu + non_fac_pe_rvu + mp_rvu), 0) as total_rvu
+        FROM pfs.rvu WHERE mod IS NULL OR mod = ''
+        GROUP BY year ORDER BY year
+    """)
+
+    combined = cf_data.join(pool_growth, on="year", how="inner")
+
+    cf_line = (
+        alt.Chart(combined.to_pandas())
+        .mark_line(point=True, color="#1f77b4")
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("conversion_factor:Q", title="Conversion Factor ($)",
+                     scale=alt.Scale(zero=False)),
+            tooltip=["year", "conversion_factor", "total_rvu"],
+        )
+    )
+
+    rvu_line = (
+        alt.Chart(combined.to_pandas())
+        .mark_line(point=True, color="#d62728", strokeDash=[4, 4])
+        .encode(
+            x=alt.X("year:O"),
+            y=alt.Y("total_rvu:Q", title="Total Unweighted RVU Pool",
+                     scale=alt.Scale(zero=False)),
+        )
+    )
+
+    cf_chart = (
+        alt.layer(cf_line, rvu_line)
+        .resolve_scale(y="independent")
+        .properties(title="Conversion Factor vs. RVU Pool Growth", width=700, height=350)
+    )
+    cf_chart
+    return
+
+
+# ── 6. Per-code dollar impact ─────────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 6. Dollar Impact on Common Services
+
+    How much does a typical office visit, imaging study, or surgical
+    procedure lose when skin sub PE RVUs increase? This table shows
+    the implied payment reduction for commonly billed codes.
+
+    The calculation: if skin sub PE grows by Δ and the total PE pool
+    is P, then each other code's PE is effectively reduced by
+    `code_PE × (Δ / P)`, and payment drops by that × CF.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(RULES, SKIN_CODES, mo, q):
+
+    latest_year = max(RULES.keys())
+    cf = RULES[latest_year].conversion_factor
+
+    # Latest skin PE delta
+    skin_delta = q(f"""
+        WITH yearly AS (
+            SELECT year, sum(non_fac_pe_rvu) as skin_pe
+            FROM pfs.rvu WHERE hcpcs IN {SKIN_CODES} AND (mod IS NULL OR mod = '')
+            GROUP BY year
+        )
+        SELECT skin_pe - LAG(skin_pe) OVER (ORDER BY year) as delta
+        FROM yearly ORDER BY year DESC LIMIT 1
+    """).item()
+
+    total_pe = q(f"""
+        SELECT sum(non_fac_pe_rvu) FROM pfs.rvu
+        WHERE year = {latest_year} AND (mod IS NULL OR mod = '')
+          AND hcpcs NOT IN {SKIN_CODES}
+    """).item()
+
+    tax_rate = skin_delta / total_pe if total_pe else 0
+
+    common_codes = q(f"""
+        SELECT hcpcs, description,
+               non_fac_pe_rvu,
+               work_rvu,
+               mp_rvu,
+               non_fac_pe_rvu + work_rvu + mp_rvu as total_rvu,
+               round(non_fac_pe_rvu * {tax_rate}, 6) as pe_rvu_loss,
+               round(non_fac_pe_rvu * {tax_rate} * {cf}, 4) as dollar_loss
+        FROM pfs.rvu
+        WHERE year = {latest_year}
+          AND (mod IS NULL OR mod = '')
+          AND hcpcs IN ('99213','99214','99215',
+                        '99203','99204','99205',
+                        '27447','27130',
+                        '43239','45380',
+                        '93000','93306',
+                        '77067','74177',
+                        '36415','85025',
+                        '90834','90837',
+                        '17000','11102')
+        ORDER BY dollar_loss
+    """)
+
+    mo.vstack([
+        mo.md(f"""
+        **CY{latest_year}** parameters:
+        - Skin sub PE delta: **{skin_delta:+.2f} RVUs**
+        - Other-code PE pool: **{total_pe:,.1f} RVUs**
+        - Implied tax rate: **{tax_rate*100:.4f}%** of each code's PE
+        - CF: **${cf}**
+        """),
+        common_codes,
+    ])
+    return
+
+
+# ── 7. Cumulative tax since 2015 ──────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 7. Cumulative Structural Tax Since 2015
+
+    The skin sub application codes' total PE RVUs have grown from
+    16.45 (2015) to 22.05 (2026). Under budget neutrality, this
+    5.60 RVU increase must come from somewhere.
+
+    This chart shows the cumulative PE RVU "withdrawn" from the
+    rest of the fee schedule by the growth in skin sub PE.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(SKIN_CODES, alt, q):
+    cumulative = q(f"""
+        WITH yearly AS (
+            SELECT year,
+                   sum(CASE WHEN hcpcs IN {SKIN_CODES} THEN non_fac_pe_rvu ELSE 0 END) as skin_pe,
+                   sum(non_fac_pe_rvu) as total_pe
+            FROM pfs.rvu WHERE mod IS NULL OR mod = ''
+            GROUP BY year
+        )
+        SELECT year,
+               round(skin_pe, 2) as skin_pe,
+               round(skin_pe - FIRST_VALUE(skin_pe) OVER (ORDER BY year), 2) as cumulative_pe_growth,
+               round(total_pe, 1) as total_pe
+        FROM yearly ORDER BY year
+    """)
+
+    cum_chart = (
+        alt.Chart(cumulative.to_pandas())
+        .mark_area(opacity=0.3, color="#d62728")
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("cumulative_pe_growth:Q",
+                     title="Cumulative Skin Sub PE Growth (RVUs above 2015 baseline)"),
+            tooltip=["year", "skin_pe", "cumulative_pe_growth", "total_pe"],
+        )
+    ) + (
+        alt.Chart(cumulative.to_pandas())
+        .mark_line(point=True, color="#d62728")
+        .encode(
+            x="year:O",
+            y="cumulative_pe_growth:Q",
+        )
+    )
+
+    cum_chart_final = cum_chart.properties(
+        title="Cumulative PE RVU Growth — Skin Sub Application Codes vs. 2015 Baseline",
+        width=700, height=300,
+    )
+    cum_chart_final
+    return
+
+
+# ── 8. Key findings ──────────────────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 8. Key Findings
+
+    1. **Structural tax is real but small per-code.** The unweighted RVU
+       share of skin sub application codes is ~0.1% of the pool.
+       Per-code compression is fractions of a cent.
+
+    2. **Volume is the amplifier.** The structural RVU analysis
+       understates the true impact because it ignores utilization.
+       OIG documented explosive volume growth in skin sub claims — when
+       frequency-weighted, these 8 codes consume a much larger share of
+       aggregate spending than their unweighted RVUs suggest.
+
+    3. **PE is the battleground.** Work RVUs for 15271–15278 have been
+       stable (13.46 total for 10 years). PE RVUs grew from 16.45 to
+       22.05 (+34%). The CY2022 clinical labor rate update was a major
+       driver.
+
+    4. **CY2026 reclassification shifts the tax.** Moving skin subs
+       from ASP + 6% (OPPS) to flat $127.28/cm² doesn't directly
+       affect the PFS budget neutrality pool — but it does change
+       the volume incentives that drive utilization of application
+       codes 15271–15278.
+
+    5. **The real cost is in the product, not the application.** For a
+       25cm² wound, the application fee (~$140) is dwarfed by the
+       product cost (often >$5,000). Budget neutrality only governs
+       the application fee; the product cost is outside the PFS.
+    """)
+    return
+
+
+if __name__ == "__main__":
+    app.run()

notebooks/skin_sub_cost_sharing.py

@@ -0,0 +1,512 @@
+import marimo
+
+__generated_with = "0.21.1"
+app = marimo.App(width="medium")
+
+
+@app.cell(hide_code=True)
+def _():
+    import marimo as mo
+
+    return (mo,)
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    # Beneficiary Cost Sharing for Skin Substitutes
+
+    Medicare Part B beneficiaries pay **20% coinsurance** on both the
+    application procedure (PFS) and the product (ASP + 6%). This
+    notebook tracks how cost sharing has changed over time and how
+    the CY2026 reclassification affects out-of-pocket exposure.
+
+    **Cost-sharing components:**
+    - **Application fee**: 20% of PFS carrier locality fee (15271–15278)
+    - **Product cost**: 20% of ASP + 6% payment limit × units
+    - **OPPS copayment**: minimum unadjusted copayment per APC (hospital outpatient)
+    - **Limiting charge**: non-participating providers may charge up to 115% of fee schedule
+
+    **CY2026 change:** All skin subs move to a flat $127.28/cm². OPPS copay
+    drops from ~$165–$366 per code to $25.43 flat. Product coinsurance
+    collapses to $25.43/unit regardless of actual ASP.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _():
+    import altair as alt
+    import polars as pl
+    from conf import connect
+    from pfs.rules import RULES
+
+    con = connect.duckdb()
+
+    def q(sql):
+        return con.execute(sql).pl()
+
+    SKIN_CODES = "('15271','15272','15273','15274','15275','15276','15277','15278')"
+
+    # Part B deductible history (published by CMS annually)
+    DEDUCTIBLES = {
+        2015: 147.00, 2016: 166.00, 2017: 183.00, 2018: 183.00,
+        2019: 185.00, 2020: 198.00, 2021: 203.00, 2022: 233.00,
+        2023: 226.00, 2024: 240.00, 2025: 257.00, 2026: 257.00,
+    }
+
+    return DEDUCTIBLES, RULES, SKIN_CODES, alt, con, pl, q
+
+
+# ── 1. Application fee coinsurance over time ─────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 1. Application Fee Coinsurance by Region
+
+    The beneficiary pays 20% of the PFS-approved amount for the
+    application procedure. This varies by locality and year.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(alt, q):
+    regions = [
+        "MANHATTAN", "REST OF FLORIDA", "REST OF TEXAS",
+        "REST OF CALIFORNIA", "SOUTH CAROLINA",
+    ]
+    region_list = ", ".join(f"'{r}'" for r in regions)
+
+    app_coinsurance = q(f"""
+        SELECT c.year, g.locality_name,
+               c.non_fac_fee,
+               round(c.non_fac_fee * 0.20, 2) as bene_coinsurance,
+               c.non_fac_limiting_charge,
+               round(c.non_fac_limiting_charge * 0.20, 2) as bene_limiting_coinsurance
+        FROM pfs.carrier_locality c
+        JOIN pfs.gpci g ON c.mac = g.mac AND c.locality = g.locality AND c.year = g.year
+        WHERE c.hcpcs = '15271'
+          AND g.locality_name IN ({region_list})
+        ORDER BY c.year, g.locality_name
+    """)
+
+    app_chart = (
+        alt.Chart(app_coinsurance.to_pandas())
+        .mark_line(point=True)
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("bene_coinsurance:Q", title="Beneficiary Coinsurance ($)"),
+            color=alt.Color("locality_name:N", title="Region"),
+            tooltip=["year", "locality_name", "non_fac_fee",
+                     "bene_coinsurance", "non_fac_limiting_charge"],
+        )
+        .properties(
+            title="15271 Application — Beneficiary 20% Coinsurance by Region",
+            width=700, height=350,
+        )
+    )
+    app_chart
+    return
+
+
+# ── 2. Product coinsurance trajectory ────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 2. Product Coinsurance Over Time (Per Unit)
+
+    The beneficiary pays 20% of the ASP + 6% payment limit **per cm²**.
+    For high-ASP products, this adds up fast — a 25cm² application of
+    a $150/cm² product means $750 in coinsurance just for the product.
+
+    The red dashed line shows the CY2026 flat-rate coinsurance:
+    20% × $127.28 = **$25.46/cm²**.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(alt, pl, q):
+    product_coins = q("""
+        SELECT quarter, hcpcs_code, short_description,
+               payment_limit,
+               round(payment_limit * 0.20, 2) as bene_per_unit
+        FROM skin_subs.asp_quarterly
+        WHERE hcpcs_code IN ('Q4101','Q4186','Q4132','Q4116','Q4100')
+        ORDER BY quarter, hcpcs_code
+    """)
+
+    flat_coins = (
+        alt.Chart(pl.DataFrame({"y": [127.28 * 0.20]}).to_pandas())
+        .mark_rule(color="red", strokeDash=[4, 4], strokeWidth=2)
+        .encode(y="y:Q")
+    )
+
+    product_lines = (
+        alt.Chart(product_coins.to_pandas())
+        .mark_line()
+        .encode(
+            x=alt.X("quarter:O", title="Quarter",
+                     axis=alt.Axis(labelAngle=-45, labelFontSize=8)),
+            y=alt.Y("bene_per_unit:Q", title="Beneficiary Coinsurance per cm² ($)"),
+            color=alt.Color("short_description:N", title="Product"),
+            tooltip=["quarter", "hcpcs_code", "short_description",
+                     "payment_limit", "bene_per_unit"],
+        )
+    )
+
+    product_chart = (
+        (product_lines + flat_coins)
+        .properties(
+            title="Product Coinsurance per Unit (20% of ASP + 6%)",
+            width=700, height=400,
+        )
+    )
+    product_chart
+    return
+
+
+# ── 3. Total episode cost sharing ────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 3. Total Episode Cost Sharing (25cm² Wound)
+
+    For a typical 25cm² wound treated in an office setting, the
+    beneficiary's total out-of-pocket is:
+
+    ```
+    total = 20% × application_fee + 20% × (ASP+6% × 25 units)
+    ```
+
+    This chart shows how that total has changed over time for
+    selected products, with the post-2026 flat-rate equivalent shown.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(alt, pl, q):
+    # Use Q1 of each year for annual comparison
+    episode_sharing = q("""
+        WITH asp_annual AS (
+            SELECT CAST(substr(quarter, 1, 4) AS INTEGER) as year,
+                   hcpcs_code, short_description,
+                   payment_limit
+            FROM skin_subs.asp_quarterly
+            WHERE substr(quarter, 6, 2) = 'Q1'
+              AND hcpcs_code IN ('Q4101','Q4186','Q4132','Q4116')
+        ),
+        app_fee AS (
+            SELECT c.year, avg(c.non_fac_fee) as avg_app_fee
+            FROM pfs.carrier_locality c
+            WHERE c.hcpcs = '15271'
+            GROUP BY c.year
+        )
+        SELECT a.year, a.hcpcs_code, a.short_description,
+               round(a.payment_limit, 2) as asp_per_unit,
+               round(f.avg_app_fee, 2) as avg_app_fee,
+               round(f.avg_app_fee * 0.20, 2) as app_coinsurance,
+               round(a.payment_limit * 25 * 0.20, 2) as product_coinsurance_25cm,
+               round(f.avg_app_fee * 0.20 + a.payment_limit * 25 * 0.20, 2) as total_bene_cost
+        FROM asp_annual a
+        LEFT JOIN app_fee f ON a.year = f.year
+        ORDER BY a.year, a.hcpcs_code
+    """)
+
+    flat_episode = 127.28 * 25 * 0.20  # $636.40 product + ~$28 application
+
+    flat_line = (
+        alt.Chart(pl.DataFrame({"y": [flat_episode]}).to_pandas())
+        .mark_rule(color="red", strokeDash=[4, 4], strokeWidth=2)
+        .encode(y="y:Q")
+    )
+
+    episode_lines = (
+        alt.Chart(episode_sharing.to_pandas())
+        .mark_line(point=True)
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("total_bene_cost:Q",
+                     title="Beneficiary Total Cost Sharing ($)"),
+            color=alt.Color("short_description:N", title="Product"),
+            tooltip=["year", "hcpcs_code", "short_description",
+                     "asp_per_unit", "app_coinsurance",
+                     "product_coinsurance_25cm", "total_bene_cost"],
+        )
+    )
+
+    episode_chart = (
+        (episode_lines + flat_line)
+        .properties(
+            title="Total Beneficiary Cost Sharing — 25cm² Wound Episode",
+            width=700, height=400,
+        )
+    )
+    episode_chart
+    return
+
+
+# ── 4. OPPS copayment collapse ───────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 4. OPPS Copayment Collapse in CY2026
+
+    In the hospital outpatient setting (HOPD), beneficiaries pay a
+    minimum unadjusted copayment per APC. For skin substitutes:
+
+    - **2021–2025**: Copays ranged from $105–$366 depending on the
+      product's APC assignment (high-cost vs. low-cost categories)
+    - **CY2026**: All products reclassified to a single flat rate
+      with copay of **$25.43** — a >85% reduction
+
+    This is the most dramatic beneficiary cost-sharing change in
+    the CY2026 reclassification.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(alt, q):
+    opps_copay = q("""
+        SELECT year,
+               count(*) as products,
+               round(avg(minimum_unadjusted_copayment), 2) as avg_copay,
+               round(min(minimum_unadjusted_copayment), 2) as min_copay,
+               round(max(minimum_unadjusted_copayment), 2) as max_copay,
+               round(max(minimum_unadjusted_copayment) -
+                     min(minimum_unadjusted_copayment), 2) as copay_spread
+        FROM opps.skin_sub_addendum_b
+        WHERE minimum_unadjusted_copayment > 0
+        GROUP BY year ORDER BY year
+    """)
+
+    copay_chart = (
+        alt.Chart(opps_copay.to_pandas())
+        .mark_bar()
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("avg_copay:Q", title="Average OPPS Copayment ($)"),
+            color=alt.condition(
+                alt.datum.year == 2026,
+                alt.value("#2ca02c"),
+                alt.value("#1f77b4"),
+            ),
+            tooltip=["year", "products", "avg_copay", "min_copay",
+                     "max_copay", "copay_spread"],
+        )
+        .properties(
+            title="OPPS Minimum Unadjusted Copayment — Skin Substitutes",
+            width=700, height=300,
+        )
+    )
+
+    mo.vstack([copay_chart, opps_copay])
+    return
+
+
+# ── 5. Winners and losers: beneficiary perspective ────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 5. Beneficiary Impact: Who Pays More, Who Pays Less?
+
+    Under the flat rate, beneficiaries using **high-ASP products**
+    (EpiFix, GrafixCore) see their per-unit coinsurance drop
+    dramatically. Those using **low-ASP products** (Apligraf at
+    ~$30/unit) see coinsurance rise from ~$6 to ~$25.
+
+    This is a wealth transfer: beneficiaries who previously used
+    expensive products benefit; those with cheaper products pay more.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(alt, q):
+    bene_impact = q("""
+        SELECT hcpcs_code, short_description,
+               payment_limit as current_asp_payment,
+               round(payment_limit * 0.20, 2) as current_bene_per_unit,
+               127.28 as flat_rate,
+               round(127.28 * 0.20, 2) as flat_bene_per_unit,
+               round(127.28 * 0.20 - payment_limit * 0.20, 2) as bene_delta_per_unit,
+               CASE
+                   WHEN payment_limit * 0.20 > 127.28 * 0.20 THEN 'bene saves'
+                   WHEN payment_limit * 0.20 < 127.28 * 0.20 THEN 'bene pays more'
+                   ELSE 'neutral'
+               END as bene_impact
+        FROM skin_subs.asp_quarterly
+        WHERE quarter = (SELECT max(quarter) FROM skin_subs.asp_quarterly)
+        ORDER BY bene_delta_per_unit
+    """)
+
+    bene_chart = (
+        alt.Chart(bene_impact.to_pandas())
+        .mark_bar()
+        .encode(
+            x=alt.X("bene_delta_per_unit:Q",
+                     title="Change in Beneficiary Coinsurance per Unit ($)"),
+            y=alt.Y("short_description:N", title="", sort="x",
+                     axis=alt.Axis(labelLimit=300)),
+            color=alt.Color("bene_impact:N",
+                            scale=alt.Scale(
+                                domain=["bene saves", "bene pays more", "neutral"],
+                                range=["#2ca02c", "#d62728", "#7f7f7f"],
+                            ),
+                            title="Impact"),
+            tooltip=["hcpcs_code", "short_description",
+                     "current_bene_per_unit", "flat_bene_per_unit",
+                     "bene_delta_per_unit"],
+        )
+        .properties(
+            title="CY2026 Beneficiary Coinsurance Change per Unit",
+            width=700,
+        )
+    )
+    bene_chart
+    return
+
+
+# ── 6. Episode scenario comparison ───────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 6. Episode Scenarios: Before vs. After CY2026
+
+    Three wound scenarios showing the complete beneficiary cost
+    breakdown before and after the flat-rate reclassification.
+    All assume office setting (POS 11), participating provider,
+    national average application fee.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(pl, q):
+    # Get average national app fee for 2025
+    avg_fee = q("""
+        SELECT round(avg(non_fac_fee), 2) as avg_fee
+        FROM pfs.carrier_locality
+        WHERE year = 2025 AND hcpcs = '15271'
+    """).item()
+
+    scenarios = pl.DataFrame({
+        "scenario": [
+            "Small wound (10cm²) — Q4101 Apligraf",
+            "Medium wound (25cm²) — Q4186 EpiFix",
+            "Large wound (50cm²) — Q4132 GrafixCore",
+        ],
+        "units": [10, 25, 50],
+        "asp_per_unit": [30.23, 151.17, 106.70],
+        "product_name": ["Apligraf", "EpiFix", "GrafixCore"],
+    }).with_columns(
+        # Pre-2026: ASP + 6%
+        (pl.col("asp_per_unit") * pl.col("units") * 0.20).round(2).alias("pre_product_coins"),
+        pl.lit(avg_fee * 0.20).round(2).alias("pre_app_coins"),
+        # Post-2026: flat $127.28
+        (pl.lit(127.28) * pl.col("units") * 0.20).round(2).alias("post_product_coins"),
+        pl.lit(avg_fee * 0.20).round(2).alias("post_app_coins"),
+    ).with_columns(
+        (pl.col("pre_product_coins") + pl.col("pre_app_coins")).alias("pre_total"),
+        (pl.col("post_product_coins") + pl.col("post_app_coins")).alias("post_total"),
+    ).with_columns(
+        (pl.col("post_total") - pl.col("pre_total")).round(2).alias("delta"),
+    )
+    scenarios
+    return
+
+
+# ── 7. Part B deductible context ─────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 7. Part B Deductible Context
+
+    Before coinsurance applies, beneficiaries must meet the annual
+    Part B deductible. A single skin substitute episode can exceed
+    the entire deductible, meaning the full coinsurance amount is
+    additional out-of-pocket cost for most beneficiaries.
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _(DEDUCTIBLES, alt, pl):
+    deductible_df = pl.DataFrame({
+        "year": list(DEDUCTIBLES.keys()),
+        "deductible": list(DEDUCTIBLES.values()),
+    })
+
+    ded_chart = (
+        alt.Chart(deductible_df.to_pandas())
+        .mark_bar(color="#ff7f0e")
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("deductible:Q", title="Annual Part B Deductible ($)"),
+            tooltip=["year", "deductible"],
+        )
+        .properties(
+            title="Medicare Part B Annual Deductible",
+            width=700, height=250,
+        )
+    )
+    ded_chart
+    return
+
+
+# ── 8. Key takeaways ─────────────────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 8. Key Takeaways
+
+    1. **OPPS copay drops >85% in CY2026.** From $105–$366 per code
+       to $25.43 flat. This is the single largest beneficiary-facing
+       change in the reclassification.
+
+    2. **High-ASP product users save significantly.** For EpiFix
+       (Q4186, ~$151/unit), coinsurance per cm² drops from ~$30 to
+       ~$25 — but for a 25cm² wound, the total product coinsurance
+       drops from ~$755 to ~$636 ($119 savings).
+
+    3. **Low-ASP product users pay more.** For Apligraf (Q4101,
+       ~$30/unit), per-unit coinsurance rises from ~$6 to ~$25.
+       A 25cm² wound goes from ~$30 to ~$636 in product coinsurance.
+
+    4. **Application fee coinsurance is stable.** The PFS component
+       (~$27–$35 for 15271 depending on locality) is a small fraction
+       of total cost sharing and relatively stable over time.
+
+    5. **Deductible is a floor, not a ceiling.** The $257 Part B
+       deductible (2025–2026) is typically exceeded by a single skin
+       sub episode, so coinsurance applies to the full amount.
+
+    6. **No Medigap/supplement analysis.** Most beneficiaries have
+       supplemental coverage (Medigap, employer, Medicaid dual) that
+       covers the 20% coinsurance. Actual out-of-pocket may be lower
+       depending on coverage type.
+    """)
+    return
+
+
+if __name__ == "__main__":
+    app.run()