一款用于评估基因组质量的新方法
一般用于评估的方法
- 二代reads 比对率； 偏向于重复序列...
- BUSCO；仅仅能检测单拷贝
最近，GEnome Biology 发表一新工具'Merqury: reference-free quality, completeness, and phasing assessment for genome assemblies',
githup: https://github.com/marbl/merqury
一款新的软件Merqury, 基于Kmer的方法进行鉴定，据文章说是目前最牛叉的，现就consensus quality value (QV) 评估进行说明，
其中图b的红色区域，用于鉴定QV，Q30等于99.9%正确率，Q40为99.99%正确率
1 软件安装
conda create -n merqury -c bioconda merqury
2 简单操作
ln -s $MERQURY/merqury.sh # Link merqury
./merqury.sh <read-db.meryl> [<mat.meryl> <pat.meryl>] <asm1.fasta> [asm2.fasta] <out>
Usage: merqury.sh <read-db.meryl> [<mat.meryl> <pat.meryl>] <asm1.fasta> [asm2.fasta] <out>
<read-db.meryl> : k-mer counts of the read set
<mat.meryl> : k-mer counts of the maternal haplotype (ex. mat.only.meryl or mat.hapmer.meryl)
<pat.meryl> : k-mer counts of the paternal haplotype (ex. pat.only.meryl or pat.hapmer.meryl)
<asm1.fasta> : Assembly fasta file (ex. pri.fasta, hap1.fasta or maternal.fasta)
[asm2.fasta] : Additional fasta file (ex. alt.fasta, hap2.fasta or paternal.fasta)
*asm1.meryl and asm2.meryl will be generated. Avoid using the same names as the hap-mer dbs
<out> : Output prefix
2.1 构建kmer的db
- （1）确定合适的kmer 大小
best_k.sh <genome_size> [tolerable_collision_rate=0.001]
# 如下
best_k.sh 3100000000
Kmer 为21
- （2）构建db （illumina 为例子）
for each read$i.fastq.gz
do
# 1. Build meryl dbs
meryl k=$k count output read$i.meryl read$i.fastq.gz
done
## 合并两个reads 的db
meryl union-sum output $genome.meryl read*.meryl
如果没有hap-mers就使用如下命令即可
$MERQURY/merqury.sh read-db.meryl asm1.fasta out_prefix
其中read-db.meryl 是根据二代reads建的db； asm1.fastat 组装基因组序列；ou t_prefix，输出结果
2.2 结果
查看QV结果
head */*.qv
==> triocanu_clr/triocanu_.qv <==
col 682359 123511626 35.1183 0.000307729
cvi 506306 122349407 36.3761 0.00023035
Both 1188665 245861033 35.6991 0.00026921
第一列：组装基因组
第二列：基因组中特有的kmer
第三列：基因组和reads中均存在的kmer
第四列：QV
第五列：错误率
查看完整度
head completeness.stats
col all 104975080 125303808 83.7764
cvi all 104809523 125303808 83.6443
both all 123134729 125303808 98.2689
计算公式：
k-mer completeness = found solid k-mers in an assembly / solid k-mers in a read set
第一列：组装基因组
第二列：all- reads set
第三列：基因组的solid kmer
第四列：reads 中的总solid kmer
第五列：完整度(%)
详细用法请查看https://github.com/marbl/merqury， 若果使用该软件，不要忘记引用：
Rhie, A., Walenz, B.P., Koren, S. et al. Merqury: reference-free quality, completeness, and phasing assessment for genome assemblies. Genome Biol 21, 245 (2020). https://doi.org/10.1186/s13059-020-02134-9
============================================
2024年1.5
meryl的一些参数学习
- 用于输出
- count-forward； 正向序列的kmer
- count-reverse：反向序列的kmer
- greater-than： 输出>指定的kmer
meryl print greater-than 50000 merylDB > repetitive_k15.txt
- union-sum: 2个数据集取并集，出现次数之和
