长读测序基因组从头组装¶

# ===== 长读测序基因组从头组装流程 =====

# ===== 方案1：HiFi 数据组装（最推荐） =====

# 1. 使用 hifiasm 组装
# -o: 输出前缀
# -t: 线程数
# --primary: 输出主要组装（非单倍型模式）
hifiasm -o asm -t 32 hifi_reads.fastq.gz

# 转换 GFA 为 FASTA
awk '/^S/{print ">"$2; print $3}' \
    asm.bp.p_ctg.gfa > primary_assembly.fa

# 2. 单倍型组装（如果有父母数据或 Hi-C）
# 方法A：用 trio binning（需要父母的短读数据）
# 先用 yak 建立父母 k-mer 数据库
yak count -t 16 -b 37 -o pat.yak paternal_illumina.fastq.gz
yak count -t 16 -b 37 -o mat.yak maternal_illumina.fastq.gz
# 再用 hifiasm trio 模式
hifiasm -o asm_trio -t 32 \
    -1 pat.yak -2 mat.yak \
    hifi_reads.fastq.gz
# 输出 hap1 和 hap2 两套组装

# 方法B：用 Hi-C 辅助单倍型分装
hifiasm -o asm_hic -t 32 \
    --h1 hic_R1.fastq.gz \
    --h2 hic_R2.fastq.gz \
    hifi_reads.fastq.gz

# ===== 方案2：Verkko 混合组装（T2T 级别） =====
# 需要 HiFi + ONT ultra-long 数据
verkko -d verkko_output \
    --hifi hifi_reads.fastq.gz \
    --nano ont_ultralong.fastq.gz \
    --threads 32

# ===== 方案3：ONT 数据组装 =====
# 使用 Flye
flye --nano-hq ont_reads.fastq.gz \
    --genome-size 3g \    # 预估基因组大小
    --out-dir flye_output \
    --threads 32

# ===== 组装质量评估 =====

# 1. 基本统计（QUAST）
quast primary_assembly.fa \
    -r GRCh38.fa \       # 参考基因组（如果有）
    -o quast_output/ \
    -t 16 \
    --large               # 大基因组模式

# 2. BUSCO 完整性评估
busco -i primary_assembly.fa \
    -l mammalia_odb10 \   # 哺乳动物保守基因集
    -o busco_output \
    -m genome \
    -c 16

# 3. Merqury 精确度评估（使用 k-mer）
# 先用 Illumina 短读建 k-mer 数据库
meryl count k=21 threads=16 \
    illumina_R1.fastq.gz illumina_R2.fastq.gz \
    output read-db.meryl

# 评估组装质量
merqury read-db.meryl \
    primary_assembly.fa \
    merqury_output

# 输出包含 QV 值（碱基准确率）和 completeness（k-mer 覆盖率）
cat merqury_output.qv
# QV >50 表示错误率 <1/100,000