Metadata-Version: 2.4
Name: seekrare
Version: 0.9.4
Summary: Rare disease candidate variant diagnosis system — three-stage VCF-to-ranking with optional Genos model and AlphaFold3 structure prediction
Author-email: Wang Zilu <1923524070@qq.com>
License: MIT
Project-URL: Homepage, https://github.com/wangz1lu/SeekRare
Project-URL: Documentation, https://github.com/wangz1lu/SeekRare#readme
Project-URL: Repository, https://github.com/wangz1lu/SeekRare
Project-URL: Issues, https://github.com/wangz1lu/SeekRare/issues
Project-URL: PyPI, https://pypi.org/project/seekrare/
Keywords: rare-disease,genomics,variant-analysis,llm,clinical-genetics,wgs,vcf,gtex,alphafold,genos
Classifier: Development Status :: 3 - Alpha
Classifier: Intended Audience :: Science/Research
Classifier: License :: OSI Approved :: MIT License
Classifier: Topic :: Scientific/Engineering :: Bio-Informatics
Classifier: Programming Language :: Python :: 3.10
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Operating System :: OS Independent
Requires-Python: >=3.10
Description-Content-Type: text/markdown
Requires-Dist: pandas>=2.0
Requires-Dist: numpy>=1.24
Requires-Dist: pyarrow>=14.0
Requires-Dist: pyyaml>=6.0
Requires-Dist: openai>=1.0
Requires-Dist: anthropic>=0.18
Requires-Dist: tqdm>=4.66
Requires-Dist: loguru>=0.7
Requires-Dist: requests>=2.31
Provides-Extra: dev
Requires-Dist: pytest>=7.4; extra == "dev"
Requires-Dist: pytest-cov>=4.1; extra == "dev"
Requires-Dist: black>=23.0; extra == "dev"
Requires-Dist: ruff>=0.1; extra == "dev"
Provides-Extra: alphafold
Requires-Dist: alphafold3; extra == "alphafold"
Provides-Extra: genomodel
Requires-Dist: torch>=2.0; extra == "genomodel"
Requires-Dist: transformers>=4.30; extra == "genomodel"
Requires-Dist: pygromaix>=0.1; extra == "genomodel"

# SeekRare

**Rare Disease Candidate Variant Diagnosis System — Three-Stage VCF-to-Ranking**
<img width="2332" height="1170" alt="seekrare_github" src="https://github.com/user-attachments/assets/7600bf6b-8852-4c6a-87d3-c878113133ce" />

---

## Table of Contents

- [Installation](#installation)
- [Architecture Overview](#architecture-overview)
- [Stage 1 — VCF Preprocessing & Basic Annotation](#stage-1--vcf-preprocessing--basic-annotation)
- [Stage 2 — Advanced Annotation (Optional, Multi-Step)](#stage-2--advanced-annotation-optional-multi-step)
- [Stage 3 — LLM Dynamic Scoring & Ranking](#stage-3--llm-dynamic-scoring--ranking)
- [Full Pipeline Usage](#full-pipeline-usage)
- [Output File Reference](#output-file-reference)

---

## Installation

```bash
pip install seekrare
```

**Reference files required for Stage 1 (user-provided):**

| File | Description |
|------|-------------|
| `ref_fasta` | Reference genome FASTA |
| `gtf_file` | Gene annotation GTF file |
| `clinvar_vcf` | ClinVar VCF.gz/csv |
| `dbsnp_vcf` | dbSNP common all chr VCF.gz (optional; used for common variant filtering) |

```bash
# GRCh38 / hg38
# ref_fasta
wget https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/405/GCF_000001405.40_GRCh38.p14/GCF_000001405.40_GRCh38.p14_genomic.fna.gz
gunzip GCF_000001405.40_GRCh38.p14_genomic.fna.gz
mv GCF_000001405.40_GRCh38.p14_genomic.fna GRCh38.fa

# gtf_file
wget https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/405/GCF_000001405.40_GRCh38.p14/GCF_000001405.40_GRCh38.p14_genomic.gtf.gz
gunzip GCF_000001405.40_GRCh38.p14_genomic.gtf.gz
mv GCF_000001405.40_GRCh38.p14_genomic.gtf GRCh38.gtf

# clinvar_vcf
wget https://ftp.ncbi.nlm.nih.gov/pub/clinvar/vcf_GRCh38/clinvar.vcf.gz
wget https://ftp.ncbi.nlm.nih.gov/pub/clinvar/vcf_GRCh38/clinvar.vcf.gz.tbi

# dbsnp_vcf
wget https://ftp.ncbi.nih.gov/snp/organisms/human_9606_b151_GRCh38p7/VCF/00-common_all.vcf.gz
wget https://ftp.ncbi.nih.gov/snp/organisms/human_9606_b151_GRCh38p7/VCF/00-common_all.vcf.gz.tbi

```

---

## Architecture Overview

### Core Design Principles

- **Stage 1 is required** — VCF normalization + basic annotation
- **Stage 2 is optional** — advanced annotation on top of Stage 1 (eQTL / SpliceVARDB / OMIM-HPO)
- **Stage 3 is required** — LLM dynamic scoring + ranking

### Mode Detection (Automatic)

| Parameters | Mode |
|-----------|------|
| `vcf_proband` only | Singleton mode |
| `vcf_proband` + `vcf_father` + `vcf_mother` | Family trio mode |

---


### Stage 1 Usage

**Option 1: via SeekRarePipeline (auto-detects singleton/family)**

```python
from seekrare import SeekRarePipeline

pipeline = SeekRarePipeline(
    vcf_proband="child.vcf.gz",
    vcf_father="father.vcf.gz",   # all three → family mode
    vcf_mother="mother.vcf.gz",
    ref_fasta="/path/to/GRCh38.fa",
    gtf_file="/path/to/genomic.gtf",
    clinvar_vcf="/path/to/clinvar.vcf.gz",
    dbsnp_vcf="/path/to/dbsnp.vcf.gz",   # optional
    work_dir="seekrare_output",
)
pipeline.stage1_preprocess()
# Output: seekrare_output/3_clinvar_annotated.csv
```

**Option 2: standalone function `run_family_preprocess`**

```python
from seekrare.preprocess import run_family_preprocess

result = run_family_preprocess(
    work_dir="seekrare_output",
    child_vcf="child.vcf.gz",
    father_vcf="father.vcf.gz",
    mother_vcf="mother.vcf.gz",
    ref_fasta="/path/to/hg38.fa",
    gtf_file="/path/to/GRCh38.gtf",
    clinvar_vcf="/path/to/clinvar.vcf.gz",
    dbsnp_vcf="/path/to/00-common_all.vcf.gz",
)
# result = {"modes": {...}, "combined": DataFrame, "output_csv": str}
```

**Option 3: individual step functions (each usable standalone)**

```python
from seekrare.preprocess import (
    stage1_vcf_to_gt_csv,       # VCF → GT CSV
    stage1_annotate_by_gtf,     # GT CSV → GTF annotation
    stage1_merge_filter_clinvar, # ClinVar annotation
    stage1_dbsnp_filter,        # dbSNP common filtering
    simplify_clinvar_csv,       # Column cleanup
)
```

---

## Stage 2 — Advanced Annotation (Optional, Multi-Step)

Stage 2 consists of multiple optional steps. Each is a standalone function — use only what you need:

```
Stage 1 output (3_clinvar_annotated.csv)
       │
  ┌────┴────────────────────────────┐
  │  step 1: stage2_eqtl_annotation  │
  │  GTEx eQTL tissue annotation    │
  └────┬────────────────────────────┘
       │ (if tissue_dir provided)
  ┌────┴────────────────────────────┐
  │  step 2: stage2_splicevardb_..│
  │  SpliceVARDB splice annotation │
  └────┬────────────────────────────┘
       │ (if splicevardb_tsv provided)
  ┌────┴────────────────────────────┐
  │  step 3: stage2_omim_hpo_...  │
  │  OMIM + HPO secondary annot.  │
  └────┴────────────────────────────┘
         ↓
  6_omim_hpo_annotated.csv (example)
```

### Step 1: GTEx eQTL Annotation (`stage2_eqtl_annotation`)

**Purpose**: LLM selects relevant tissues from GTEx based on symptoms; eQTL data merged from parquets.

```python
from seekrare import stage2_eqtl_annotation

stage2_eqtl_annotation(
    stage1_csv="seekrare_output/3_clinvar_annotated.csv",
    tissue_dir="/path/to/GTEx_v11_eQTL_parquets",
    symptoms="thalassemia",
    output_csv="seekrare_output/4_eqtl_annotated.csv",
    llm_model="deepseek-v4-flash",
    api_key="sk-xxxxxxxx",
    base_url="https://api.deepseek.com",
)
```

**New columns**: `eqtl_gene`, `eqtl_slope`, `eqtl_pval`, `eqtl_tissue`, `n_eqtl_tissues`

---

### Step 2: SpliceVARDB Annotation (`stage2_splicevardb_annotation`)

**Purpose**: Match variants by `hg38` column (format `chr-pos-ref-alt`) to SpliceVARDB, annotate `classification`.

```python
from seekrare import stage2_splicevardb_annotation

stage2_splicevardb_annotation(
    input_csv="seekrare_output/4_eqtl_annotated.csv",
    splicevardb_tsv="/path/to/splicevardb.download.tsv",
    output_csv="seekrare_output/5_splicevardb_annotated.csv",
)
```

**New column**: `splicevardb` (values: `Splice-altering` / `Low-frequency` / `Normal` / `Conflicting`)

---

### Step 3: OMIM + HPO Secondary Annotation (`stage2_omim_hpo_annotation`)

**Purpose**:
- `OMIM`: gene_name → MIM numbers from `genemap2.txt` (existing values preserved)
- `HPO`: append HPO tags from `phenotype.hpoa` by OMIM (existing values preserved)
- `diseasename`: supplement from `phenotype.hpoa` disease_name (existing values preserved)

```python
from seekrare import stage2_omim_hpo_annotation

stage2_omim_hpo_annotation(
    input_csv="seekrare_output/5_splicevardb_annotated.csv",
    genemap2_path="/path/to/2022_05_05-genemap2.txt",
    mimtitles_path="/path/to/mimTitles.txt",
    phenotype_hpoa_path="/path/to/phenotype.hpoa",
    output_csv="seekrare_output/6_omim_hpo_annotated.csv",
)
```

**Modified columns**: `OMIM` (fills blanks), `HPO` (appends new tags), `diseasename` (fills blanks)

---

### Via SeekRarePipeline

```python
pipeline = SeekRarePipeline(
    vcf_proband="child.vcf.gz",
    work_dir="seekrare_output",
    gtex_tissue_dir="/path/to/GTEx_v11_eQTL_parquets",
    splicevardb_tsv="/path/to/splicevardb.download.tsv",
    genemap2_path="/path/to/genemap2.txt",
    mimtitles_path="/path/to/mimTitles.txt",
    phenotype_hpoa_path="/path/to/phenotype.hpoa",
)
pipeline.stage1_preprocess()

pipeline.stage2_eqtl_annotation(symptoms="thalassemia")
pipeline.stage2_splicevardb_annotation()
pipeline.stage2_omim_hpo_annotation()
```

---

## Stage 3 — LLM Dynamic Scoring & Ranking

### Scoring System

**Static columns (built-in maps, no LLM needed):**

| Column | Mapping Rule |
|--------|-------------|
| `feature_type` | CDS=1.0, exon=0.9, gene=0.7, start_codon=0.8, stop_codon=0.8, transcript=0.5 |
| `significance` | `Benign/Likely_benign` → worst case (min) = -1.0; `Pathogenic` → 1.0; `/`-split takes worst end; empty → -0.5 |
| `clinvarstar` | 0–5 directly mapped |
| `eqtl_tissue` | has content=0.5, empty=0 |
| `splicevardb` | `Splice-altering`=1.0, `Low-frequency`=0.6, `Conflicting`=0.5, `Normal`=0.2 |

**Dynamic columns (LLM scores each unique value based on symptoms):**

| Column | Description |
|--------|-------------|
| `gene_name` | LLM scores each gene by relevance to symptoms |
| `diseasename` | LLM scores each disease name by relevance to symptoms |
| `HPO` | HP:xxxx tags scored by semantic relevance to symptoms |
| `OMIM` | OMIM:xxxxx scored by relevance to symptoms |
| `Orphanet` | Orphanet:xxxxx scored by relevance to symptoms |
| `inheritance_mode` | LLM decides which of de_novo / recessive / xlinked fits best given disease characteristics |
| `MC` | SO:xxxx tags scored by functional impact |

**Consistency bonus**: if both `gene_name`>0.6 AND `diseasename`>0.6 → extra +0.1. Rewards variants where the gene and disease name both strongly match the suspected diagnosis.

**Default weight distribution** (LLM may adjust based on symptoms):

```json
{
  "gene_name": 0.15, "HPO": 0.15, "OMIM": 0.15,
  "Orphanet": 0.08, "diseasename": 0.12,
  "inheritance_mode": 0.10, "MC": 0.05,
  "feature_type": 0.10, "significance": 0.10,
  "clinvarstar": 0.025, "eqtl_tissue": 0.025, "splicevardb": 0.025
}
```

### Usage

**Option 1: standalone function**

```python
from seekrare import stage3_score_and_rank

top = stage3_score_and_rank(
    csv_path="seekrare_output/6_omim_hpo_annotated.csv",
    symptoms="thalassemia",
    top_k=50,
    output_csv="seekrare_output/stage3_ranked.csv",
    llm_model="deepseek-v4-flash",
    api_key="sk-xxxxxxxx",
    base_url="https://api.deepseek.com",
)
print(top[["CPRA", "gene_name", "seekrare_score", "rank"]].head(10))
```

**Option 2: via SeekRarePipeline**

```python
pipeline.stage3_score_and_rank(symptoms="thalassemia")
```

**Option 3: use Stage3Scorer directly**

```python
from seekrare.scoring import Stage3Scorer

scorer = Stage3Scorer(
    csv_path="seekrare_output/6_omim_hpo_annotated.csv",
    symptoms="thalassemia",
    top_k=50,
    llm_model="deepseek-v4-flash",
    api_key="sk-xxxxxxxx",
)
df = scorer.run()
scorer.save(df, "seekrare_output/stage3_ranked.csv")
```

### Stage 3 Output Columns

```
All original columns + seekrare_score (weighted total) + rank (1-based)
```

---

## Stage 4 — Genos / AlphaFold3 (Optional)

### Stage 4A: Genos Model Analysis

**Purpose**: Run Genos model prediction on Stage 3 ranked top-K sites.

```python
from seekrare import stage4_genos_analysis

stage4_genos_analysis(
    sites="top:10",    # take top 10 from stage3_ranked.csv
    stage3_csv="seekrare_output/stage3_ranked.csv",
    genome_fa="/path/to/hg38.fa",
    model_path="/path/to/Genos-1.2B",
    output_dir="seekrare_output/genos_result",
)
```

**Note**: Genos model weights must be pre-deployed at `model_path`.

---

### Stage 4B: AlphaFold3 Structure Prediction

```python
from seekrare import stage4_alphafold_prediction

stage4_alphafold_prediction(
    csv_path="seekrare_output/stage3_ranked.csv",
    ref_fasta="/path/to/hg38.fa",
    gtf_file="/path/to/genomic.gtf",
    top_n=5,
    output_dir="seekrare_output/alphafold_results",
)
```

---

## Full Pipeline Usage

```python
from seekrare import SeekRarePipeline

# ── Config ─────────────────────────────────────────────
pipeline = SeekRarePipeline(
    # Stage 1
    vcf_proband="/path/to/proband.vcf.gz",
    vcf_father="/path/to/father.vcf.gz",
    vcf_mother="/path/to/mother.vcf.gz",
    ref_fasta="/path/to/hg38.fa",
    gtf_file="/path/to/genomic.gtf",
    clinvar_vcf="/path/to/clinvar.vcf.gz",
    dbsnp_vcf="/path/to/dbsnp.vcf.gz",
    work_dir="seekrare_output",
    # Stage 2 (all optional)
    gtex_tissue_dir="/path/to/GTEx_v11_eQTL_parquets",
    splicevardb_tsv="/path/to/splicevardb.download.tsv",
    genemap2_path="/path/to/genemap2.txt",
    mimtitles_path="/path/to/mimTitles.txt",
    phenotype_hpoa_path="/path/to/phenotype.hpoa",
    # Stage 3 LLM
    llm_model="deepseek-v4-flash",
    api_key="sk-xxxxxxxx",
    base_url="https://api.deepseek.com",
)

# ── Stage 1 ─────────────────────────────────────────────
pipeline.stage1_preprocess()

# ── Stage 2 (multi-step optional) ───────────────────────
pipeline.stage2_eqtl_annotation(symptoms="thalassemia")
pipeline.stage2_splicevardb_annotation()
pipeline.stage2_omim_hpo_annotation()

# ── Stage 3 ─────────────────────────────────────────────
result = pipeline.stage3_score_and_rank(symptoms="thalassemia")
print(result[["CPRA", "gene_name", "seekrare_score", "rank"]].head(20))

# ── Stage 4 (optional) ───────────────────────────────────
pipeline.stage4_genos_analysis(
    sites="top:10",
    genome_fa="/path/to/hg38.fa",
    model_path="/path/to/Genos-1.2B",
    output_dir="seekrare_output/genos_result",
)
```

---

## Output File Reference

| File | Stage | Description |
|------|-------|-------------|
| `3_clinvar_annotated.csv` | End of Stage 1 | Family/singleton merged result with `inheritance_mode` |
| `4_eqtl_annotated.csv` | End of Stage 2 Step 1 | After GTEx eQTL annotation |
| `5_splicevardb_annotated.csv` | End of Stage 2 Step 2 | After SpliceVARDB annotation |
| `6_omim_hpo_annotated.csv` | End of Stage 2 Step 3 | After OMIM + HPO secondary annotation |
| `stage3_ranked.csv` | End of Stage 3 | Sorted results with `seekrare_score` and `rank` |
| `genos_result/` | Stage 4A | Genos model prediction output directory |
| `alphafold_results/` | Stage 4B | AlphaFold3 structure prediction output directory |

---

## Dependencies

**Core** (`pip install seekrare`):
- pandas, numpy, pyarrow, pyyaml
- openai, anthropic
- tqdm, loguru, requests

**Optional groups**:
- `[alphafold]` — alphafold3
- `[genomodel]` — torch, transformers, pygromaix
- `[dev]` — pytest, black, ruff
