ExpansionHunter — 短串联重复序列扩展检测工具¶

# 方法1：conda安装（v5.0.0，推荐）
conda create -n expansionhunter python=3.9
conda activate expansionhunter
conda install -c bioconda expansionhunter   # 安装最新版（v5.0.0）

# 验证安装
ExpansionHunter --version    # 查看版本号

# 方法2：下载预编译二进制（Linux/macOS）
# 访问 https://github.com/Illumina/ExpansionHunter/releases 下载
wget https://github.com/Illumina/ExpansionHunter/releases/download/v5.0.0/ExpansionHunter-v5.0.0-linux_x86_64.tar.gz
tar xzf ExpansionHunter-v5.0.0-linux_x86_64.tar.bz2
./ExpansionHunter-v5.0.0-linux_x86_64/bin/ExpansionHunter --version

# 获取变异目录文件（STR位点定义）
# 官方提供hg19/hg38的STR变异目录
# 位于 ExpansionHunter 安装目录的 variant_catalog/ 子目录下
ls ExpansionHunter-v5.0.0-linux_x86_64/variant_catalog/
# hg19/  hg38/  - 各含 variant_catalog.json（常见致病STR位点）

# 最常用：用官方变异目录检测全部已知致病STR
ExpansionHunter \
    --reads sample.bam \              # 输入BAM或CRAM文件（需排序和索引）
    --reference reference.fasta \    # 参考基因组（需.fai索引）
    --variant-catalog variant_catalog.json \  # STR位点定义文件
    --output-prefix sample_EH \       # 输出文件前缀
    --threads 8                       # 线程数（v5.0支持多线程）

# 流式分析（默认，适合少量位点<1000）
ExpansionHunter \
    --reads sample.bam \
    --reference reference.fasta \
    --variant-catalog variant_catalog.json \
    --output-prefix sample_EH \
    --analysis-mode streaming          # 默认模式

# 寻路分析（大型目录，>10k位点时用，速度更快）
ExpansionHunter \
    --reads sample.bam \
    --reference reference.fasta \
    --variant-catalog large_catalog.json \
    --output-prefix sample_EH \
    --analysis-mode seeking \          # 寻路模式，适合大目录
    --threads 16                       # 多线程加速

# REViewer生成每个STR位点的reads可视化图（SVG格式）
# 先安装REViewer
conda install -c bioconda reviewer

# 生成可视化（需要ExpansionHunter的BAM输出）
REViewer \
    --reads sample_EH_realigned.bam \  # ExpansionHunter输出的BAM
    --vcf sample_EH.vcf \             # ExpansionHunter输出的VCF
    --reference reference.fasta \
    --catalog variant_catalog.json \
    --locus HTT \                      # 指定查看HTT（亨廷顿）位点
    --output-prefix HTT_review

# 完整流程：检测神经退行性疾病相关STR扩增

REF="hg38.fasta"
CATALOG="ExpansionHunter-v5.0.0/variant_catalog/hg38/variant_catalog.json"

# 1. 运行ExpansionHunter
ExpansionHunter \
    --reads patient_sample.bam \
    --reference $REF \
    --variant-catalog $CATALOG \
    --output-prefix patient_STR \
    --sex female \               # 指定样本性别（影响X连锁位点如FMR1）
    --threads 16

# 2. 查看输出文件
ls patient_STR.*
# patient_STR.vcf           - VCF格式结果（每个STR位点的基因型）
# patient_STR.json          - JSON格式详细结果
# patient_STR_realigned.bam - 重新比对的reads（供REViewer可视化）

# 3. 查看关键致病位点结果
# 常见致病STR及其正常/致病阈值：
# HTT   (亨廷顿病):     正常 ≤36 CAG，致病 ≥40
# C9ORF72 (ALS/FTD):    正常 ≤30 GGGGCC，致病 ≥200
# FMR1  (脆性X综合征):  正常 ≤45 CGG，致病 ≥200
# ATXN3 (SCA3脊髓小脑): 正常 ≤44 CAG，致病 ≥55

# 提取HTT位点结果
grep "HTT" patient_STR.vcf | grep -v "^#"

# 4. 解读VCF中的GT字段
# 示例结果：
# HTT  GT:REPCN  1/1:17/19   - 两条染色体分别17和19个CAG重复（正常）
# HTT  GT:REPCN  1/2:17/45   - 一条17个，一条45个CAG重复（致病！）

# 5. 提取扩增位点（重复次数大于某阈值）
python3 << 'EOF'
import json
with open("patient_STR.json") as f:   # 读取JSON结果
    data = json.load(f)               # 解析JSON
for locus in data["LocusResults"]:   # 遍历所有位点
    for variant in locus["Variants"].values():  # 查看每个变异
        if "RepeatCount" in variant:
            rc = variant["RepeatCount"]
            if max(rc) > 50:         # 超过50次重复则输出（示例阈值）
                print(locus["LocusId"], rc)  # 输出位点名和重复次数
EOF