Zamora-PCT Score

Executable Code

"""
Zamora-PCT Score: Bayesian-Derived Clinical Score for Infection vs Flare in SLE
Based on Bayesian bivariate meta-analysis of procalcitonin diagnostic accuracy

FHE-compatible implementation for RheumaScore platform.
"""

import numpy as np
import math

# ══════════════════════════════════════════════════════════════════
# SCORE METADATA
# ══════════════════════════════════════════════════════════════════

ZAMORA_PCT_DESCRIPTORS = [
    {
        "name": "pct_elevated",
        "label": "PCT ≥ 0.5 ng/mL",
        "weight": 5,
        "category": "biomarker",
        "direction": "infection",
        "evidence": "Bayesian meta-analysis: LR+ 5.07 (95% CrI 2.82-8.37), PMID: 22605405",
    },
    {
        "name": "crp_elevated",
        "label": "CRP > 50 mg/L",
        "weight": 3,
        "category": "biomarker",
        "direction": "infection",
        "evidence": "LR+ ~2.5 for infection in SLE, PMID: 35786213",
    },
    {
        "name": "wbc_elevated",
        "label": "WBC > 11,000/μL",
        "weight": 3,
        "category": "laboratory",
        "direction": "infection",
        "evidence": "LR+ ~3.0, PMID: 33199770, 39143483",
    },
    {
        "name": "prolonged_fever",
        "label": "Fever ≥ 38.5°C for ≥3 days",
        "weight": 2,
        "category": "clinical",
        "direction": "infection",
        "evidence": "LR+ ~2.0, PMID: 35786213",
    },
    {
        "name": "low_complement",
        "label": "Low C3 or C4",
        "weight": -3,
        "category": "immunological",
        "direction": "flare",
        "evidence": "LR+ 0.4 (favors flare), PMID: 27744359",
    },
    {
        "name": "rising_anti_dsdna",
        "label": "Rising anti-dsDNA antibodies",
        "weight": -4,
        "category": "immunological",
        "direction": "flare",
        "evidence": "LR+ 0.3 (favors flare), PMID: 27744359",
    },
    {
        "name": "active_urine_sediment",
        "label": "Active urine sediment (non-infectious)",
        "weight": -2,
        "category": "laboratory",
        "direction": "flare",
        "evidence": "Favors nephritis flare, LR+ 0.5",
    },
]

ZAMORA_PCT_WEIGHTS = np.array(
    [d["weight"] for d in ZAMORA_PCT_DESCRIPTORS], dtype=np.int8
)
N_ZAMORA_PCT = len(ZAMORA_PCT_DESCRIPTORS)

SCORE_META = {
    "zamora_pct": {
        "name": "Zamora-PCT Score",
        "full_name": "Zamora-PCT Bayesian Score for Infection vs Flare in SLE",
        "disease": "SLE",
        "category": "differential_diagnosis",
        "descriptors": ZAMORA_PCT_DESCRIPTORS,
        "n_inputs": N_ZAMORA_PCT,
        "input_type": "binary",
        "input_format": f"Array of {N_ZAMORA_PCT} binary integers (0=absent, 1=present)",
        "score_range": [-9, 13],
        "version": "1.0",
        "methodology": "Bayesian bivariate meta-analysis (Reitsma/Zamora model), MCMC 50k iterations",
        "bayesian_posterior": {
            "pooled_sensitivity": {"median": 0.742, "ci95": [0.588, 0.860]},
            "pooled_specificity": {"median": 0.854, "ci95": [0.749, 0.911]},
            "positive_lr": {"median": 5.07, "ci95": [2.82, 8.37]},
            "negative_lr": {"median": 0.30, "ci95": [0.16, 0.49]},
            "diagnostic_or": {"median": 16.88, "ci95": [6.84, 41.15]},
            "n_studies": 10,
            "n_patients": 604,
        },
        "references": [
            "Zamora J et al. BMC Med Res Methodol 2006;6:31. PMID:16836745",
            "Wu JY et al. Arthritis Rheum 2012;64:3034-42. PMID:22605405",
        ],
    }
}


# ══════════════════════════════════════════════════════════════════
# COMPUTE FUNCTION (FHE-compatible: integer arithmetic only)
# ══════════════════════════════════════════════════════════════════


def compute_zamora_pct(indicators: list[int]) -> int:
    """
    Compute Zamora-PCT score from binary indicators.

    Args:
        indicators: List of 7 binary integers (0 or 1):
            [0] pct_elevated       - PCT ≥ 0.5 ng/mL
            [1] crp_elevated       - CRP > 50 mg/L
            [2] wbc_elevated       - WBC > 11,000/μL
            [3] prolonged_fever    - Fever ≥38.5°C ≥3 days
            [4] low_complement     - Low C3 or C4
            [5] rising_anti_dsdna  - Rising anti-dsDNA
            [6] active_urine_sediment - Active urine sediment

    Returns:
        Integer score from -9 to +13.
        Positive = infection likely; Negative = flare likely.
    """
    assert len(indicators) == N_ZAMORA_PCT, f"Expected {N_ZAMORA_PCT} inputs"
    score = 0
    for i in range(N_ZAMORA_PCT):
        score += indicators[i] * int(ZAMORA_PCT_WEIGHTS[i])
    return score


# ══════════════════════════════════════════════════════════════════
# INTERPRET FUNCTION
# ══════════════════════════════════════════════════════════════════

INTERPRETATION_THRESHOLDS = [
    (-9, -3, "low", "Flare Likely", "infection_probability_lt_20"),
    (-2, 2, "moderate", "Indeterminate", "infection_probability_20_40"),
    (3, 6, "high", "Moderate Infection Risk", "infection_probability_40_80"),
    (7, 13, "very_high", "High Infection Risk", "infection_probability_gt_80"),
]


def interpret_zamora_pct(score: int, pretest_prevalence: float = 0.40) -> dict:
    """
    Interpret Zamora-PCT score with Bayesian probability estimation.

    Args:
        score: Integer score from compute_zamora_pct()
        pretest_prevalence: Prior probability of infection (default 0.40)

    Returns:
        Dict with interpretation, risk category, and post-test probability.
    """
    # Bayesian update: post_odds = pre_odds * exp(score/3)
    pre_odds = pretest_prevalence / (1 - pretest_prevalence)
    post_odds = pre_odds * math.exp(score / 3.0)
    post_prob = post_odds / (1 + post_odds)

    # Determine category
    risk_level = "indeterminate"
    risk_label = "Indeterminate"
    for lo, hi, level, label, _ in INTERPRETATION_THRESHOLDS:
        if lo <= score <= hi:
            risk_level = level
            risk_label = label
            break

    # Clinical recommendations
    if score >= 7:
        recommendation = (
            "High probability of bacterial infection. "
            "Obtain blood/urine cultures, start empiric antibiotics. "
            "Consider CT chest/abdomen if source unclear."
        )
    elif score >= 3:
        recommendation = (
            "Moderate probability of infection. "
            "Obtain cultures, consider empiric antibiotics while awaiting results. "
            "Monitor closely for clinical deterioration."
        )
    elif score >= -2:
        recommendation = (
            "Indeterminate. Both infection and flare possible. "
            "Full workup including cultures AND immunological markers. "
            "Consider treating both until diagnosis clarified."
        )
    else:
        recommendation = (
            "Autoimmune flare most likely. "
            "Adjust immunosuppression as indicated. "
            "Continue monitoring for occult infection. "
            "Repeat PCT in 24-48h if clinical concern persists."
        )

    return {
        "score": score,
        "score_range": "−9 to +13",
        "risk_level": risk_level,
        "risk_label": risk_label,
        "post_test_probability_infection": round(post_prob, 3),
        "pretest_prevalence_used": pretest_prevalence,
        "recommendation": recommendation,
        "methodology": "Bayesian bivariate meta-analysis (Zamora/Reitsma), k=10 studies, n=604",
        "pct_meta_analysis": {
            "pooled_sensitivity": 0.742,
            "pooled_specificity": 0.854,
            "positive_lr": 5.07,
            "negative_lr": 0.30,
            "auc_hsroc": 0.87,
        },
    }


# ══════════════════════════════════════════════════════════════════
# FHE CIRCUIT DEFINITION (for Zama Concrete)
# ══════════════════════════════════════════════════════════════════


def zamora_pct_fhe_circuit(indicators):
    """
    FHE-compatible circuit for Zamora-PCT score.
    All operations are integer additions compatible with Concrete FHE.

    Input: 7 encrypted binary values (0 or 1)
    Output: encrypted integer score (shifted by +9 to keep non-negative: range 0-22)
    """
    # Weights: [+5, +3, +3, +2, -3, -4, -2]
    # Shift by +9 so output is always non-negative (FHE-friendly)
    score = 9  # offset
    score += indicators[0] * 5   # PCT
    score += indicators[1] * 3   # CRP
    score += indicators[2] * 3   # WBC
    score += indicators[3] * 2   # Fever
    # For negative weights, subtract via complement:
    # -3 * x = 3 * (1-x) - 3, but simpler: just subtract
    score -= indicators[4] * 3   # Low complement (negative = flare)
    score -= indicators[5] * 4   # Rising anti-dsDNA (negative = flare)
    score -= indicators[6] * 2   # Active sediment (negative = flare)
    return score  # range: 0 (=-9 real) to 22 (=+13 real)


# ══════════════════════════════════════════════════════════════════
# SELF-TEST
# ══════════════════════════════════════════════════════════════════

if __name__ == "__main__":
    print("=" * 60)
    print("Zamora-PCT Score — Self-Test")
    print("=" * 60)

    # Test 1: Clear infection pattern
    case1 = [1, 1, 1, 1, 0, 0, 0]  # PCT+, CRP+, WBC+, fever+
    s1 = compute_zamora_pct(case1)
    r1 = interpret_zamora_pct(s1)
    print(f"\nCase 1 (infection): score={s1}, prob={r1['post_test_probability_infection']:.1%}")
    print(f"  → {r1['risk_label']}: {r1['recommendation'][:80]}...")

    # Test 2: Clear flare pattern
    case2 = [0, 0, 0, 0, 1, 1, 1]  # low C3, rising dsDNA, active sediment
    s2 = compute_zamora_pct(case2)
    r2 = interpret_zamora_pct(s2)
    print(f"\nCase 2 (flare): score={s2}, prob={r2['post_test_probability_infection']:.1%}")
    print(f"  → {r2['risk_label']}: {r2['recommendation'][:80]}...")

    # Test 3: Mixed/indeterminate
    case3 = [1, 0, 0, 0, 1, 0, 0]  # PCT+ but low complement
    s3 = compute_zamora_pct(case3)
    r3 = interpret_zamora_pct(s3)
    print(f"\nCase 3 (mixed): score={s3}, prob={r3['post_test_probability_infection']:.1%}")
    print(f"  → {r3['risk_label']}: {r3['recommendation'][:80]}...")

    # Test 4: FHE circuit
    fhe_result = zamora_pct_fhe_circuit(case1)
    real_score = fhe_result - 9
    print(f"\nFHE circuit test: raw={fhe_result}, real_score={real_score} (expected {s1})")
    assert real_score == s1, "FHE circuit mismatch!"

    print("\n✅ All tests passed.")

Demo Output

Case 1 (infection): score=13, prob=98.1%
Case 2 (flare): score=-9, prob=3.2%
Case 3 (mixed): score=2, prob=56.5%
FHE circuit test: passed