feat: skin sub pricing notebook — RVU, locality fees, and ASP over time

New marimo notebook tracking skin substitute economics across payment
systems and years:

1. Application code RVUs (15271-15278) over time
2. Locality-based PFS fees by region (MAC-aware carrier join)
3. Regional comparison for 15271 across 10 high-volume areas
4. Quarterly ASP payment limits for top 15 products with $127.28
   flat-rate reference line
5. OPPS payment method timeline (pass-through → high/low → flat rate)
6. Total episode cost: product ASP + application fee for 25cm² wound
7. Flat-rate impact analysis: winners and losers under CY2026
8. Summary table: fee range spread across all localities by year

Data sources: pfs.rvu, pfs.carrier_locality, pfs.gpci,
skin_subs.asp_quarterly, opps.skin_sub_addendum_b

refs #246

notebooks/skin_sub_pricing.py

@@ -0,0 +1,502 @@
+import marimo
+
+__generated_with = "0.20.2"
+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("""
+    # Skin Substitute Pricing Over Time
+
+    Tracks skin substitute application code RVUs and locality-based
+    payment amounts (PFS) alongside quarterly ASP drug pricing, layered
+    over the OPPS payment method timeline (pass-through → high/low →
+    flat rate).
+
+    **Data sources:**
+    - `pfs.rvu` — Work, PE, MP RVUs for application codes 15271–15278
+    - `pfs.carrier_locality` — CMS-published locality fees (ground truth)
+    - `pfs.gpci` — Geographic practice cost indices by MAC/locality
+    - `skin_subs.asp_quarterly` — Quarterly ASP + 6% payment limits
+    - `opps.skin_sub_addendum_b` — OPPS status indicators and APC rates
+    """)
+    return
+
+
+@app.cell(hide_code=True)
+def _():
+    import altair as alt
+    import polars as pl
+    from conf import connect
+
+    con = connect.duckdb()
+
+    def q(sql):
+        return con.execute(sql).pl()
+
+    return alt, con, pl, q
+
+
+# ── 1. Application code RVUs over time ──────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 1. Application Code RVUs Over Time
+
+    CPT 15271–15278 are the procedure codes for applying skin
+    substitutes. Their RVUs determine the surgeon's payment
+    independent of the product cost.
+    """)
+    return
+
+
+@app.cell
+def _(alt, q):
+    rvu_ts = q("""
+        SELECT year, hcpcs,
+               CASE hcpcs
+                   WHEN '15271' THEN 'Trunk/limbs <100cm² (initial)'
+                   WHEN '15272' THEN 'Trunk/limbs <100cm² (add-on)'
+                   WHEN '15273' THEN 'Trunk/limbs ≥100cm² (initial)'
+                   WHEN '15274' THEN 'Trunk/limbs ≥100cm² (add-on)'
+                   WHEN '15275' THEN 'Face/hands/feet <100cm² (initial)'
+                   WHEN '15276' THEN 'Face/hands/feet <100cm² (add-on)'
+                   WHEN '15277' THEN 'Face/hands/feet ≥100cm² (initial)'
+                   WHEN '15278' THEN 'Face/hands/feet ≥100cm² (add-on)'
+               END as description,
+               work_rvu, non_fac_pe_rvu, fac_pe_rvu, mp_rvu,
+               work_rvu + non_fac_pe_rvu + mp_rvu as non_fac_total,
+               work_rvu + fac_pe_rvu + mp_rvu as fac_total
+        FROM pfs.rvu
+        WHERE hcpcs IN ('15271','15272','15273','15274',
+                        '15275','15276','15277','15278')
+        ORDER BY year, hcpcs
+    """)
+
+    rvu_chart = (
+        alt.Chart(rvu_ts.to_pandas())
+        .mark_line(point=True)
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("non_fac_total:Q", title="Total Non-Facility RVUs"),
+            color=alt.Color("description:N", title="Code"),
+            tooltip=["year", "hcpcs", "description",
+                     "work_rvu", "non_fac_pe_rvu", "mp_rvu", "non_fac_total"],
+        )
+        .properties(title="Application Code RVUs (Non-Facility)", width=700, height=400)
+    )
+    rvu_chart
+    return (rvu_ts,)
+
+
+# ── 2. Locality-based payment for application codes ─────────────────
+
+
+@app.cell(hide_code=True)
+def _(con, mo):
+    _regions = con.execute("""
+        SELECT DISTINCT locality_name
+        FROM pfs.gpci WHERE year = 2025
+        ORDER BY locality_name
+    """).fetchall()
+
+    region_picker = mo.ui.dropdown(
+        options={r[0]: r[0] for r in _regions},
+        value="MANHATTAN",
+        label="Region",
+    )
+    mo.md(f"""
+    ## 2. Application Code Payment by Region
+
+    CMS-published carrier locality fees for skin sub application codes.
+    Select a region to see how payments vary geographically over time.
+
+    {region_picker}
+    """)
+    return (region_picker,)
+
+
+@app.cell
+def _(alt, pl, q, region_picker):
+    _region = region_picker.value
+
+    carrier_ts = q(f"""
+        SELECT c.year, c.hcpcs,
+               CASE c.hcpcs
+                   WHEN '15271' THEN 'Trunk <100cm²'
+                   WHEN '15272' THEN 'Trunk <100cm² add-on'
+                   WHEN '15273' THEN 'Trunk ≥100cm²'
+                   WHEN '15274' THEN 'Trunk ≥100cm² add-on'
+                   WHEN '15275' THEN 'Face/hands <100cm²'
+                   WHEN '15276' THEN 'Face/hands <100cm² add-on'
+                   WHEN '15277' THEN 'Face/hands ≥100cm²'
+                   WHEN '15278' THEN 'Face/hands ≥100cm² add-on'
+               END as description,
+               c.non_fac_fee, c.fac_fee,
+               g.locality_name
+        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 IN ('15271','15272','15273','15274',
+                          '15275','15276','15277','15278')
+          AND g.locality_name = '{_region}'
+        ORDER BY c.year, c.hcpcs
+    """)
+
+    carrier_chart = (
+        alt.Chart(carrier_ts.to_pandas())
+        .mark_line(point=True)
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("non_fac_fee:Q", title="Non-Facility Fee ($)"),
+            color=alt.Color("description:N", title="Code"),
+            tooltip=["year", "hcpcs", "description", "non_fac_fee",
+                     "fac_fee", "locality_name"],
+        )
+        .properties(
+            title=f"Application Code Fees — {_region}",
+            width=700, height=400,
+        )
+    )
+    carrier_chart
+    return (carrier_ts,)
+
+
+# ── 3. Regional comparison for a single code ────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 3. Regional Comparison — 15271 (Primary Application)
+
+    Non-facility fee for 15271 (trunk/limbs, <100cm², initial
+    application) across selected high-volume regions.
+    """)
+    return
+
+
+@app.cell
+def _(alt, q):
+    regions = [
+        "MANHATTAN", "ALASKA*", "REST OF FLORIDA", "REST OF TEXAS",
+        "REST OF CALIFORNIA", "CHICAGO", "DETROIT", "ATLANTA",
+        "HAWAII, GUAM", "SOUTH CAROLINA",
+    ]
+    region_list = ", ".join(f"'{r}'" for r in regions)
+
+    regional_15271 = q(f"""
+        SELECT c.year, g.locality_name,
+               c.non_fac_fee
+        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
+    """)
+
+    regional_chart = (
+        alt.Chart(regional_15271.to_pandas())
+        .mark_line(point=True)
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("non_fac_fee:Q", title="Non-Facility Fee ($)"),
+            color=alt.Color("locality_name:N", title="Region"),
+            tooltip=["year", "locality_name", "non_fac_fee"],
+        )
+        .properties(
+            title="15271 Non-Facility Fee — Regional Comparison",
+            width=700, height=400,
+        )
+    )
+    regional_chart
+    return
+
+
+# ── 4. ASP product pricing over time ────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 4. ASP Drug Pricing Over Time
+
+    Quarterly ASP + 6% payment limits for skin substitute products.
+    The flat horizontal line at **$127.28** marks the CY2026 flat rate —
+    products above it lose revenue, products below it gain.
+    """)
+    return
+
+
+@app.cell
+def _(alt, pl, q):
+    # Top products by volume (most quarters of data)
+    top_products = q("""
+        SELECT hcpcs_code, short_description,
+               count(*) as quarters
+        FROM skin_subs.asp_quarterly
+        GROUP BY hcpcs_code, short_description
+        HAVING count(*) >= 8
+        ORDER BY quarters DESC
+        LIMIT 15
+    """)
+
+    top_codes = top_products.select("hcpcs_code").to_series().to_list()
+    code_list = ", ".join(f"'{c}'" for c in top_codes)
+
+    asp_ts = q(f"""
+        SELECT quarter, hcpcs_code, short_description,
+               payment_limit, asp_per_unit
+        FROM skin_subs.asp_quarterly
+        WHERE hcpcs_code IN ({code_list})
+        ORDER BY quarter, hcpcs_code
+    """)
+
+    flat_rate_rule = (
+        alt.Chart(pl.DataFrame({"y": [127.28]}).to_pandas())
+        .mark_rule(color="red", strokeDash=[4, 4], strokeWidth=2)
+        .encode(y="y:Q")
+    )
+
+    asp_lines = (
+        alt.Chart(asp_ts.to_pandas())
+        .mark_line()
+        .encode(
+            x=alt.X("quarter:O", title="Quarter",
+                     axis=alt.Axis(labelAngle=-45, labelFontSize=8)),
+            y=alt.Y("payment_limit:Q", title="Payment Limit (ASP + 6%) $",
+                     scale=alt.Scale(domainMax=800)),
+            color=alt.Color("short_description:N", title="Product"),
+            tooltip=["quarter", "hcpcs_code", "short_description",
+                     "payment_limit", "asp_per_unit"],
+        )
+    )
+
+    asp_chart = (
+        (asp_lines + flat_rate_rule)
+        .properties(
+            title="ASP Quarterly Payment Limits — Top 15 Products",
+            width=700, height=450,
+        )
+    )
+    asp_chart
+    return
+
+
+# ── 5. OPPS payment method timeline ─────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 5. OPPS Payment Method Timeline
+
+    How CMS has paid for skin substitutes in the outpatient setting:
+
+    | Period | Method | Detail |
+    |--------|--------|--------|
+    | Pre-2023 | Pass-through | ASP + 6% per unit |
+    | 2023–2025 | High/low split | High-cost: pass-through; Low-cost: packaged into APC |
+    | 2026+ | Flat rate | $127.28/cm² for all products |
+
+    The table below shows OPPS status indicators and payment rates
+    for skin sub codes over time.
+    """)
+    return
+
+
+@app.cell
+def _(q):
+    opps_ts = q("""
+        SELECT year, hcpcs, short_description,
+               status_indicator, apc, payment_rate
+        FROM opps.skin_sub_addendum_b
+        WHERE payment_rate IS NOT NULL
+          AND payment_rate > 0
+        ORDER BY year, hcpcs
+    """)
+    opps_ts
+    return (opps_ts,)
+
+
+# ── 6. Combined view: ASP + application fee ─────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 6. Total Episode Cost: Product + Application
+
+    For a typical 25cm² wound, the total Medicare payment is:
+
+    - **Product cost**: ASP + 6% payment limit × 25 units
+    - **Application fee**: 15271 non-facility fee (PFS, locality-specific)
+    - **Post-2026**: $127.28 × 25 = $3,182 flat + application fee
+
+    This chart shows the combined cost trend for selected products.
+    """)
+    return
+
+
+@app.cell
+def _(alt, pl, q, region_picker):
+    _region = region_picker.value
+
+    # Get 15271 fee for the selected region, all years
+    app_fee = q(f"""
+        SELECT c.year, c.non_fac_fee as application_fee
+        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 = '{_region}'
+    """)
+
+    # Get annual ASP for top products (use Q1 of each year)
+    episode_data = q("""
+        SELECT
+            CAST(substr(quarter, 1, 4) AS INTEGER) as year,
+            hcpcs_code, short_description,
+            payment_limit,
+            payment_limit * 25 as product_cost_25cm2
+        FROM skin_subs.asp_quarterly
+        WHERE substr(quarter, 6, 2) = 'Q1'
+          AND hcpcs_code IN ('Q4101','Q4186','Q4132','Q4116','Q4100')
+        ORDER BY year, hcpcs_code
+    """)
+
+    # Join to get total episode cost
+    combined = (
+        episode_data.join(app_fee, on="year", how="left")
+        .with_columns(
+            (pl.col("product_cost_25cm2") + pl.col("application_fee").fill_null(0))
+            .alias("total_episode_cost")
+        )
+    )
+
+    flat_line = (
+        alt.Chart(pl.DataFrame({"y": [127.28 * 25]}).to_pandas())
+        .mark_rule(color="red", strokeDash=[4, 4], strokeWidth=2)
+        .encode(y="y:Q")
+    )
+
+    episode_lines = (
+        alt.Chart(combined.to_pandas())
+        .mark_line(point=True)
+        .encode(
+            x=alt.X("year:O", title="Year"),
+            y=alt.Y("total_episode_cost:Q", title="Total Episode Cost ($)"),
+            color=alt.Color("short_description:N", title="Product"),
+            tooltip=["year", "hcpcs_code", "short_description",
+                     "product_cost_25cm2", "application_fee",
+                     "total_episode_cost"],
+        )
+    )
+
+    episode_chart = (
+        (episode_lines + flat_line)
+        .properties(
+            title=f"Total Episode Cost (25cm² wound) — {_region}",
+            width=700, height=400,
+        )
+    )
+    episode_chart
+    return
+
+
+# ── 7. Flat-rate impact analysis ────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 7. Flat-Rate Impact: Winners and Losers
+
+    Products with ASP + 6% above $127.28 lose revenue under the 2026
+    flat rate; those below gain. This chart shows the latest quarter's
+    payment limit vs. the flat rate.
+    """)
+    return
+
+
+@app.cell
+def _(alt, pl, q):
+    impact = q("""
+        SELECT hcpcs_code, short_description,
+               payment_limit,
+               127.28 as flat_rate,
+               payment_limit - 127.28 as delta,
+               CASE
+                   WHEN payment_limit > 127.28 THEN 'loses'
+                   WHEN payment_limit < 127.28 THEN 'gains'
+                   ELSE 'neutral'
+               END as impact
+        FROM skin_subs.asp_quarterly
+        WHERE quarter = (SELECT max(quarter) FROM skin_subs.asp_quarterly)
+        ORDER BY delta DESC
+    """)
+
+    impact_chart = (
+        alt.Chart(impact.to_pandas())
+        .mark_bar()
+        .encode(
+            x=alt.X("delta:Q", title="Payment Limit − $127.28 Flat Rate"),
+            y=alt.Y("short_description:N", title="", sort="-x",
+                     axis=alt.Axis(labelLimit=300)),
+            color=alt.Color("impact:N",
+                            scale=alt.Scale(
+                                domain=["loses", "gains", "neutral"],
+                                range=["#d62728", "#2ca02c", "#7f7f7f"],
+                            ),
+                            title="Impact"),
+            tooltip=["hcpcs_code", "short_description",
+                     "payment_limit", "flat_rate", "delta"],
+        )
+        .properties(title="CY2026 Flat-Rate Impact by Product", width=700)
+    )
+    impact_chart
+    return
+
+
+# ── 8. Summary table ────────────────────────────────────────────────
+
+
+@app.cell(hide_code=True)
+def _(mo):
+    mo.md("""
+    ## 8. Summary: Application Code Fee Range by Year
+
+    Min, median, and max non-facility fee across all MAC/locality
+    combinations for application code 15271.
+    """)
+    return
+
+
+@app.cell
+def _(q):
+    fee_summary = q("""
+        SELECT c.year,
+               count(*) as localities,
+               round(min(c.non_fac_fee), 2) as min_fee,
+               round(percentile_cont(0.5) WITHIN GROUP (ORDER BY c.non_fac_fee), 2) as median_fee,
+               round(max(c.non_fac_fee), 2) as max_fee,
+               round(max(c.non_fac_fee) - min(c.non_fac_fee), 2) as fee_spread
+        FROM pfs.carrier_locality c
+        WHERE c.hcpcs = '15271'
+        GROUP BY c.year
+        ORDER BY c.year
+    """)
+    fee_summary
+    return
+
+
+if __name__ == "__main__":
+    app.run()