---
name: hcq-retina
description: Hydroxychloroquine retinal toxicity risk stratification with Monte Carlo uncertainty estimation based on AAO 2016/2020 screening guidelines.
authors: Erick Adrián Zamora Tehozol, DNAI, RheumaAI
version: 1.0.0
tags: [hydroxychloroquine, retinal-toxicity, screening, OCT, AAO-guidelines, rheumatology, SLE, monte-carlo, DeSci, RheumaAI]
---

# HCQ-RETINA

**Hydroxychloroquine Retinal Toxicity Risk Stratification with Monte Carlo Uncertainty Estimation**

## Purpose

Hydroxychloroquine (HCQ) is a cornerstone therapy in SLE and RA but carries cumulative risk of irreversible retinal toxicity (bull's eye maculopathy). The AAO 2016 revised guidelines emphasize screening based on real body weight dosing and cumulative exposure. HCQ-RETINA quantifies this risk as a composite score (0–100) across 8 domains and provides evidence-based screening recommendations.

## Clinical Problem

- HCQ retinal toxicity prevalence: ~7.5% after 5 years of use (Melles 2020)
- Risk is dose- and duration-dependent, markedly increasing after 5 years
- Toxicity is IRREVERSIBLE once established — early detection is critical
- AAO recommends ≤5.0 mg/kg/day based on real body weight (not ideal)
- CKD, tamoxifen, and pre-existing macular disease amplify risk

## Domains (8 weighted)

| Domain | Weight | Key Reference |
|--------|--------|---------------|
| Cumulative dose (g) | 0.25 | Melles & Jorge 2020 |
| Duration (years) | 0.20 | AAO 2016, Melles 2020 |
| Daily dose/kg | 0.20 | AAO 2016 (≤5.0 mg/kg/day) |
| Renal function (eGFR) | 0.10 | Yusuf 2017, pharmacokinetics |
| Tamoxifen use | 0.08 | Marmor 2016 |
| Macular disease | 0.07 | AAO 2016 |
| Age | 0.05 | Browning 2014 |
| CYP2D6 PM status | 0.05 | Petri 2020 |

## Risk Categories

- **LOW** (<15): Baseline exam, annual after 5 years
- **MODERATE** (15-34): Annual SD-OCT + 10-2 VF
- **HIGH** (35-59): Screen every 6 months, consider dose reduction
- **VERY HIGH** (≥60): Urgent ophthalmology referral, consider discontinuation

## Usage

```bash
python3 hcq_retina.py
```

## References

1. Marmor MF et al. Revised recommendations on screening for chloroquine and hydroxychloroquine retinopathy. Ophthalmology 2016;123:1386-94.
2. Melles RB, Jorge AM. Prevalence of hydroxychloroquine retinopathy. JAMA Ophthalmol 2020;138(4):e200370.
3. Petri M et al. Hydroxychloroquine blood levels predict toxicity. Arthritis Rheumatol 2020;72(11):1894-1901.
4. Browning DJ. Hydroxychloroquine and chloroquine retinopathy: screening for drug toxicity. Am J Ophthalmol 2014;158(6):1207-12.
5. Yusuf IH et al. Hydroxychloroquine retinopathy. Eye 2017;31:828-845.


## Executable Code

```python
#!/usr/bin/env python3
"""
HCQ-RETINA: Hydroxychloroquine Retinal Toxicity Risk Stratification
with Monte Carlo Uncertainty Estimation

Authors: Erick Adrián Zamora Tehozol, DNAI, RheumaAI
License: MIT

Computes a composite retinal toxicity risk score (0-100) for patients on
hydroxychloroquine (HCQ), based on AAO 2016/2020 screening guidelines,
Marmor et al. 2016, Melles & Jorge 2020 prevalence data, and
pharmacokinetic principles.

Domains:
  1. Cumulative dose (g) — strongest predictor per Melles 2020
  2. Duration of use (years) — exponential risk after 5 years
  3. Daily dose per real body weight (mg/kg/day) — AAO threshold ≤5.0
  4. Renal function (eGFR) — HCQ renally cleared, CKD increases tissue levels
  5. Tamoxifen co-administration — doubles risk (Marmor 2016)
  6. Macular disease — pre-existing retinal vulnerability
  7. Age — older patients have higher risk
  8. CYP2D6 poor metabolizer status — reduced clearance

References:
  - AAO Screening Guidelines 2016 (Marmor MF et al. Ophthalmology 2016;123:1386-94)
  - Melles RB, Jorge AM. JAMA Ophthalmol 2020;138(4):e200370
  - Marmor MF. Am J Ophthalmol 2017;177:xv-xvi
  - Petri M et al. Arthritis Rheumatol 2020;72(11):1894-1901
  - Browning DJ. Am J Ophthalmol 2014;158(6):1207-12
  - Yusuf IH et al. Eye 2017;31:828-845
"""

import json
import math
import random
import sys
from dataclasses import dataclass, field
from typing import Optional


@dataclass
class HCQPatient:
    """Patient data for HCQ retinal toxicity risk assessment."""
    cumulative_dose_g: float          # Total lifetime HCQ in grams
    duration_years: float             # Years on HCQ
    daily_dose_mg: float              # Current daily dose (mg)
    real_body_weight_kg: float        # Actual body weight (kg)
    egfr: float = 90.0               # eGFR mL/min/1.73m² (CKD-EPI)
    tamoxifen: bool = False           # Concurrent tamoxifen use
    macular_disease: bool = False     # Pre-existing macular pathology
    age: int = 50                     # Patient age
    cyp2d6_poor_metabolizer: bool = False  # CYP2D6 PM status


@dataclass
class DomainScore:
    """Individual domain score with details."""
    name: str
    raw_value: float
    score: float        # 0-100 within domain
    weight: float
    weighted: float     # score * weight
    detail: str


@dataclass
class HCQResult:
    """Complete risk assessment result."""
    composite_score: float
    risk_category: str
    domains: list
    ci_lower: float = 0.0
    ci_upper: float = 0.0
    screening_recommendation: str = ""
    monitoring_frequency: str = ""
    notes: list = field(default_factory=list)


# Domain weights (sum = 1.0)
WEIGHTS = {
    "cumulative_dose":   0.25,
    "duration":          0.20,
    "daily_dose_kg":     0.20,
    "renal_function":    0.10,
    "tamoxifen":         0.08,
    "macular_disease":   0.07,
    "age":               0.05,
    "cyp2d6":            0.05,
}


def score_cumulative_dose(dose_g: float) -> tuple:
    """
    Cumulative dose is the strongest predictor.
    Melles 2020: <100g ~2%, 100-200g ~2%, 200-400g ~5%, >600g ~20%, >1000g ~40%+
    """
    if dose_g < 100:
        s = (dose_g / 100) * 10
    elif dose_g < 200:
        s = 10 + ((dose_g - 100) / 100) * 10
    elif dose_g < 400:
        s = 20 + ((dose_g - 200) / 200) * 20
    elif dose_g < 600:
        s = 40 + ((dose_g - 400) / 200) * 20
    elif dose_g < 1000:
        s = 60 + ((dose_g - 600) / 400) * 25
    else:
        s = min(85 + ((dose_g - 1000) / 500) * 15, 100)
    detail = f"Cumulative HCQ: {dose_g:.0f}g"
    return min(s, 100), detail


def score_duration(years: float) -> tuple:
    """
    AAO 2016: risk increases markedly after 5 years.
    Melles 2020: ~2% at 5y, ~20% at 20y.
    Sigmoid growth centered at 10 years.
    """
    if years <= 1:
        s = years * 3
    elif years <= 5:
        s = 3 + (years - 1) * 3.5
    elif years <= 10:
        s = 17 + (years - 5) * 8
    elif years <= 20:
        s = 57 + (years - 10) * 3.5
    else:
        s = min(92 + (years - 20) * 1.0, 100)
    detail = f"Duration: {years:.1f} years"
    return min(s, 100), detail


def score_daily_dose_kg(daily_mg: float, weight_kg: float) -> tuple:
    """
    AAO 2016 threshold: ≤5.0 mg/kg/day real body weight.
    Risk escalates above 5.0, critical above 6.5.
    """
    dose_per_kg = daily_mg / max(weight_kg, 30)
    if dose_per_kg <= 4.0:
        s = (dose_per_kg / 4.0) * 15
    elif dose_per_kg <= 5.0:
        s = 15 + ((dose_per_kg - 4.0) / 1.0) * 25
    elif dose_per_kg <= 6.5:
        s = 40 + ((dose_per_kg - 5.0) / 1.5) * 35
    else:
        s = min(75 + ((dose_per_kg - 6.5) / 2.0) * 25, 100)
    detail = f"Daily dose: {daily_mg:.0f}mg / {weight_kg:.0f}kg = {dose_per_kg:.2f} mg/kg/day (AAO threshold: ≤5.0)"
    return min(s, 100), detail


def score_renal(egfr: float) -> tuple:
    """
    HCQ is ~40-50% renally excreted. CKD increases tissue accumulation.
    eGFR <30 is high risk, 30-60 moderate, >60 low.
    """
    if egfr >= 90:
        s = 0
    elif egfr >= 60:
        s = ((90 - egfr) / 30) * 20
    elif egfr >= 30:
        s = 20 + ((60 - egfr) / 30) * 40
    elif egfr >= 15:
        s = 60 + ((30 - egfr) / 15) * 30
    else:
        s = 90 + ((15 - egfr) / 15) * 10
    detail = f"eGFR: {egfr:.0f} mL/min/1.73m²"
    return min(s, 100), detail


def score_tamoxifen(on_tamoxifen: bool) -> tuple:
    """Tamoxifen co-use doubles retinal toxicity risk (Marmor 2016)."""
    s = 80 if on_tamoxifen else 0
    detail = f"Tamoxifen: {'Yes — doubles retinal risk' if on_tamoxifen else 'No'}"
    return s, detail


def score_macular(has_macular: bool) -> tuple:
    """Pre-existing macular disease increases vulnerability."""
    s = 70 if has_macular else 0
    detail = f"Macular disease: {'Yes — increased vulnerability' if has_macular else 'No'}"
    return s, detail


def score_age(age: int) -> tuple:
    """Older age associated with higher risk (pharmacokinetic changes, retinal aging)."""
    if age < 40:
        s = 0
    elif age < 60:
        s = ((age - 40) / 20) * 30
    elif age < 75:
        s = 30 + ((age - 60) / 15) * 40
    else:
        s = min(70 + ((age - 75) / 15) * 30, 100)
    detail = f"Age: {age} years"
    return min(s, 100), detail


def score_cyp2d6(is_pm: bool) -> tuple:
    """CYP2D6 poor metabolizers have reduced HCQ clearance."""
    s = 75 if is_pm else 0
    detail = f"CYP2D6 PM: {'Yes — reduced clearance' if is_pm else 'No/Unknown'}"
    return s, detail


def compute_risk(patient: HCQPatient, n_simulations: int = 5000, seed: int = 42) -> HCQResult:
    """Compute composite retinal toxicity risk score with Monte Carlo CI."""
    
    # Score each domain
    domains = []
    
    scorers = [
        ("cumulative_dose", lambda: score_cumulative_dose(patient.cumulative_dose_g)),
        ("duration", lambda: score_duration(patient.duration_years)),
        ("daily_dose_kg", lambda: score_daily_dose_kg(patient.daily_dose_mg, patient.real_body_weight_kg)),
        ("renal_function", lambda: score_renal(patient.egfr)),
        ("tamoxifen", lambda: score_tamoxifen(patient.tamoxifen)),
        ("macular_disease", lambda: score_macular(patient.macular_disease)),
        ("age", lambda: score_age(patient.age)),
        ("cyp2d6", lambda: score_cyp2d6(patient.cyp2d6_poor_metabolizer)),
    ]
    
    composite = 0.0
    for name, scorer in scorers:
        raw_score, detail = scorer()
        w = WEIGHTS[name]
        weighted = raw_score * w
        composite += weighted
        domains.append(DomainScore(
            name=name,
            raw_value=raw_score,
            score=round(raw_score, 1),
            weight=w,
            weighted=round(weighted, 1),
            detail=detail,
        ))
    
    composite = round(min(composite, 100), 1)
    
    # Monte Carlo for uncertainty estimation
    rng = random.Random(seed)
    mc_scores = []
    for _ in range(n_simulations):
        mc_composite = 0.0
        noise_cum = patient.cumulative_dose_g * (1 + rng.gauss(0, 0.05))
        noise_dur = patient.duration_years * (1 + rng.gauss(0, 0.03))
        noise_dose = patient.daily_dose_mg * (1 + rng.gauss(0, 0.02))
        noise_egfr = patient.egfr * (1 + rng.gauss(0, 0.08))
        noise_age = patient.age + rng.gauss(0, 0.5)
        
        s1, _ = score_cumulative_dose(max(0, noise_cum))
        s2, _ = score_duration(max(0, noise_dur))
        s3, _ = score_daily_dose_kg(max(0, noise_dose), patient.real_body_weight_kg)
        s4, _ = score_renal(max(5, noise_egfr))
        s5, _ = score_tamoxifen(patient.tamoxifen)
        s6, _ = score_macular(patient.macular_disease)
        s7, _ = score_age(max(18, int(noise_age)))
        s8, _ = score_cyp2d6(patient.cyp2d6_poor_metabolizer)
        
        vals = [s1, s2, s3, s4, s5, s6, s7, s8]
        keys = list(WEIGHTS.keys())
        for i, k in enumerate(keys):
            mc_composite += vals[i] * WEIGHTS[k]
        mc_scores.append(min(mc_composite, 100))
    
    mc_scores.sort()
    ci_lower = round(mc_scores[int(0.025 * n_simulations)], 1)
    ci_upper = round(mc_scores[int(0.975 * n_simulations)], 1)
    
    # Risk category
    if composite < 15:
        cat = "LOW"
        screening = "Baseline exam at start or within first year. Annual OCT/VF after 5 years."
        freq = "Annual after 5 years of use"
    elif composite < 35:
        cat = "MODERATE"
        screening = "Annual screening with SD-OCT and 10-2 VF starting now."
        freq = "Annual"
    elif composite < 60:
        cat = "HIGH"
        screening = "Screen every 6 months with SD-OCT, mfERG, and 10-2 VF. Consider dose reduction."
        freq = "Every 6 months"
    else:
        cat = "VERY HIGH"
        screening = "URGENT: Ophthalmology referral. Consider discontinuation. SD-OCT + mfERG + FAF immediately."
        freq = "Every 3-6 months or discontinue"
    
    # Clinical notes
    notes = []
    dose_per_kg = patient.daily_dose_mg / max(patient.real_body_weight_kg, 30)
    if dose_per_kg > 5.0:
        notes.append(f"⚠️ Daily dose {dose_per_kg:.1f} mg/kg exceeds AAO maximum of 5.0 mg/kg/day — reduce dose")
    if patient.duration_years > 5 and composite < 35:
        notes.append("Duration >5 years — annual screening recommended per AAO 2016")
    if patient.tamoxifen:
        notes.append("⚠️ Tamoxifen co-use doubles retinal toxicity risk (Marmor 2016)")
    if patient.egfr < 60:
        notes.append("⚠️ CKD (eGFR<60) — HCQ accumulation risk, consider dose adjustment")
    if patient.cumulative_dose_g > 1000:
        notes.append("⚠️ Cumulative dose >1000g — very high prevalence of toxicity (Melles 2020)")
    
    return HCQResult(
        composite_score=composite,
        risk_category=cat,
        domains=[{
            "name": d.name,
            "score": d.score,
            "weight": d.weight,
            "weighted": d.weighted,
            "detail": d.detail,
        } for d in domains],
        ci_lower=ci_lower,
        ci_upper=ci_upper,
        screening_recommendation=screening,
        monitoring_frequency=freq,
        notes=notes,
    )


def print_result(result: HCQResult, label: str = ""):
    """Pretty-print a risk assessment result."""
    if label:
        print(f"\n{'='*60}")
        print(f"  {label}")
        print(f"{'='*60}")
    
    print(f"\n  COMPOSITE SCORE: {result.composite_score}/100  [{result.risk_category}]")
    print(f"  95% CI: [{result.ci_lower}, {result.ci_upper}]")
    print(f"\n  Screening: {result.screening_recommendation}")
    print(f"  Frequency: {result.monitoring_frequency}")
    
    print(f"\n  {'Domain':<20} {'Score':>6} {'Weight':>7} {'Weighted':>9}")
    print(f"  {'-'*42}")
    for d in result.domains:
        print(f"  {d['name']:<20} {d['score']:>6.1f} {d['weight']:>7.2f} {d['weighted']:>9.1f}")
    
    if result.notes:
        print(f"\n  Clinical Notes:")
        for n in result.notes:
            print(f"    {n}")
    print()


def demo():
    """Run three clinical scenarios."""
    
    # Scenario 1: Low risk — young SLE patient, early in treatment
    p1 = HCQPatient(
        cumulative_dose_g=73,      # ~1 year at 200mg/day
        duration_years=1.0,
        daily_dose_mg=200,
        real_body_weight_kg=65,
        egfr=105,
        age=32,
    )
    r1 = compute_risk(p1)
    print_result(r1, "Scenario 1: Young SLE, 1 year HCQ 200mg/day, 65kg")
    
    # Scenario 2: High risk — long-term RA patient, high cumulative dose
    p2 = HCQPatient(
        cumulative_dose_g=730,     # ~10 years at 200mg/day
        duration_years=10,
        daily_dose_mg=400,
        real_body_weight_kg=55,    # Small patient → high mg/kg
        egfr=52,                   # CKD stage 3a
        tamoxifen=False,
        macular_disease=False,
        age=68,
    )
    r2 = compute_risk(p2)
    print_result(r2, "Scenario 2: RA 10y, HCQ 400mg/day, 55kg, CKD3a, age 68")
    
    # Scenario 3: Very high risk — long duration, tamoxifen, CKD
    p3 = HCQPatient(
        cumulative_dose_g=1200,    # ~16 years at 200mg/day
        duration_years=16,
        daily_dose_mg=400,
        real_body_weight_kg=50,
        egfr=28,                   # CKD stage 4
        tamoxifen=True,
        macular_disease=True,
        age=72,
        cyp2d6_poor_metabolizer=True,
    )
    r3 = compute_risk(p3)
    print_result(r3, "Scenario 3: SLE 16y, HCQ 400mg, tamoxifen, CKD4, macular dz, CYP2D6 PM")
    
    return r1, r2, r3


if __name__ == "__main__":
    demo()

```