stack/notebooks/skin_subs_explorer.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 Substitutes Analysis

    Predatory pricing, distribution, and utilization analysis of skin
    substitute / cellular tissue products (CTPs) in the Medicare program.

    **Hypothesis:** Certain distributors, products, geographies, and care
    settings are softer or harder targets for predatory use and pricing.

    Navigate sections using the tabs below.
    """)
    return


@app.cell(hide_code=True)
def _():
    import altair as alt
    import polars as pl

    from conf import connect

    con = connect.duckdb()
    return alt, con, pl


@app.cell(hide_code=True)
def _(mo):
    tabs = mo.ui.tabs(
        {
            "ASP Pricing": mo.md("Loading..."),
            "Market": mo.md("Loading..."),
            "Geographic": mo.md("Loading..."),
            "Setting": mo.md("Loading..."),
            "Anomaly": mo.md("Loading..."),
            "Evidence": mo.md("Loading..."),
        }
    )
    tabs
    return


@app.cell(hide_code=True)
def _(con, mo, pl):
    asp = pl.from_pandas(
        con.execute("""
            SELECT quarter, hcpcs_code, product_name, manufacturer, category,
                   asp_per_unit, payment_limit, entry_quarter, total_quarters,
                   qoq_pct_change, cumulative_pct_change,
                   flat_rate_impact, anomalous_spike
            FROM skin_subs.asp_trajectories
            ORDER BY quarter
        """).fetchdf()
    )

    products = sorted(asp["hcpcs_code"].unique().to_list())
    product_select = mo.ui.multiselect(
        products, label="Products", value=["Q4101", "Q4186", "Q4121"]
    )
    product_select
    return asp, product_select


@app.cell(hide_code=True)
def _(alt, asp, mo, pl, product_select):
    mo.stop(not product_select.value, mo.md("Select products above."))
    filtered = asp.filter(pl.col("hcpcs_code").is_in(product_select.value))
    labels = (
        filtered.select(["hcpcs_code", "product_name"])
        .unique()
        .with_columns(
            (pl.col("hcpcs_code") + " " + pl.col("product_name").fill_null("")).alias(
                "label"
            )
        )
    )
    chart_data = filtered.join(labels, on="hcpcs_code").to_pandas()

    chart = (
        alt.Chart(chart_data)
        .mark_line(point=True)
        .encode(
            x=alt.X("quarter:N", title="Quarter", axis=alt.Axis(labelAngle=-45)),
            y=alt.Y("payment_limit:Q", title="Payment Limit ($/cm²)"),
            color=alt.Color("label:N", title="Product"),
            tooltip=[
                "quarter",
                "hcpcs_code",
                "product_name",
                "manufacturer",
                "asp_per_unit",
                "payment_limit",
                "qoq_pct_change",
            ],
        )
        .properties(width=700, height=400, title="ASP Payment Limit Over Time")
    )

    flat_rate_line = (
        alt.Chart(chart_data.head(1))
        .mark_rule(color="red", strokeDash=[4, 4])
        .encode(y=alt.datum(127.28))
    )

    flat_label = (
        alt.Chart(chart_data.head(1))
        .mark_text(align="left", dx=5, dy=-5, color="red", fontSize=10)
        .encode(y=alt.datum(127.28), text=alt.datum("$127.28 flat rate (Jan 2026)"))
    )

    chart + flat_rate_line + flat_label
    return


@app.cell(hide_code=True)
def _(alt, con, mo, pl):
    mfr = pl.from_pandas(
        con.execute("""
            SELECT manufacturer, product_count, active_products,
                   total_revenue, avg_asp, categories
            FROM skin_subs.manufacturers
            WHERE total_revenue IS NOT NULL
            ORDER BY total_revenue DESC LIMIT 15
        """).fetchdf()
    )

    mfr_chart = (
        alt.Chart(mfr.to_pandas())
        .mark_bar()
        .encode(
            x=alt.X("total_revenue:Q", title="Total Revenue ($)"),
            y=alt.Y("manufacturer:N", sort="-x", title=None),
            color=alt.Color(
                "avg_asp:Q", scale=alt.Scale(scheme="reds"), title="Avg ASP ($/cm²)"
            ),
            tooltip=[
                "manufacturer",
                "product_count",
                "active_products",
                "total_revenue",
                "avg_asp",
                "categories",
            ],
        )
        .properties(width=600, height=350, title="Top Manufacturers by Revenue")
    )

    # Category breakdown
    cat = pl.from_pandas(
        con.execute("""
            SELECT category, count(*) as products,
                   round(avg(latest_asp), 2) as avg_asp,
                   sum(claim_lines) as claims
            FROM skin_subs.products_enriched
            WHERE category IS NOT NULL
            GROUP BY category ORDER BY products DESC LIMIT 12
        """).fetchdf()
    )

    cat_chart = (
        alt.Chart(cat.to_pandas())
        .mark_bar()
        .encode(
            x=alt.X("products:Q", title="Product Count"),
            y=alt.Y("category:N", sort="-x", title=None),
            color=alt.Color(
                "avg_asp:Q", scale=alt.Scale(scheme="blues"), title="Avg ASP"
            ),
            tooltip=["category", "products", "avg_asp", "claims"],
        )
        .properties(width=600, height=300, title="Products by Category")
    )

    mo.vstack([mfr_chart, cat_chart])
    return


@app.cell(hide_code=True)
def _(alt, con, mo, pl):
    geo = pl.from_pandas(
        con.execute("""
            SELECT state, mac_jurisdiction, claim_lines, unique_benes,
                   unique_providers, total_paid, paid_per_bene,
                   pct_office, pct_podiatry
            FROM skin_subs.geographic_summary
            ORDER BY total_paid DESC
        """).fetchdf()
    )

    geo_chart = (
        alt.Chart(geo.to_pandas())
        .mark_bar()
        .encode(
            x=alt.X("paid_per_bene:Q", title="Paid per Beneficiary ($)"),
            y=alt.Y("state:N", sort="-x", title=None),
            color=alt.Color("mac_jurisdiction:N", title="MAC Jurisdiction"),
            tooltip=[
                "state",
                "mac_jurisdiction",
                "unique_benes",
                "total_paid",
                "paid_per_bene",
                "pct_office",
                "pct_podiatry",
            ],
        )
        .properties(width=600, height=400, title="Spend per Beneficiary by State")
    )

    mac = pl.from_pandas(
        con.execute(
            "SELECT * FROM skin_subs.mac_summary ORDER BY total_paid DESC"
        ).fetchdf()
    )

    mac_chart = (
        alt.Chart(mac.to_pandas())
        .mark_bar()
        .encode(
            x=alt.X("total_paid:Q", title="Total Paid ($)"),
            y=alt.Y("mac_jurisdiction:N", sort="-x", title=None),
            color=alt.Color(
                "avg_pct_podiatry:Q",
                scale=alt.Scale(scheme="oranges"),
                title="% Podiatry",
            ),
            tooltip=[
                "mac_jurisdiction",
                "states",
                "claim_lines",
                "total_paid",
                "avg_paid_per_bene",
                "avg_pct_office",
                "avg_pct_podiatry",
            ],
        )
        .properties(width=600, height=250, title="Total Paid by MAC Jurisdiction")
    )

    mo.vstack([geo_chart, mac_chart])
    return


@app.cell(hide_code=True)
def _(alt, con, mo, pl):
    setting = pl.from_pandas(
        con.execute("""
            SELECT setting, provider_specialty,
                   sum(claim_lines) as lines, round(sum(total_paid), 2) as paid
            FROM skin_subs.setting_analysis
            GROUP BY setting, provider_specialty
            ORDER BY paid DESC
        """).fetchdf()
    )

    setting_chart = (
        alt.Chart(setting.to_pandas())
        .mark_bar()
        .encode(
            x=alt.X("paid:Q", title="Total Paid ($)"),
            y=alt.Y("provider_specialty:N", sort="-x", title=None),
            color=alt.Color("setting:N", title="Setting"),
            tooltip=["setting", "provider_specialty", "lines", "paid"],
        )
        .properties(width=600, height=300, title="Spend by Specialty and Setting")
    )

    setting_chart
    return


@app.cell(hide_code=True)
def _(alt, con, mo, pl):
    risk = pl.from_pandas(
        con.execute("""
            SELECT rendering_npi, provider_specialty, state, primary_setting,
                   total_paid, paid_per_patient, unique_patients,
                   volume_zscore, intensity_zscore, composite_risk,
                   risk_tier, fraud_pattern_match, benford_chi2
            FROM skin_subs.provider_risk_profile
            ORDER BY composite_risk DESC
        """).fetchdf()
    )

    scatter = (
        alt.Chart(risk.to_pandas())
        .mark_circle(size=60)
        .encode(
            x=alt.X("volume_zscore:Q", title="Volume Z-Score (vs specialty peers)"),
            y=alt.Y("intensity_zscore:Q", title="Intensity Z-Score ($/patient)"),
            color=alt.Color(
                "risk_tier:N",
                scale=alt.Scale(
                    domain=["critical", "high", "moderate", "low"],
                    range=["#d62728", "#ff7f0e", "#bcbd22", "#2ca02c"],
                ),
                title="Risk Tier",
            ),
            size=alt.Size("total_paid:Q", title="Total Paid"),
            shape=alt.Shape("primary_setting:N", title="Setting"),
            tooltip=[
                "rendering_npi",
                "provider_specialty",
                "state",
                "primary_setting",
                "total_paid",
                "paid_per_patient",
                "composite_risk",
                "risk_tier",
                "fraud_pattern_match",
            ],
        )
        .properties(width=650, height=450, title="Provider Risk: Volume vs Intensity")
    )

    # Risk tier distribution
    tier = risk.group_by("risk_tier").agg(
        pl.col("rendering_npi").count().alias("n"),
        pl.col("total_paid").mean().round(0).alias("avg_paid"),
    )

    tier_table = mo.ui.table(tier.to_pandas())

    mo.vstack(
        [scatter, mo.md("### Risk Tier Distribution"), tier_table]
    )
    return


@app.cell(hide_code=True)
def _(con, mo, pl):
    ev = pl.from_pandas(
        con.execute("""
            SELECT source_tag, count(*) as items,
                   count(*) FILTER (WHERE is_snowball) as snowball
            FROM skin_subs.evidence_base
            GROUP BY source_tag ORDER BY items DESC
        """).fetchdf()
    )

    study = pl.from_pandas(
        con.execute("""
            SELECT pub_type, count(*) as n,
                   count(*) FILTER (WHERE product_count > 0) as with_products,
                   count(*) FILTER (WHERE sample_size IS NOT NULL) as with_sample,
                   round(avg(sample_size), 0) as avg_sample
            FROM skin_subs.study_characteristics
            GROUP BY pub_type ORDER BY n DESC
        """).fetchdf()
    )

    mo.vstack(
        [
            mo.md("### Evidence Base by Source"),
            mo.ui.table(ev.to_pandas()),
            mo.md("### Study Characteristics by Type"),
            mo.ui.table(study.to_pandas()),
            mo.md(f"**Total evidence items:** {ev['items'].sum():,}"),
        ]
    )
    return


if __name__ == "__main__":
    app.run()