组装"预实验" | 物种倍性分析

文摘   2024-12-07 09:55   北京  


在自然演化过程中,多倍化和全基因组复制事件常常被被认为是物进化的推动力,它们为生物进化提供了原始的遗传材料,同时也导致了物种的倍性差异。倍性直接影响基因组的大小和复杂性,不同倍性的物种需要不同的组装策略,对于多倍体物种,可能需要采用特殊的算法和工具来处理重复序列和同源区域。因此我们在进行基因组评估之前,需要先做一下倍性分析。

01

倍性分析有哪些方法?

1.实验方法
压片法

通过促进细胞有丝分裂,使用秋水仙素使细胞分裂停止在有丝分裂中期,然后固定、酶解、染色、烘烤等方法将细胞染色体固定在载玻片上,使用核型分析显微镜观察分析染色体的数目以及形态;

荧光原位杂交
将荧光素标记的特异性探针与待测细胞中的目标DNA进行杂交,洗涤后在荧光显微镜下观察目标DNA的位置、大小、形态、数量等;
流式细胞:
细胞样品固定后加入DAPI荧光染料,通过流式细胞仪检测,与已知倍性的对照样品比较计算后可以获得DNA的倍性关系。
2. 测序数据分析
Smudgeplot:
从k-mer数据库中提取杂合k-mer对,然后训练这些杂合k-mer对。通过比较k-mer对覆盖度的总数CovA+CovB和相对覆盖度CovB/(CovA + CovB)统计杂合k-mers对的数量,解析基因组结构,从而估计物种的倍性;
GenomeScope:
主要基于k-mer进行基因组大小的评估。GenomeScope2.0添加了参数倍性Ploidy,可以评估多倍体的基因组特征。四倍体共有两种可能的拓扑结构,代表着同源四倍体和异源四倍体,每种拓扑包含三种杂合基因型和一种纯合基因型,共有五种基因型。(五倍体有五种可能的拓扑,六倍体有十六种)。根据结果中杂合基因型的分布模式可以区分异源四倍体和同源四倍体。

02


什么是k-mer?


k-mer是指包含在一段序列中的长度为k的子串

现有一段序列,长度是15个碱基,把k的长度设置成5,这意味着需要从这段序列上,每隔一个碱基取一个长度为5个碱基的序列片段。那么序列长度为15的序列可以取11个5mer。如图:

注意:

K需为奇数,防止正反链混淆;

针对长度为K的K-mer,对于 A,T,C,G四种碱基类型,一共能产生的K-mer种类数为4^k个

03

jellyfish+Smudgeplot


1.获取k-mer频数分布表(假设样本名为A)

ls /Data/A* | awk '{ print "gzip -dc " $0}' > A_generate.file
jellyfish count -t 20 -C -m 21 -s 3G -g A_generate.file -G 2 -o A.sif
# jellyfish count:k-mer计数# -t 指定线程数# -C 对DNA正负链都进行统计,表示考虑DNA正义与反义链。如果是双端测序reads,需要这个参数。# -m k-mer长度设置为21bp.# -s 3G 存储用的hash表大小为3G,这个参数识别单位M(Mbp)和G(Gbp)。如果设置不够大,会生成多个hash文件。最好设置的值大于总的独特的(distinct)k-mer数。如果基因组大小为G,每个reads有一个错误,总共有n条reads,则该值可以设置为[(G + n)/0.8]。# -g --generator=path 记录产生fast[aq]命令的文件# -G 同时运行的数目# -o 指定输出文件的名字,hash格式储存的k-mer频数文件

2.统计k-mer频率

jellyfish histo -v -t 20 -h 10000000 -o A.histo A.sif
# 结果最后是两列,第一列是x,表示的是出现的次数;第二列是y,表示出现x次数的kmer数目。# -v 显示详细信息# -t 线程数# -h x的最大值,默认是10000# -o 输出文件的名字。

3. 提取杂合k-mer对
# 估计覆盖阈值上(U)下(L)限U=$(smudgeplot.py cutoff A.histo U)L=$(smudgeplot.py cutoff A.histo L)
# 计算k-mer pairsjellyfish dump -c -L $L -U $U A.sif | smudgeplot.py hetkmers -o kmer_pairs# dump:将Jellyfish计数后的k-mer数据转储(dump)出来。# -c 参数表示输出覆盖度(coverage)信息# -L $L 和 -U $U 参数分别表示覆盖度的下限(lower bound)和上限(upper bound),确定的合理阈值。# smudgeplot.py hetkmers 寻找杂合k-mer pairs
4. 倍性评估
# 绘图smudgeplot.py plot kmer_pairs_coverages.tsv -o my_genome
# 结果图片:my_genome_smudgeplot.png# 具体数值:my_genome_summary_table.tsv

04

KMC+Smudgeplot


1.用KMC计算k-mer频率,生成k-mer直方图(假设样本名为A)

mkdir tmpls /Data/A* | awk '{ print "gzip -dc " $0}' > A_generate.file
kmc -k21 -t16 -m64 -ci1 -cs10000 A_generate.file kmcdb tmp #计算k-mer频率# -k21:k-mer长度设置为21# -t 指定线程数# -m64:内存64G,设置使用RAM的大致数量,范围1-1024。# -ci1 -cs10000:统计k-mer coverages覆盖度范围在[1-10000]。# A_generate.file:保存了输入文件列表的文件名为A_generate.file# kmcdb:KMC数据库的输出文件名前缀# tmp:临时目录
kmc_tools transform kmcdb histogram A.histo -cx10000 #生成k-mer频数直方表A.histo和k-mer直方图# -cx10000:储存在直方图文件中counter的最大值。
2. 提取阈值范围的k-mers,计算k-mer pairs
# 估计覆盖阈值上(U)下(L)限U=$(smudgeplot.py cutoff A.histo U)L=$(smudgeplot.py cutoff A.histo L)
# 用kmc_tools提取k-merskmc_tools transform kmcdb -ci"$L" -cx"$U" reduce kmcdb_L"$L"_U"$U"# 用KMC的smudge_pairs计算k-mer pairssmudge_pairs kmcdb_L"$L"_U"$U" kmcdb_L"$L"_U"$U"_coverages.tsv kmcdb_L"$L"_U"$U"_pairs.tsv > kmcdb_L"$L"_U"$U"_familysizes.tsv
3. 倍性评估
# 绘图smudgeplot.py plot kmcdb_L"$L"_U"$U"_coverages.tsv -o my_genome
# 结果图片:my_genome_smudgeplot.png# 具体数值:my_genome_summary_table.tsv



05


结果解读


在当前文件夹中会存在2个图片,三个文件。其中最重要的结果图片是 my_genome_smudgeplot.png,横轴是杂合k-mer对的深度占总深度的比例,纵轴是所有k-mer对的总深度,二者的交点即代表了不同的杂合结构。交点的亮度代表了落入其中的k-mer数。下图表示AAB可能性最高,所以推测是三倍体。


参考资料:
https://www.nature.com/articles/s41467-020-14998-3
https://github.com/KamilSJaron/smudgeplot?tab=readme-ov-file
https://zhuanlan.zhihu.com/p/529851006


扫描下方二维码了解基因组组装课程:



生信课堂
生信笔记
 最新文章