stack/notebooks/skin_sub_pricing.py

import marimo

__generated_with = "0.23.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("""
    # 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()
    # OPPS reference data lives in the DuckLake lakehouse (M5, #514);
    # the monolith's opps schema is a deprecated mirror.
    lake = connect.ducklake()

    def q(sql):
        return con.execute(sql).pl()

    def ql(sql):
        """Query the lake (OPPS reference tables)."""
        return lake.execute(sql).pl()

    return alt, con, lake, pl, q, ql


@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(hide_code=True)
def _(alt, ql):
    rvu_ts = ql("""
        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


@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(hide_code=True)
def _(alt, 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


@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(hide_code=True)
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


@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(hide_code=True)
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


@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(hide_code=True)
def _(ql):
    opps_ts = ql("""
        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


@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(hide_code=True)
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


@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(hide_code=True)
def _(alt, 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


@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(hide_code=True)
def _(ql):
    fee_summary = ql("""
        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()