"""Reproduce ProteinIQ's Medicare price-per-gram analysis of monoclonal antibodies (October 6, 2026).

Input: raw/october-2026-medicare-part-b-payment-limit-file.csv, the Section 508 CSV inside CMS's
"October 2026 Medicare Part B Payment Limit Files" ZIP (final, posted September 17, 2026).

Run from this directory: python analyze.py   (Python 3.13, standard library only)
Writes medicare-price-per-gram.csv, biosimilar-comparison.csv and summary.json.
"""
import csv
import json
import re
import statistics as st

RAW = 'raw/october-2026-medicare-part-b-payment-limit-file.csv'

C, I, N, B, E, O = 'Cancer', 'Immune and inflammatory', 'Neurology', 'Blood and rare disease', 'Eye', 'Other'

# Reference (originator) product codes: HCPCS code -> (molecule, primary therapy area, format).
# Formats: mAb = conventional full-length antibody; ADC = antibody-drug conjugate;
# bispecific = bispecific antibody or bispecific T-cell engager; fragment = Fab, Fab' or scFv.
REFERENCE = {
    'J0174': ('lecanemab', N, 'mAb'), 'J0175': ('donanemab', N, 'mAb'), 'J0179': ('brolucizumab', E, 'fragment'),
    'J0202': ('alemtuzumab', N, 'mAb'), 'J0480': ('basiliximab', I, 'mAb'), 'J0490': ('belimumab', I, 'mAb'),
    'J0491': ('anifrolumab', I, 'mAb'), 'J0517': ('benralizumab', I, 'mAb'), 'J0584': ('burosumab', B, 'mAb'),
    'J0593': ('lanadelumab', I, 'mAb'), 'J0638': ('canakinumab', I, 'mAb'), 'J0717': ('certolizumab pegol', I, 'fragment'),
    'J0791': ('crizanlizumab', B, 'mAb'), 'J0897': ('denosumab', O, 'mAb'), 'J1289': ('narsoplimab', B, 'mAb'),
    'J1299': ('eculizumab', B, 'mAb'), 'J1302': ('sutimlimab', B, 'mAb'), 'J1303': ('ravulizumab', B, 'mAb'),
    'J1305': ('evinacumab', O, 'mAb'), 'J1307': ('crovalimab', B, 'mAb'), 'J1323': ('elranatamab', C, 'bispecific'),
    'J1326': ('zolbetuximab', C, 'mAb'), 'J1602': ('golimumab', I, 'mAb'), 'J1628': ('guselkumab', I, 'mAb'),
    'J1745': ('infliximab', I, 'mAb'), 'J1746': ('ibalizumab', O, 'mAb'), 'J1747': ('spesolimab', I, 'mAb'),
    'J1823': ('inebilizumab', N, 'mAb'), 'J2182': ('mepolizumab', I, 'mAb'), 'J2267': ('mirikizumab', I, 'mAb'),
    'J2323': ('natalizumab', N, 'mAb'), 'J2327': ('risankizumab', I, 'mAb'), 'J2329': ('ublituximab', N, 'mAb'),
    'J2350': ('ocrelizumab', N, 'mAb'), 'J2356': ('tezepelumab', I, 'mAb'), 'J2357': ('omalizumab', I, 'mAb'),
    'J2361': ('depemokimab', I, 'mAb'), 'J2777': ('faricimab', E, 'bispecific'), 'J2778': ('ranibizumab', E, 'fragment'),
    'J2786': ('reslizumab', I, 'mAb'), 'J2860': ('siltuximab', B, 'mAb'), 'J3032': ('eptinezumab', N, 'mAb'),
    'J3055': ('talquetamab', C, 'bispecific'), 'J3111': ('romosozumab', O, 'mAb'), 'J3241': ('teprotumumab', E, 'mAb'),
    'J3245': ('tildrakizumab', I, 'mAb'), 'J3247': ('secukinumab', I, 'mAb'), 'J3262': ('tocilizumab', I, 'mAb'),
    'J3263': ('toripalimab', C, 'mAb'), 'J3358': ('ustekinumab', I, 'mAb'), 'J3380': ('vedolizumab', I, 'mAb'),
    'J7170': ('emicizumab', B, 'bispecific'), 'J7173': ('concizumab', B, 'mAb'),
    'J9011': ('datopotamab deruxtecan', C, 'ADC'), 'J9022': ('atezolizumab', C, 'mAb'), 'J9023': ('avelumab', C, 'mAb'),
    'J9026': ('tarlatamab', C, 'bispecific'), 'J9035': ('bevacizumab', C, 'mAb'), 'J9038': ('axatilimab', I, 'mAb'),
    'J9039': ('blinatumomab', C, 'bispecific'), 'J9042': ('brentuximab vedotin', C, 'ADC'),
    'J9053': ('belantamab mafodotin', C, 'ADC'), 'J9055': ('cetuximab', C, 'mAb'), 'J9061': ('amivantamab', C, 'bispecific'),
    'J9119': ('cemiplimab', C, 'mAb'), 'J9145': ('daratumumab', C, 'mAb'), 'J9173': ('durvalumab', C, 'mAb'),
    'J9176': ('elotuzumab', C, 'mAb'), 'J9177': ('enfortumab vedotin', C, 'ADC'),
    'J9203': ('gemtuzumab ozogamicin', C, 'ADC'), 'J9204': ('mogamulizumab', C, 'mAb'), 'J9210': ('emapalumab', B, 'mAb'),
    'J9227': ('isatuximab', C, 'mAb'), 'J9228': ('ipilimumab', C, 'mAb'), 'J9229': ('inotuzumab ozogamicin', C, 'ADC'),
    'J9256': ('nipocalimab', N, 'mAb'), 'J9271': ('pembrolizumab', C, 'mAb'), 'J9272': ('dostarlimab', C, 'mAb'),
    'J9273': ('tisotumab vedotin', C, 'ADC'), 'J9275': ('cosibelimab', C, 'mAb'), 'J9276': ('zanidatamab', C, 'bispecific'),
    'J9286': ('glofitamab', C, 'bispecific'), 'J9299': ('nivolumab', C, 'mAb'), 'J9301': ('obinutuzumab', C, 'mAb'),
    'J9303': ('panitumumab', C, 'mAb'), 'J9306': ('pertuzumab', C, 'mAb'), 'J9308': ('ramucirumab', C, 'mAb'),
    'J9309': ('polatuzumab vedotin', C, 'ADC'), 'J9312': ('rituximab', C, 'mAb'),
    'J9317': ('sacituzumab govitecan', C, 'ADC'), 'J9321': ('epcoritamab', C, 'bispecific'),
    'J9326': ('telisotuzumab vedotin', C, 'ADC'), 'J9329': ('tislelizumab', C, 'mAb'),
    'J9333': ('rozanolixizumab', N, 'mAb'), 'J9345': ('retifanlimab', C, 'mAb'), 'J9347': ('tremelimumab', C, 'mAb'),
    'J9348': ('naxitamab', C, 'mAb'), 'J9349': ('tafasitamab', C, 'mAb'), 'J9350': ('mosunetuzumab', C, 'bispecific'),
    'J9353': ('margetuximab', C, 'mAb'), 'J9354': ('ado-trastuzumab emtansine', C, 'ADC'),
    'J9355': ('trastuzumab', C, 'mAb'), 'J9358': ('fam-trastuzumab deruxtecan', C, 'ADC'),
    'J9359': ('loncastuximab tesirine', C, 'ADC'), 'J9380': ('teclistamab', C, 'bispecific'),
    'J9381': ('teplizumab', O, 'mAb'), 'J9382': ('zenocutuzumab', C, 'bispecific'),
    'J9601': ('linvoseltamab', C, 'bispecific'),
}

# Biosimilar codes for each reference code. Hospital-specific tocilizumab codes (Q0234, Q0237, Q0238)
# and the COVID-19 tocilizumab code (Q0249) are excluded.
BIOSIMILARS = {
    'J9312': ['Q5115', 'Q5119', 'Q5123'],
    'J9035': ['Q5107', 'Q5118', 'Q5126', 'Q5129', 'Q5160'],
    'J9355': ['Q5112', 'Q5113', 'Q5114', 'Q5116', 'Q5117', 'Q5146'],
    'J3358': ['Q5098', 'Q5099', 'Q5100', 'Q5138', 'Q9997', 'Q9998', 'Q9999'],
    'J0897': ['Q5136', 'Q5157', 'Q5158', 'Q5161', 'Q5162', 'Q5167'],
    'J1299': ['Q5151', 'Q5152'],
    'J3262': ['Q5133', 'Q5135', 'Q5156'],
    'J1745': ['Q5103', 'Q5104', 'Q5121'],
    'J2778': ['Q5124', 'Q5128', 'Q5168'],
    'J2323': ['Q5134'],
}
# Excluded from the biosimilar summary: ranibizumab (the reference payment limit is below every
# biosimilar's) and natalizumab (one biosimilar, paid at parity). Both remain in the CSV.
SUMMARY_EXCLUDE = {'J2778', 'J2323'}


def load():
    rows = list(csv.reader(open(RAW, encoding='latin-1')))
    header = next(i for i, r in enumerate(rows) if r and r[0] == 'HCPCS Code')
    return {r[0]: r for r in rows[header + 1:] if r and r[0]}


def mg_per_unit(unit):
    m = re.fullmatch(r'\s*([\d.]+)\s*(MG|MCG)\s*', unit, re.I)
    if not m:
        raise ValueError(f'Not a mass billing unit: {unit!r}')
    value = float(m.group(1))
    return value / 1000 if m.group(2).upper() == 'MCG' else value


def per_gram(row):
    return float(row[3]) / mg_per_unit(row[2]) * 1000


def summarize(values):
    q = st.quantiles(values, n=4, method='inclusive')
    return {'n': len(values), 'median': round(st.median(values)), 'q1': round(q[0]), 'q3': round(q[2]),
            'min': round(min(values)), 'max': round(max(values))}


def main():
    data = load()
    products = []
    for code, (molecule, area, fmt) in REFERENCE.items():
        row = data[code]
        products.append({'hcpcs': code, 'short_description': row[1], 'billing_unit': row[2],
                         'payment_limit_usd': float(row[3]), 'molecule': molecule, 'therapy_area': area,
                         'format': fmt, 'usd_per_gram': round(per_gram(row), 2)})
    products.sort(key=lambda p: p['usd_per_gram'])
    with open('medicare-price-per-gram.csv', 'w', newline='', encoding='utf-8') as f:
        writer = csv.DictWriter(f, fieldnames=list(products[0]))
        writer.writeheader()
        writer.writerows(products)

    biosimilar_rows, biosimilar_summary = [], []
    for ref, codes in BIOSIMILARS.items():
        ref_value = per_gram(data[ref])
        values = []
        for code in codes:
            value = per_gram(data[code])
            values.append(value)
            biosimilar_rows.append({'molecule': REFERENCE[ref][0], 'reference_hcpcs': ref, 'biosimilar_hcpcs': code,
                                    'short_description': data[code][1], 'billing_unit': data[code][2],
                                    'payment_limit_usd': float(data[code][3]), 'usd_per_gram': round(value, 2),
                                    'reference_usd_per_gram': round(ref_value, 2)})
        if ref not in SUMMARY_EXCLUDE:
            median = st.median(values)
            biosimilar_summary.append({'molecule': REFERENCE[ref][0], 'reference_usd_per_gram': round(ref_value),
                                       'biosimilars': len(values), 'biosimilar_median_usd_per_gram': round(median),
                                       'median_discount_percent': round(100 * (1 - median / ref_value))})
    with open('biosimilar-comparison.csv', 'w', newline='', encoding='utf-8') as f:
        writer = csv.DictWriter(f, fieldnames=list(biosimilar_rows[0]))
        writer.writeheader()
        writer.writerows(biosimilar_rows)

    values = [p['usd_per_gram'] for p in products]
    summary = {
        'source': 'CMS October 2026 Medicare Part B payment limit file (2Q26 ASP data)',
        'all_antibodies': summarize(values),
        'by_format': {fmt: summarize([p['usd_per_gram'] for p in products if p['format'] == fmt])
                      for fmt in ['mAb', 'ADC', 'bispecific', 'fragment']},
        'conventional_by_area': {area: summarize([p['usd_per_gram'] for p in products
                                                  if p['format'] == 'mAb' and p['therapy_area'] == area])
                                 for area in [C, I, N, B, O]},
        'count_below_10000': sum(v < 10_000 for v in values),
        'count_above_100000': sum(v > 100_000 for v in values),
        'count_above_1000000': sum(v > 1_000_000 for v in values),
        'biosimilars': biosimilar_summary,
    }
    with open('summary.json', 'w', encoding='utf-8') as f:
        json.dump(summary, f, indent=2)
    print(json.dumps(summary, indent=2))


if __name__ == '__main__':
    main()