Skip to content

Repository files navigation

SeizeVar

Mechanism-aware variant interpretation for monogenic epilepsy

Not just "is this variant pathogenic?" — but "which direction does it break the channel?"

English 简体中文

Accepted License: MIT Code DOI Data DOI Python 3.10+ PyTorch


TL;DR

SeizeVar is a 49-gene epilepsy variant-interpretation framework that answers a question mainstream predictors do not: direction of effect. It pairs a consensus pathogenicity layer (Random Forest + ESM-2 LoRA cross-attention) with a gain-versus-loss-of-function mechanism classifier and a deterministic sodium-channel direction rule, all scored under a leave-one-out dynamic-evidence design that re-reads the live ClinVar reference at inference time.

Scale 45,143 variants end-to-end — 4,576 train, 11,274 validation, 29,293 prospectively scored VUS
Key result Mechanism direction at AUROC 0.736 under leave-one-gene-out, where AlphaMissense scores 0.428, REVEL 0.502, MetaRNN 0.401
Precision Two-head consensus reaches 78 % PPV, roughly double the 29–41 % of single competitor tools
Output 4,708 consensus Likely-Pathogenic VUS calls, 1,500 carrying a definite mechanism direction
Status Accepted at Clinical Genetics (Wiley), MS 6385742

The problem

Missense pathogenicity prediction is close to saturated: several tools now exceed AUROC 0.97 on expert-curated ClinVar. Yet 29,293 canonical-missense ClinVar entries across this 49-gene epilepsy panel still sit at "uncertain significance." Raw discrimination is no longer the bottleneck.

What clinicians actually need next is direction. In monogenic epilepsy this is not academic: SCN1A loss-of-function variants are commonly worsened by sodium-channel blockers, whereas SCN8A and neonatal SCN2A gain-of-function variants are recognised indications for high-dose blockade. Same gene family, opposite molecular consequence, opposite clinical implication.

General-purpose predictors are structurally unable to help here. They return a single pathogenicity probability, and when that probability is repurposed to separate gain-of-function from loss-of-function it performs at or below chance — AlphaMissense 0.428, MetaRNN 0.401. SeizeVar treats direction as its own supervised task, and adds a second missing capability: evidence that updates as ClinVar grows, without retraining.

Why SeizeVar?

Capability Typical variant predictor SeizeVar
Pathogenicity score Yes Yes — Random Forest + ESM-2 LoRA consensus
Mechanism direction (GoF vs LoF) No Yes — a dedicated supervised classifier
Evidence that updates with ClinVar No — frozen at training time Yes — five features recomputed at inference
Epilepsy-gene specific Rarely Yes — 49-gene ClinGen-curated panel
Precision on flagged candidates 29–41 % PPV 78 % PPV via two-head consensus gating
Released prospective predictions Rarely 29,293 scored VUS, 1,500 with a direction call

The two-layer design is deliberate: pathogenicity and mechanism are separate biological questions supported by different evidence types, and a single score cannot simultaneously answer "is this variant harmful?" and "which direction does it push channel function?"

Scope note. Earlier versions of this project described the output as a prescribing recommendation. Following peer review, the framing is now explicitly a molecular-mechanism prediction. Direction of effect is one input a clinician may weigh; it is not a treatment recommendation, and this repository does not issue one.


How it works

flowchart TB
    subgraph S1["① Data — 01_data/"]
        A["ClinVar April 2026<br/>GRCh38, 49-gene panel"] --> B["Strict HGVSp missense filter<br/>42,862 pass, 53.8% discarded"]
        B --> C["Split by evidence tier"]
        C --> D1["Train 4,576<br/>6 deduplicated sources"]
        C --> D2["6 validation cohorts<br/>11,274 variants"]
        C --> D3["VUS pool<br/>29,293"]
        D1 & D2 & D3 --> E["N×N leakage audit<br/>zero variant overlap · KCNQ2 held out"]
    end

    subgraph S2["② Features — 02_features/"]
        E --> F["39 features in 4 tiers"]
        F --> F1["gene-level ×8"]
        F --> F2["residue-level ×20<br/>AlphaFold v4 structure"]
        F --> F3["substitution-level ×6"]
        F --> F4["dynamic-evidence ×5<br/>leave-one-out enforced"]
        F1 & F2 & F3 & F4 --> FM["39-dim matrix<br/>z-scored on train statistics"]
    end

    subgraph S3["③ Models — 03_models/"]
        FM --> G1["Random Forest<br/>gene-grouped 5-fold CV"]
        FM --> G2["ESM-2 t30 + LoRA<br/>cross-attention fusion"]
        FM --> G3["Mechanism head v5<br/>XGBoost, GoF vs LoF"]
    end

    subgraph S4["④ Decision layer — 05_vus_application/"]
        G1 & G2 --> H["Consensus gate<br/>both heads must agree"]
        G3 --> I["Direction call<br/>P·GoF thresholds"]
        H --> J["LP / LB / Uncertain"]
        I --> K["Na⁺-channel direction rule<br/>6 genes, deterministic"]
        J & K --> L["29,293 scored VUS<br/>4,708 Likely-Pathogenic"]
    end

    subgraph S5["⑤ Benchmark — 06_competitors/"]
        L --> M["9 tools × 6 cohorts<br/>AlphaMissense · REVEL · MetaRNN<br/>CADD · VEST4 · SIFT · PolyPhen-2 · LoGoFunc"]
    end

    style S1 fill:#e8f4f8,stroke:#4a90a4
    style S2 fill:#f0f4e8,stroke:#7a9a4a
    style S3 fill:#f8e8f0,stroke:#a44a7a
    style S4 fill:#fff4e0,stroke:#c49a3a
    style S5 fill:#eeeeee,stroke:#888888
Loading

The pathogenicity head: language model meets tabular features

The ESM-2 head does something more interesting than concatenating an embedding to a feature vector. The mutant-residue embedding is the query, and each of the 39 hand-engineered features becomes its own key/value token with a learned per-feature identity embedding. The model therefore learns which biological evidence to attend to for this particular variant, and the attention weights are directly readable.

flowchart LR
    A["Protein sequence<br/>+ mutant position"] --> B["ESM-2 t30<br/>150M params · FROZEN"]
    B -.->|"LoRA r=8, α=16<br/>last 2 layers, q/k/v/out"| B
    B --> C["Mutant-residue<br/>embedding 640-d"]
    C --> D["Linear → 128-d<br/>QUERY"]

    E["39 biological<br/>features"] --> F["Per-feature Linear 1→128<br/>+ learned feature-ID embedding"]
    F --> G["39 tokens × 128-d<br/>KEYS / VALUES"]

    D --> H["Multi-head cross-attention<br/>4 heads · residual · LayerNorm"]
    G --> H
    H --> I["MLP 128→64→1"]
    I --> J["Pathogenicity logit"]
    H --> K["Attention weights [1×39]<br/>→ interpretability"]

    style B fill:#e8f0f8,stroke:#4a7aa4
    style H fill:#f8e8e8,stroke:#a44a4a
    style K fill:#f0f8e8,stroke:#6a9a4a
Loading

Only ~132 K parameters are trainable — 0.09 % of the backbone (~82 K of those are the LoRA adapters themselves). The interpretability payoff is concrete: the model places 27.6 % of its attention on channel-topology features for ion-channel calls versus 20.4 % for non-channel calls, i.e. it learned the relevant biology rather than memorising gene identity.

Dynamic evidence: predictions that improve without retraining

Five of the 39 features are computed against the live ClinVar reference at inference time — for example the spatial distance to the nearest pathogenic variant (SPDV) in the AlphaFold structure. During training these are strictly leave-one-out (the query variant, and for allelic series its own substitution, are removed from the reference), which is verified empirically: 0 of 3,153 training P/LP variants have SPDV = 0.

The consequence is that a new ClinVar submission raises a neighbouring VUS's score with no retraining. Measured across six historical ClinVar snapshots (2015–2026) with architecture held constant, median SPDV fell from 26 Å to 7 Å and Random-Forest recall@0.7 rose from 19 % to 34 %, reaching near-parity with frozen AlphaMissense (36 %) by 2024. These five features are 13 % of the feature count but carry 35.9 % of Gini importance on ion-channel genes.


Results

Mechanism direction — the core contribution

Evaluated on 415 variants across 22 panel genes with patch-clamp / TEVC ground truth, under honest leave-one-gene-out (every test variant scored by a model that never saw any variant from its gene):

Method AUROC Coverage
SeizeVar mechanism head v5 0.736 [0.685, 0.785] 100 %
ESM-2 zero-shot LLR 0.693 100 %
MetaRNN 0.623 27 %
AlphaMissense 0.601 89 %
CADD / SIFT / REVEL / PolyPhen-2 / VEST4 0.55 – 0.61 ~25 %

Head-to-head against LoGoFunc, the one dedicated GoF/LoF tool, on the 329 coordinate-available variants: SeizeVar 0.770 [0.719, 0.819] vs LoGoFunc 0.760 [0.707, 0.817] — statistically indistinguishable, with SeizeVar covering all 22 panel genes.

On honest numbers. An earlier version of this work reported mechanism AUROC 0.905. A self-audit found that 86 % of the evaluation set had been seen in training. The contamination-free leave-one-gene-out estimate is 0.736, and that is the number this repository and the accepted paper report. Similarly, LoGoFunc's apparent 0.821 was traced to an 85-variant training-overlap artifact and corrected to 0.760. Superseded figures are documented in the revision history rather than quietly dropped.

Three independent confirmations of the mechanism head, none of them ClinVar-derived:

  • Retrospective drug response (n = 110): concordance between predicted mechanism and sodium-channel-blocker response was 43 % (37/86) when aligned versus 4 % (1/24) when not — OR 17.37 [2.24, 134.52], p = 1.5 × 10⁻⁴.
  • KCNQ2 deep mutational scanning (n = 301, gene fully held out of training): Spearman ρ = +0.378, p = 1.1 × 10⁻¹¹.
  • TSC2 deep mutational scanning (n = 8,166): 100 % LoF-leaning, max P(GoF) = 0.205 — correct for a tumour-suppressor.

Per-gene direction emerges without gene identity as an input feature: LoF-leaning for SCN1A/SCN1B, GoF-leaning for SCN2A/SCN3A/SCN8A.

Pathogenicity — six cohorts, nine tools

Val-ClinGen Val-Temporal Reclassified-VUS KCNQ2-DMS TSC2-MAVE SCN2A-Ephys
SeizeVar Hybrid (deployed) 0.989 0.990 0.965 0.898 0.855 1.000
SeizeVar-ESM 0.956 0.958 0.820 0.835 0.748 1.000
SeizeVar-RF 0.949 0.973 0.776 0.841 0.819 1.000
AlphaMissense 0.976 0.965 0.881 0.889 0.797 1.000
REVEL 0.975 0.958 0.886 0.822 0.774 0.909
MetaRNN 0.988 0.980 0.887 0.855 0.804 1.000
CADD 0.919 0.879 0.869 0.682 0.676 1.000
VEST4 0.956 0.949 0.850 0.845 0.786 1.000
PolyPhen-2 0.946 0.896 0.756 0.634 0.627 0.545

Read this table honestly: SeizeVar's individual pathogenicity heads do not beat AlphaMissense, REVEL or MetaRNN — the ESM head trails AlphaMissense by 0.02–0.06 on four of six cohorts. The density-aware hybrid ensemble reaches or exceeds the best single component on 6 of 6 cohorts, but the contribution of this work is mechanism direction, dynamic evidence and precision-oriented integration, not pathogenicity-AUROC dominance. The paper says so explicitly.

A methodological finding worth its own paper

Every tool appeared to lose 0.09–0.14 AUROC on the Reclassified-VUS cohort — the textbook signature of concept drift. It is an artifact. 31 of the 210 variants are nonsense; missense-only competitors silently return NaN on all of them, and only SeizeVar-ESM scores the complete set. On the missense subset (n = 176), SeizeVar-ESM reaches AUROC 0.959 — first of nine predictors, and identical to its Val-Temporal performance with no decline whatsoever.

The 210-variant Reclassified-VUS benchmark is released with paired-core indices and per-tool scores so others can avoid the same trap.

Robustness against the "ClinVar retrieval" critique

A reasonable objection to any ClinVar-trained model is that it retrieves rather than predicts. Four tests:

Stress test Result
Strip 6 most gene-correlated constraint features Δ = −0.039 (0.728 → 0.689)
Reduce to purely variant-level 26 features 0.675, still well above chance
Entirely unseen gene (KCNQ2 held out of all training) 0.796 [0.602, 0.943]
Leave-one-gene-family-out NMDA 0.702 · K⁺ 0.655 · Na⁺ 0.640
Remove all 5 dynamic-evidence features Δ = −0.017 (−0.044 on variants absent from ClinVar)
Retrain at 2018 cutoff (n = 742), test on 2024+ 0.943 [0.936, 0.950] — only Δ = +0.028 from 8 years of ClinVar growth
Genes with < 20 training P/LP variants (n = 372) 0.928 [0.894, 0.956] — only Δ = +0.05 from dense genes

Prospective application

Scoring all 29,293 panel VUS produced 4,708 consensus Likely-Pathogenic (16.1 %), 12,980 Likely-Benign (44.3 %) and 11,605 that remain uncertain. Highest yield: SCN1A 656 (51.7 % of its VUS pool), SCN2A 403, SCN8A 300, POLG 294, CACNA1A 257.

Of those, 1,582 fall in the six sodium-channel genes and 1,500 carry a definite mechanism direction — 679 LoF-leaning, 821 GoF-leaning. On the Reclassified-VUS benchmark the consensus rule achieved 78 % PPV (29 of 37 flagged candidates correct) versus AlphaMissense 41 %, REVEL 40 %, MetaRNN 29 %.

Scope. The direction call is a prediction of molecular mechanism, not a prescribing recommendation. A 25-cell threshold-sensitivity sweep (n = 6,129) shows < 10 % category change in the central region and no per-gene direction reversal at any tested cutoff.


Repository layout

seizevar/
├── 01_data/code/                ① ClinVar extraction, splits, augmentation, N×N leakage audit
├── 02_features/code/            ② 39-feature matrix — gene / residue / substitution / dynamic-evidence
├── 03_models/
│   ├── pathogenicity/code/      ③ Random-Forest + ESM-2 LoRA training and held-out inference
│   └── mechanism/code/          ③ Gain-vs-loss-of-function head
├── 05_vus_application/code/     ④ Prospective scoring of the 29,293-variant VUS pool
├── 06_competitors/code/         ⑤ AlphaMissense / REVEL / MetaRNN / dbNSFP / LoGoFunc benchmarking
├── data_release/                Curated prediction tables + MANIFEST with SHA-1 hashes
├── figures/scripts/             Paper figures (Fig 1–6, Fig S1–S7, response figures)
├── infra/                       Colab bridge for the GPU training path
└── notebooks/                   Data and feature pipeline walkthroughs

The 26 numbered revision sub-analyses under submit_clinical_genetics/revision1_2026-05-23/analyses/ are tracked with code and outputs — bootstrap CIs, calibration, hybrid ensemble, topology ablations, threshold sweeps, the external functional cohort, and the prospective ClinVar simulation, each mapped to a specific reviewer comment.


Installation

git clone https://github.com/pipi0616/seizevar.git
cd seizevar
python -m venv .venv && source .venv/bin/activate
pip install -r requirements.txt

Python 3.10+. The ESM-2 fine-tune is GPU-bound; everything else runs comfortably on CPU. Training and inference were developed on Apple Silicon (MPS) with a Colab A100/V100 path in infra/. Scoring 29 K variants with the ESM-LoRA head takes ~30–45 min on M-series MPS, ~5 min on an A100.

Running the pipeline

Stages are numbered and run in order; each reads the previous stage's outputs/. Scripts resolve paths relative to __file__, so the working directory does not matter.

# ① Extract ClinVar, build splits, prove no leakage
python 01_data/code/01_extract_clinvar.py
python 01_data/code/02_build_splits.py
python 01_data/code/03_augment_train.py
python 01_data/code/04_build_extra_valsets.py
python 01_data/code/05_leakage_audit.py          # emits PASS/WARN/FAIL summary.md

# ② Compute the 39-feature matrix
python 02_features/code/compute_features.py
python 02_features/code/audit_features.py

# ③ Train
python 03_models/pathogenicity/code/train_rf_pathogenicity.py
python 03_models/pathogenicity/code/train_esm_lora.py        # GPU recommended
python submit_clinical_genetics/revision1_2026-05-23/analyses/11_external_functional_cohort/code/r12a_train_v5.py   # deployed mechanism head

# ④ Score the prospective VUS pool
python 05_vus_application/code/predict_vus_rf.py
python 05_vus_application/code/predict_vus_esm.py
python 05_vus_application/code/merge_vus_predictions.py

# ⑤ Benchmark against nine competitors
python 06_competitors/code/fetch_am_for_vus.py
python 06_competitors/code/extract_logofunc.py               # requires tabix
python 06_competitors/code/compute_auroc.py
python 06_competitors/code/mechanism_benchmark.py

bash infra/pack_for_colab.sh builds the Colab bundle. python data_release/build_release.py regenerates the curated tables.

Data and model availability

This repository holds source code and small reference tables only. Inputs (~454 MB) and trained weights (~628 MB, the ESM-LoRA checkpoint alone is 594 MB) are on Zenodo — unpack both at the repository root before running.

Location
Source code archive 10.5281/zenodo.19912082
Data + trained weights 10.5281/zenodo.19912547 (CC BY 4.0, 762 MB)
Per-variant predictions data_release/predictions_vus_all.csv — 29,293 rows
Reclassified-VUS benchmark data_release/ — 210 variants with per-tool scores

Upstream inputs are subject to their own terms: ClinVar, gnomAD, UniProt, AlphaFold v4, dbNSFP v4.5, MaveDB, and SCN2A electrophysiology from Northwestern PRISM (CC-BY 4.0).

Reproducibility notes

Every headline number is bound to a code output path in an internal 728-line numbers registry, verified by a linter before submission. All seeds fixed (random_state=42); headline AUROCs carry 2,000-iteration bootstrap CIs. Cross-validation is gene-grouped throughout, escalating to leave-one-gene-out and leave-one-gene-family-out for the mechanism task. The leakage audit emits a machine-checkable PASS/WARN/FAIL report covering variant-level overlap, whole-gene hold-out integrity, internal duplicates and position-level overlap.

Negative results are shipped, not hidden: ESM mutant-minus-wild-type delta embeddings were tested for the mechanism head, gave no gain, and are documented as a null.

Citation

Ye S, Chen P. From pathogenicity to mechanism: a variant interpretation framework for monogenic epilepsy. Clinical Genetics, 2026. (Accepted, MS 6385742)

@article{ye2026seizevar,
  title   = {From pathogenicity to mechanism: a variant interpretation
             framework for monogenic epilepsy},
  author  = {Ye, Shaopei and Chen, Peng},
  journal = {Clinical Genetics},
  year    = {2026},
  note    = {Accepted, MS 6385742}
}

Authors

Shaopei Ye — Universitat de Barcelona, Barcelona, Spain · ORCID 0009-0006-3598-0469 Peng Chen — Dept. of Otorhinolaryngology, The Second People's Hospital of Meishan, Sichuan, China · ORCID 0009-0001-1105-6638

License

MIT — see LICENSE. Prediction tables in data_release/ are CC BY 4.0.

Related work

  • HEARS-LM — hearing-loss-specialized protein language model
  • OtoVCE — LLM evidence layer for hereditary hearing loss
  • QAFI — quantitative functional-impact prediction


🌏 中文版本 / Chinese Version

点击下方展开完整中文文档 · Click below to expand the full Chinese documentation

👉 SeizeVar 中文说明(点击展开 / Click to expand)


SeizeVar 中文说明

面向单基因遗传性癫痫的机制感知型变异解读框架

不只回答"这个变异致病吗",而是回答"它把离子通道往哪个方向破坏了"

English 简体中文


一句话概括

SeizeVar 是一个覆盖 49 个癫痫相关基因的变异解读框架,回答了主流预测工具无法回答的问题:效应方向。它将致病性共识层(随机森林 + ESM-2 LoRA 交叉注意力)与功能获得/功能缺失(GoF/LoF)机制分类器、钠通道方向判定规则结合,并在留一法(leave-one-out)动态证据设计下打分——推理时实时读取最新的 ClinVar 参考集。

规模 端到端处理 45,143 个变异——训练 4,576、验证 11,274、前瞻性打分 VUS 29,293
核心结果 机制方向在留一基因法下 AUROC 0.736,而 AlphaMissense 仅 0.428、REVEL 0.502、MetaRNN 0.401
精确度 双头共识规则达到 78% 阳性预测值(PPV),约为单一工具(29–41%)的两倍
产出 4,708 个共识"可能致病"VUS 判定,其中 1,500 个带明确机制方向
状态 已被 Clinical Genetics(Wiley)接收,稿号 MS 6385742

研究问题

错义变异致病性预测已接近饱和:多个工具在专家审校的 ClinVer 数据上 AUROC 超过 0.97。然而在这个 49 基因癫痫panel中,仍有 29,293 条经典错义 ClinVar 记录停留在"意义未明"。单纯的判别能力已不再是瓶颈。

临床真正需要的下一步是方向。在单基因癫痫中这绝非学术问题:SCN1A 功能缺失变异通常会被钠通道阻滞剂加重,而 SCN8A 和新生儿期 SCN2A 功能获得变异恰恰是高剂量阻滞剂的公认适应证。同一基因家族,相反的分子后果,相反的临床含义。

通用预测器在结构上无法提供帮助。它们只输出单一致病概率;当这个概率被挪用来区分功能获得与功能缺失时,其表现处于或低于随机水平——AlphaMissense 0.428,MetaRNN 0.401。SeizeVar 把"方向"当作独立的监督学习任务,并补上第二项缺失能力:证据随 ClinVar 增长而自动更新,无需重新训练

为什么需要 SeizeVar?

能力 常规变异预测器 SeizeVar
致病性打分 有——随机森林 + ESM-2 LoRA 双头共识
机制方向(功能获得 vs 缺失) 有——独立的监督分类器
证据随 ClinVar 更新 无——训练时即冻结 有——5 个特征在推理时重新计算
针对癫痫基因优化 少见 有——49 基因 ClinGen 审校 panel
候选变异的精确度 阳性预测值 29–41% 78% 阳性预测值(双头共识门控)
发布前瞻性预测结果 少见 29,293 个已打分 VUS,其中 1,500 个带方向判定

两层设计是有意为之:致病性与机制是两个由不同类型证据支撑的独立生物学问题,单一分数无法同时回答"这个变异有害吗"和"它把通道功能往哪个方向推"。

适用范围说明。 本项目早期版本曾把输出描述为处方建议。经同行评审后,其定位已明确为分子机制预测。效应方向只是临床医生可以纳入考量的一项输入,它本身不是治疗建议,本仓库也不提供治疗建议。


工作原理

flowchart TB
    subgraph S1["① 数据 — 01_data/"]
        A["ClinVar 2026年4月版<br/>GRCh38 · 49 基因 panel"] --> B["严格 HGVSp 错义过滤<br/>保留 42,862 · 剔除 53.8%"]
        B --> C["按证据等级划分"]
        C --> D1["训练集 4,576<br/>6 个来源去重融合"]
        C --> D2["6 个验证队列<br/>共 11,274 变异"]
        C --> D3["VUS 池<br/>29,293"]
        D1 & D2 & D3 --> E["N×N 泄漏审计<br/>零变异重叠 · KCNQ2 整基因留出"]
    end

    subgraph S2["② 特征 — 02_features/"]
        E --> F["39 个特征 · 4 个层级"]
        F --> F1["基因层 ×8"]
        F --> F2["残基层 ×20<br/>AlphaFold v4 结构"]
        F --> F3["替换层 ×6"]
        F --> F4["动态证据层 ×5<br/>强制留一法"]
        F1 & F2 & F3 & F4 --> FM["39 维特征矩阵<br/>按训练集统计量标准化"]
    end

    subgraph S3["③ 模型 — 03_models/"]
        FM --> G1["随机森林<br/>按基因分组 5 折交叉验证"]
        FM --> G2["ESM-2 t30 + LoRA<br/>交叉注意力融合"]
        FM --> G3["机制头 v5<br/>XGBoost · GoF vs LoF"]
    end

    subgraph S4["④ 决策层 — 05_vus_application/"]
        G1 & G2 --> H["共识门控<br/>两个头必须一致"]
        G3 --> I["方向判定<br/>P·GoF 阈值"]
        H --> J["可能致病 / 可能良性 / 未定"]
        I --> K["钠通道方向规则<br/>6 个基因 · 确定性"]
        J & K --> L["29,293 个 VUS 完成打分<br/>4,708 判为可能致病"]
    end

    subgraph S5["⑤ 基准测试 — 06_competitors/"]
        L --> M["9 个工具 × 6 个队列<br/>AlphaMissense · REVEL · MetaRNN<br/>CADD · VEST4 · SIFT · PolyPhen-2 · LoGoFunc"]
    end

    style S1 fill:#e8f4f8,stroke:#4a90a4
    style S2 fill:#f0f4e8,stroke:#7a9a4a
    style S3 fill:#f8e8f0,stroke:#a44a7a
    style S4 fill:#fff4e0,stroke:#c49a3a
    style S5 fill:#eeeeee,stroke:#888888
Loading

致病性预测头:让语言模型与表格特征对话

ESM-2 预测头的设计比"把嵌入向量拼接到特征向量后面"更有意思。突变残基的嵌入作为 Query,而 39 个人工设计的生物学特征各自成为一个 Key/Value token,并配有可学习的"特征身份"嵌入。因此模型学到的是针对当前这个变异,应该关注哪一类生物学证据,且注意力权重可以直接读出。

flowchart LR
    A["蛋白序列<br/>+ 突变位点"] --> B["ESM-2 t30<br/>150M 参数 · 冻结"]
    B -.->|"LoRA r=8, α=16<br/>最后 2 层 q/k/v/out"| B
    B --> C["突变残基嵌入<br/>640 维"]
    C --> D["线性层 → 128 维<br/>QUERY"]

    E["39 个生物学特征"] --> F["逐特征线性层 1→128<br/>+ 可学习特征身份嵌入"]
    F --> G["39 个 token × 128 维<br/>KEYS / VALUES"]

    D --> H["多头交叉注意力<br/>4 头 · 残差 · LayerNorm"]
    G --> H
    H --> I["MLP 128→64→1"]
    I --> J["致病性 logit"]
    H --> K["注意力权重 [1×39]<br/>→ 可解释性"]

    style B fill:#e8f0f8,stroke:#4a7aa4
    style H fill:#f8e8e8,stroke:#a44a4a
    style K fill:#f0f8e8,stroke:#6a9a4a
Loading

可训练参数仅约 13.2 万,占骨干网络的 0.09%(其中约 8.2 万是 LoRA 适配器本身)。可解释性的收益是具体的:模型在离子通道基因判定上把 27.6% 的注意力放在通道拓扑特征上,而非通道基因仅 20.4%——说明它学到的是相关生物学,而不是死记基因身份。

动态证据:无需重训练即可自我改进

39 个特征中有 5 个是在推理时针对最新 ClinVar 参考集实时计算的——例如在 AlphaFold 结构中到最近致病变异的空间距离(SPDV)。训练阶段这些特征严格执行留一法(查询变异本身、以及等位基因系列中它自己的替换形式,都从参考集中移除),并经过实证验证:3,153 个训练集致病变异中,SPDV = 0 的数量为 0

其结果是:一条新的 ClinVar 提交会直接抬高邻近 VUS 的分数,完全不需要重新训练。在架构固定不变的前提下,跨 6 个历史 ClinVar 快照(2015–2026)测量显示,SPDV 中位数从 26 Å 降到 7 Å,随机森林在阈值 0.7 处的召回率从 19% 升到 34%,到 2024 年已接近冻结版 AlphaMissense(36%)。这 5 个特征只占特征总数的 13%,却在离子通道基因上承担了 35.9% 的 Gini 重要性


结果

机制方向——本工作的核心贡献

在 22 个 panel 基因、415 个具有膜片钳/双电极电压钳(TEVC)金标准的变异上,采用严格留一基因法(每个测试变异都由一个从未见过该基因任何变异的模型打分):

方法 AUROC 覆盖率
SeizeVar 机制头 v5 0.736 [0.685, 0.785] 100%
ESM-2 零样本 LLR 0.693 100%
MetaRNN 0.623 27%
AlphaMissense 0.601 89%
CADD / SIFT / REVEL / PolyPhen-2 / VEST4 0.55 – 0.61 约 25%

与唯一的专用 GoF/LoF 工具 LoGoFunc 在 329 个有基因组坐标的变异上正面对比:SeizeVar 0.770 [0.719, 0.819] vs LoGoFunc 0.760 [0.707, 0.817]——统计上无差异,但 SeizeVar 覆盖全部 22 个 panel 基因。

关于诚实数字的说明。 本工作早期版本曾报告机制 AUROC 为 0.905。自查发现评估集中 86% 的变异其实已在训练集中出现过。去除污染后的留一基因法估计值是 0.736,这也是本仓库和已接收论文所报告的数字。同样,LoGoFunc 表面上的 0.821 被追溯为 85 个变异的训练重叠假象,修正为 0.760。被推翻的数字保留在修订记录中,而非悄悄删除。

机制头还获得三项独立验证,均非来自 ClinVar:

  • 回顾性药物反应(n = 110):当预测机制与钠通道阻滞剂反应方向一致时,符合率为 43%(37/86);不一致时仅 4%(1/24)——OR 17.37 [2.24, 134.52],p = 1.5 × 10⁻⁴。
  • KCNQ2 深度突变扫描(n = 301,该基因完全未参与训练):Spearman ρ = +0.378,p = 1.1 × 10⁻¹¹。
  • TSC2 深度突变扫描(n = 8,166):100% 判为 LoF 倾向,最大 P(GoF) = 0.205——对一个抑癌基因而言完全正确。

值得注意的是,基因身份并未作为输入特征,但每个基因的方向倾向仍自发涌现:SCN1A/SCN1B 偏 LoF,SCN2A/SCN3A/SCN8A 偏 GoF。

致病性——6 个队列、9 个工具

Val-ClinGen Val-Temporal 重分类 VUS KCNQ2-DMS TSC2-MAVE SCN2A 电生理
SeizeVar 混合模型(部署版) 0.989 0.990 0.965 0.898 0.855 1.000
SeizeVar-ESM 0.956 0.958 0.820 0.835 0.748 1.000
SeizeVar-RF 0.949 0.973 0.776 0.841 0.819 1.000
AlphaMissense 0.976 0.965 0.881 0.889 0.797 1.000
REVEL 0.975 0.958 0.886 0.822 0.774 0.909
MetaRNN 0.988 0.980 0.887 0.855 0.804 1.000
CADD 0.919 0.879 0.869 0.682 0.676 1.000
VEST4 0.956 0.949 0.850 0.845 0.786 1.000
PolyPhen-2 0.946 0.896 0.756 0.634 0.627 0.545

请诚实地阅读这张表: SeizeVar 的单个致病性预测头并未超越 AlphaMissense、REVEL 或 MetaRNN——ESM 头在 6 个队列中的 4 个上落后 AlphaMissense 0.02–0.06。密度感知的混合集成在 6/6 个队列上达到或超过最佳单一组件,但本工作的贡献在于机制方向、动态证据和以精确度为导向的整合,而非致病性 AUROC 的领先。论文中对此有明确表述。

一项本身即可独立成文的方法学发现

所有工具在"重分类 VUS"队列上都出现了 0.09–0.14 的 AUROC 下降——这是概念漂移的教科书式特征。但它其实是假象。210 个变异中有 31 个是无义变异;仅支持错义的竞品工具在这些变异上静默返回 NaN,只有 SeizeVar-ESM 能对完整集合打分。在错义子集(n = 176)上,SeizeVar-ESM 达到 AUROC 0.959——9 个预测器中排名第一,且与其 Val-Temporal 表现完全一致,没有任何下降。

我们发布了这个 210 变异的 Reclassified-VUS 基准集,附带配对核心索引和各工具得分,供他人避开同一陷阱。

针对"只是在检索 ClinVar"质疑的稳健性验证

对任何基于 ClinVar 训练的模型,一个合理的质疑是:它是在检索而非预测。四类压力测试:

压力测试 结果
剔除 6 个与基因最相关的约束特征 Δ = −0.039(0.728 → 0.689)
仅保留纯变异层面的 26 个特征 0.675,仍远高于随机
完全未见过的基因(KCNQ2 全程留出) 0.796 [0.602, 0.943]
留一基因家族法 NMDA 0.702 · 钾通道 0.655 · 钠通道 0.640
移除全部 5 个动态证据特征 Δ = −0.017(在 ClinVar 中不存在的变异上为 −0.044)
用 2018 年截点重训(n = 742),测 2024+ 0.943 [0.936, 0.950]——8 年 ClinVar 增长仅带来 Δ = +0.028
训练致病变异 < 20 个的基因(n = 372) 0.928 [0.894, 0.956]——与高密度基因差距仅 Δ = +0.05

前瞻性应用

对全部 29,293 个 panel VUS 打分后,得到 4,708 个共识"可能致病"(16.1%)、12,980 个"可能良性"(44.3%),另有 11,605 个仍未定。产出最高的基因:SCN1A 656(占其 VUS 池的 51.7%)、SCN2A 403、SCN8A 300、POLG 294、CACNA1A 257。

其中 1,582 个落在 6 个钠通道基因内,1,500 个带有明确机制方向——679 个 LoF 倾向、821 个 GoF 倾向。在 Reclassified-VUS 基准上,共识规则达到 78% PPV(37 个候选中 29 个正确),而 AlphaMissense 为 41%、REVEL 40%、MetaRNN 29%。

适用范围声明。 方向判定是对分子机制的预测,不是处方建议。一项 25 格阈值敏感性扫描(n = 6,129)显示中心区域类别变动 < 10%,且在所有测试阈值下均未出现任何基因的方向翻转


仓库结构

seizevar/
├── 01_data/code/                ① ClinVar 提取、数据划分、增广、N×N 泄漏审计
├── 02_features/code/            ② 39 特征矩阵——基因 / 残基 / 替换 / 动态证据
├── 03_models/
│   ├── pathogenicity/code/      ③ 随机森林 + ESM-2 LoRA 训练与留出集推理
│   └── mechanism/code/          ③ 功能获得 vs 功能缺失机制头
├── 05_vus_application/code/     ④ 对 29,293 个 VUS 池的前瞻性打分
├── 06_competitors/code/         ⑤ AlphaMissense / REVEL / MetaRNN / dbNSFP / LoGoFunc 基准测试
├── data_release/                精选预测表 + 带 SHA-1 校验的 MANIFEST
├── figures/scripts/             论文图表(图 1–6、图 S1–S7、审稿回应图)
├── infra/                       GPU 训练路径的 Colab 桥接
└── notebooks/                   数据与特征流程演示

submit_clinical_genetics/revision1_2026-05-23/analyses/ 下的 26 个编号修订子分析连同代码输出一并纳入版本控制——自助法置信区间、校准度、混合集成、拓扑消融、阈值扫描、外部功能队列、前瞻性 ClinVar 模拟,每一项都对应到具体的审稿意见。


安装

git clone https://github.com/pipi0616/seizevar.git
cd seizevar
python -m venv .venv && source .venv/bin/activate
pip install -r requirements.txt

需要 Python 3.10+。ESM-2 微调依赖 GPU,其余部分在 CPU 上运行流畅。开发环境为 Apple Silicon(MPS),并在 infra/ 中提供 Colab A100/V100 路径。用 ESM-LoRA 头对 2.9 万个变异打分,在 M 系列 MPS 上约需 30–45 分钟,A100 上约 5 分钟

运行流程

各阶段按编号顺序执行,每一步读取上一步的 outputs/。脚本基于 __file__ 解析路径,因此工作目录无关紧要。

# ① 提取 ClinVar、构建划分、验证无泄漏
python 01_data/code/01_extract_clinvar.py
python 01_data/code/02_build_splits.py
python 01_data/code/03_augment_train.py
python 01_data/code/04_build_extra_valsets.py
python 01_data/code/05_leakage_audit.py          # 输出 PASS/WARN/FAIL 审计报告

# ② 计算 39 维特征矩阵
python 02_features/code/compute_features.py
python 02_features/code/audit_features.py

# ③ 训练
python 03_models/pathogenicity/code/train_rf_pathogenicity.py
python 03_models/pathogenicity/code/train_esm_lora.py        # 建议使用 GPU
python submit_clinical_genetics/revision1_2026-05-23/analyses/11_external_functional_cohort/code/r12a_train_v5.py   # 部署版机制头

# ④ 对前瞻性 VUS 池打分
python 05_vus_application/code/predict_vus_rf.py
python 05_vus_application/code/predict_vus_esm.py
python 05_vus_application/code/merge_vus_predictions.py

# ⑤ 与 9 个竞品工具基准对比
python 06_competitors/code/fetch_am_for_vus.py
python 06_competitors/code/extract_logofunc.py               # 需要 tabix
python 06_competitors/code/compute_auroc.py
python 06_competitors/code/mechanism_benchmark.py

bash infra/pack_for_colab.sh 构建 Colab 打包文件;python data_release/build_release.py 重新生成精选结果表。

数据与模型获取

本仓库仅包含源代码和小型参考表。输入数据(约 454 MB)与训练好的权重(约 628 MB,其中 ESM-LoRA 检查点单独就有 594 MB)存放在 Zenodo——运行前请将两者解压到仓库根目录。

内容 位置
源代码存档 10.5281/zenodo.19912082
数据 + 训练权重 10.5281/zenodo.19912547(CC BY 4.0,762 MB)
逐变异预测结果 data_release/predictions_vus_all.csv——29,293 行
Reclassified-VUS 基准集 data_release/——210 个变异,含各工具得分

上游数据各自遵循其原始许可条款:ClinVar、gnomAD、UniProt、AlphaFold v4、dbNSFP v4.5、MaveDB,以及来自 Northwestern PRISM 的 SCN2A 电生理数据(CC-BY 4.0)。

可复现性说明

论文中每一个关键数字都绑定到具体的代码输出路径,由一份 728 行的内部"数字登记表"管理,并在投稿前经脚本校验。所有随机种子固定(random_state=42);关键 AUROC 均附 2,000 次自助法置信区间。交叉验证全程按基因分组,机制任务进一步升级为留一基因法和留一基因家族法。泄漏审计输出机器可校验的 PASS/WARN/FAIL 报告,覆盖变异层重叠、整基因留出完整性、内部重复和位点层重叠。

阴性结果同样公开而非隐藏:机制头曾尝试引入 ESM"突变减野生型"差值嵌入,未带来任何增益,作为阴性结果如实记录。

引用

Ye S, Chen P. From pathogenicity to mechanism: a variant interpretation framework for monogenic epilepsy. Clinical Genetics, 2026.(已接收,稿号 MS 6385742)

@article{ye2026seizevar,
  title   = {From pathogenicity to mechanism: a variant interpretation
             framework for monogenic epilepsy},
  author  = {Ye, Shaopei and Chen, Peng},
  journal = {Clinical Genetics},
  year    = {2026},
  note    = {Accepted, MS 6385742}
}

作者

叶少培(Shaopei Ye) — 巴塞罗那大学(Universitat de Barcelona),西班牙巴塞罗那 · ORCID 0009-0006-3598-0469 陈鹏(Peng Chen) — 四川省眉山市第二人民医院耳鼻咽喉科 · ORCID 0009-0001-1105-6638

许可证

MIT——详见 LICENSEdata_release/ 中的预测结果表采用 CC BY 4.0。

相关项目

  • HEARS-LM — 听力损失专用蛋白质语言模型
  • OtoVCE — 遗传性听力损失的大模型证据抽取层
  • QAFI — 变异功能影响的定量预测

About

Mechanism-aware variant interpretation pipeline for monogenic epilepsy: RF + ESM-2 LoRA pathogenicity heads, gain-vs-loss-of-function mechanism classifier, sodium-channel prescribing rule, leave-one-out dynamic-evidence framework.

Topics

Resources

Stars

Watchers

Forks

Releases

Packages

Contributors

Languages