一般从连续性、完整性、准确性三个方面进行评估。
N50等统计
assembly-stats genome.fasta > genome.N50_stats
CG含量等统计
seqkit fx2tab -g -H -l -n -i genome.fasta > genome.all_stats二代短 reads 比对
# 构建基因组index
bwa index -p genome genome.fasta
# 比对并排序
bwa mem -t 8 genome ecoli_R1.fastq.gz ecoli_R2.fastq.gz | samtools sort - -o aln_short.bam
# 统计比对率
samtools flagstats aln_short.bam > genome.short.flagstas
# 计算每条序列覆盖度
samtools coverage aln_short.bam > genome.short.coverage
# 汇总获得全基因组覆盖度
awk 'NR > 1 {sum+=$3;rd+=$4; cov += $5; dp+=$7*$5 }END{print \
"Total: "sum"\nCovbase: "cov"\ncoverage: "cov/sum"\nmeandepth: "dp/sum}' genome.short.coverage > genome.short.coverage_all
三代长 reads 比对
# minimap2比对
minimap2 -t 4 -ax map-ont --secondary=no genome.fasta \
reads.fastq.gz | samtools sort - -o aln_ont.bam
# 比对率统计
samtools flagstats aln_ont.bam > genome.ont.flagstas
# 覆盖度统计
samtools coverage aln_ont.bam > genome.ont.coverage
awk 'NR > 1 {sum+=$3; rd +=$4; cov += $5; dp+=$7*$5 }END{print \
"Total: "sum"\nCovbase: "cov"\ncoverage: "cov/sum"\nmeandepth: "dp/sum}'LAI
LAI 指标基于基因组的 LTR 进行完整性评估。总LTR-RT含量小于5%，完整LTR-RT含量小于0.1%时，LAI的估计是不准确。
EDTA中会自动计算
awk '{print $1}' Ath.fasta > genome.fa
# 进行index构建
gt suffixerator \
-db genome.fa \
-indexname genome.fa \
-tis -suf -lcp -des -ssp -sds -dna
# 运行ltrharvest进行LTR预测
ltrharvest \
-index genome.fa \
-minlenltr 100 -maxlenltr 7000 -mindistltr 1000 -maxdistltr 15000 \
-mintsd 4 -maxtsd 20 -motif TGCA -motifmis 1 \
-similar 85 -vic 60 -seed 20 -seqids yes > genome.fa.harvest.scn
# 运行LTR_FINDER 进行LTR预测
LTR_FINDER_parallel_wy \
-seq genome.fa -threads 10 -harvest_out -size 5000000 -time 500 \
-finder_para "-w 2 -C -D 15000 -d 1000 -L 7000 -l 100 -p 20 -M 0.85 "
# 合并结果
cat genome.fa.harvest.scn genome.fa.finder.combine.scn >genome.fa.rawLTR.scn
# 使用LTR_retriever进行整合，并计算LAI
LTR_retriever \
-genome genome.fa \
-inharvest genome.fa.rawLTR.scn \ -threads 10| Category | LAI | Examples |
|---|---|---|
| Draft | 0 ≤ LAI <10 | Apple (v1.0), Cacao (v1.0) |
| Reference | 10 ≤ LAI <20 | Arabidopsis (TAIR10), Grape (12X) |
| Gold | 20 ≤ LAI | Rice (MSUv7), Maize (B73 v4) |
merqury
merqury 软件基于测序数据 kmer 和组装结果 kmer 的一致性进行基因组完整性和 QV 计算，从而对完整性进行评估。
# 获取K值
best_k.sh
# kmer统计
meryl \
k=16 memory=24G threads=10 \
count R1.fastq.gz output R1_out.meryl s
meryl \
k=16 memory=24G threads=10 \
count R2.fastq.gz output R2_out.meryl
# kmer合并
meryl k=16 union-sum output all.meryl R1_out.meryl R2_out.meryl
# 运行merquery
merqury.sh all.meryl genome.fasta merqury_out结果文件：merqury_out.qv为基因组 qv 结果，merqury_out.completeness.stats为基因组完整性统计
BUSCO 评估
# 列出数据库(需要联网)
busco --list-datasets
# 自动下载数据库并运行busco
busco -m geno -i genome.fa -o bacteria_odb10 -l bacteria_odb10 -c 10
# 使用本地busco数据库运行
# db下载地址:
# https://busco-data.ezlab.org/v5/data/lineages/ mkdir buscoDB
cd buscoDB
wget https://busco-data.ezlab.org/v5/data/lineages/\
bacteria_odb10.2020-03-06.tar.gz
tar -zxvf bacteria_odb10.2020-03-06.tar.gz
wget https://busco-data.ezlab.org/v5/data/lineages/\
gammaproteobacteria_odb10.2021-02-23.tar.gz
tar -zxvf gammaproteobacteria_odb10.2021-02-23.tar.gz
cd ../
# 使用本地busco数据库运行busco
busco -m geno \
-i genome.fa \
-o bacteria_offline \
-l ./buscoDB/bacteria_odb10 \ -c 10 \
--offline
busco -m geno \ -i genome.fa \
-o gammaproteobacteria_offline \
-l ./buscoDB/gammaproteobacteria_odb10 \
-c 10 \ --offline
# 统计、比较、画图
mkdir BUSCO_summaries
cp ./*/short_summary.*.txt BUSCO_summaries/
generate_plot.py \
-wd BUSCO_summaries/