(全文5058字)【推荐】用Smudgeplot评估物种倍性后,使用水母+GenomeScope1.0组合进行二倍体物种的基因组调查,使用CombineKMC+GenomeScope2.0用于多倍体物种的基因组调查。一、基因组调查1.1.GenomesurveyGenomesurvey是指对基因组特征的评价,一般是指通过k-mer对二代测序数据进行分析,得到基因组大小、杂合度、重复序列比、GC含量等基因组信息。Genomesurvey介绍在没有参考基因组的情况下,利用低、中深度短片段测序获得的reads进行基因组调查(genomesurvey),初步评估基因组特征,包括基因组大小(genomesize)、杂合度(heterozygosity)(杂合性)、重复序列比、GC含量等,从而为后续的基因组测序、组装和结构注释方案提供参考。基因组调查的目的基因组调查的主要目的是获取两个方面的信息,一是基因组的大小,二是基因组的复杂程度。Genomesize(基因组大小)因为测序成本是根据测序量计算的,基因组越大,测序成本越高。基因组复杂性基因组越复杂(杂合度越高,重复序列比例越高),测序和组装就越困难。基因组普查实践为了准确估计基因组信息,建议测序深度为50X(即估计基因组大小的50倍),最小不应低于25X。需要使用同一个体进行基因组调查和后续的基因组denovo测序,因为有些物种个体的基因组特征差异较大,可以利用同一个体调查的Illumina数据进行组装。估计的基因组大小可以通过现有研究粗略判断,包括流式细胞术研究,相关物种研究,也建议到C值数据库网站查询植物。基因组调查通常使用k-mer分析进行。1.2.基因组复杂度判断基因组复杂度方案的标准包括:基因组大小、倍性、杂合度、重复序列比、GC含量等。一般来说,基因组越大,重复序列所占比例越高;GC含量异常低或高,重复序列比例高;多倍体基因组的杂合度高于二倍体。判断基因组的复杂程度可以参考以下经验标准:简单基因组:单倍体;或纯合二倍体;或杂合度小于0.5%,重复序列小于50%,GC含量为35%-65%的二倍体。复杂基因组:杂合度在0.5%~1.2%之间的二倍体,或重复序列高于50%,或GC含量异常(<35%或>65%),或多倍体。高度复杂的基因组:杂合度>1.2%;或重复序列>65%。2.k-mer分析k-mer分析可用于生物信息学的许多方面。本博客中的k-mer分析特指用于基因组调查的k-mer分析方法。2.1.K-mer相关概念Monomericunit(mer)单体单位,单位为nt或bp。在双链核酸中通常用作一个单位,100merDNA相当于每条链中有100nt,因此整条链为100bp。k-mer的概念在生物信息学中,k-mer是指序列中包含的长度为k的子串。一段长度为L的核酸序列,以一个碱基为台阶滑动,总共可以产生(L-K+1)个k-mers。此外,该核酸的反向互补序列也可用于再生k-mer。图1.k-mer示例。图片来源:https://cloud.tencent.com/dev...2.2.k-mer分析步骤通过切割二代测序的reads对k-mers统计k-mers的总数和每个k-mer的频率绘制k-mer的频率分布图来确定预期的深度k-mer根据k-mer频率分布的主峰peak(即主峰对应的k-mer频率)。根据k-mers的预期深度和k-mers总数估计基因组大小。从低频k-mers估计数据错误率,并校正基因组大小的估计。从k-mer的其他峰估计K-mer杂合性和重复序列比率。2.3.k-mer原理应用k-mer分析的前提是测序reads在基因组上是随机分布的。先定义几个变量,方便解释原理:Genomesize:Greadreadlength:TotalnumberofLreads:nk-merlength:K2.3.1。Basedepthdistribution单个read序列覆盖某个碱基的概率:$L/G$。由于$L/G$较小,n较大,因此各碱基的覆盖深度服从泊松分布。那么每个碱基的覆盖深度的期望是:$d=(L/G)*n$。2.3.2.k-mer深度分布一个大小为G的基因组可以产生大约G个k-mer类型。假设一个基因组产生的K是唯一的,那么从一个大小为G的基因组中,可以得到$(G-K+1)$个不同的k-mers。一般来说,基因组大小G以百Mb或Gb为单位,远大于K和1,所以K和1可以忽略不计,约等于G种k-mers。单个读取序列完全覆盖某个k-mer的概率:$(L-K+1)/G$。k-mer类型总数约为G。单次read可产生的k-mer类型数量为$(L-K+1)$。假设reads随机分布在基因组上,单个read序列完全覆盖某个k-mer的概率为:$(L-K+1)/G$又因为$(L-K+1)/G$很小,n很大,每个k-mer的覆盖深度服从泊松分布。每个k-mer的预期覆盖深度为:$D=((L-K+1)/G)*n$。由此可得,基因组大小为:$G=(L-K+1)*n/D$。2.3.3.通过计算k-mers总数N和k-mers的预期深度D来计算基因组大小。对测序读数进行k-mer分割以获得k-mers的总数N。计算所有分裂的k-mers并绘制频率分布图。理想情况下(不考虑测序错误、序列重复、杂合序列),k-mer的频率分布服从泊松分布,频率分布峰值处的k-mer频率可以作为k的期望深度D-聚体。从下面的公式可以看出,基因组大小G=N/D。k-mers的总数$N=(L-K+1)*n$k-mers的期望深度$D=((L-K+1)/G)*n$只能通过k来计算聚体分析基因组大小G可以根据k聚体的总数N和k聚体的预期深度D计算得出。2.3.4.基因组调查在不考虑测序错误、序列重复和杂合序列的情况下,k-mer的深度分布服从泊松分布。但实际情况是三者都存在,因此需要计算错误率、重复序列比例和杂合度,并根据计算结果修正对基因组大小的估计。错误率测序错误:一般认为低频k-mers(K=1,2...)是由测序错误引起的。去除低频k-mers并计算错误率以纠正基因组大小估计。一般将拐点前的低频k-mer作为误差剔除。重复序列比例基因组中重复序列的存在会导致对应的k-mer频率高于预期的k-mer深度D。在k-mer频率分布图中,重复序列对应的k-mer会出现频率大于主峰的重复峰。根据主峰右侧的重复峰,可以估计基因组的重复序列比。通过模型计算一个阈值。例如,阈值是主峰的1.6倍。理论上,单拷贝序列(非重复序列)出现在1.6倍主峰右侧的概率很低,1.6倍主峰之后的可以视为重复k-mer。统计重复k-mers的总数m,这样可以计算:重复序列的大小R=重复k-mers总数m/k-mers的期望深度D。R/G是重复序列与基因组的比率。杂合子基因组中存在的杂合子序列会对应两种类型的k-mers,并使得这些k-mers出现的频率低于预期的k-mers深度D。假设两个k-mer的数目大致相等,则k-mer频率分布图中杂合序列对应的k-mer会出现频率约为主峰一半的杂合峰。可以根据主峰左侧的杂合性峰来估计基因组的杂合性。假设基因组中的杂合率为h(每个碱基为杂合点的概率),则一个k-mer为纯合k-mer的概率为$P1=(1-h)^k$,则是杂合子的概率$P2=1-P1=1-(1-h)^k$。假设生成的k-mer种类总数为M,其中k-mer种类属于单拷贝序列(非重复序列)。在二倍体中,单拷贝序列(非重复序列)的同源区域将产生U×P2的杂合k-mer。这个数可以通过计算得到:只要统计非重复k-mers的总数为M,M-U就是产生的杂合k-mers数。那么$M-U=U×P2=U×(1-(1-h)^k)$。通过计算,可以通过M和U的值得到基因组的杂合率h。在实际应用过程中,在对基因组的错误率、杂合度和重复序列比进行估计后,对基因组的估计重新校正大小以获得基因组调查的结果。图2.k-mer分析(软件GenomeScope)结果示例2.4。k-mer应用很多分析使用k-mer处理方法,测序得到的reads截取k-mer后用于分析。如评估基因组特征、组装基因组、物种样本污染评估等评估基因组特征(genomesurvey)包括评估基因组大小(size)、杂合度、重复序列比等基因组组装基因组中使用k-mer的目的assembly是去除低频k-mers,提高组装结果的准确性。评估基因组大小(size)在不考虑测序误差和序列重复性的情况下,k-mer的深度分布服从泊松分布,深度分布曲线的峰值可以作为预期的测序深度。k-mer评价基因组大小的公式在k-mer原理部分有解释:$$基因组大小G=总k-mer数N/k-mer预期深度D$$根据k-mer评价基因组杂合度频率分布图基因组的杂合度峰可以估计基因组的杂合度。估计基因组中重复序列的比例根据k-mer频率分布图的重复峰,可以估计基因组中重复序列的比例。评估物种样本污染根据k-mer频率分布图,可以评估测序样本的污染情况。如果k-mer频率分布图中有两个明显的峰,但两个峰的横坐标不是双重关系,则可能是DNA污染造成的。因为一个物种杂合峰的频率大约是主峰的一半,重复峰的频率大约是主峰的两倍。如果DNA中有两个物种的样本(即有污染),预计会有两个对应的峰,而两个峰的横坐标不是两倍的关系,因此可以初步确定测序样本被污染。至于是什么物种造成的污染,大家可以通过爆破nr库来做一个简单的判断。除了通过k-mer频率分布评价DNA样本的污染程度外,还可以通过GC含量分布图判断图中是否存在多个密度集中的簇。2.5.k-mer的优点和局限性提高了准确率。二代测序准确率达到99.9%。10,000)还是会对分析产生很大的影响。由于测序错误的随机性,在cuttingreads产生的k-mers中,大部分测序错误产生的k-mers是被测序物种中不存在的k-mers,因此只出现一次(或很少出现几次)times),如果去除这些低频k-mers,就更容易去除测序错误,从而使分析(genomesurvey,assembledgenome)结果更加可靠。不适用于过于复杂的基因组k-mer分析,适用于分析具有较大比例独特主峰区域的简单基因组。当基因组的杂合度很高或重复序列的比例很大时,其作用可能使k-mer分析无法正确估计基因组大小。2.6.k-mer的选择2.6.1k-mer大小的选择K应足够大,使k-mer能够定位到基因组中的唯一位置。过大的k-mer会降低以低频k-mer为代表的错误碱基被剔除的概率(增加错误率),同时也会降低k-mer深度(使得k-mer频率分布的峰值不明显)很明显),并且大的k-mer会增加计算资源的使用。基因组调查一般选择17,21比较常见。2.6.2选择k-mer的考虑k-mer只能是奇数?设置k-mer为奇数的目的是为了防止通过k-mer组装时正反链混淆。偶数个k-mer的反向互补序列往往与自身相同,这样组装k-mer时正负链的组装就会混淆,而奇数个则不存在这个问题编号的k-mer。k-mer的长度代表可能的k-mer种类的数量(4的K次方),k-mer片段越长,mapping的种类特异性越强。基因组越大,所需的k-mer越长。当基因组中有很多重复序列时,可以使用更大的k-mer来跨越高重复区域,从而获得更准确、更完整的基因组草图;由于reads上碱基错误率的存在,选择更长的k-mer会带来更高的错误率,这也可以通过增加测序深度来补偿。如果用于组装基因组,为了获得更完整的基因组,应尽可能使用更长的k-mer进行组装。3.k-mer分析软件3.1.k-mer分析软件介绍k-mer分析分为两个步骤:k-mer频率统计和基因组特征评估。此外,Smudgeplot还可以使用k-mer分析来评估物种的倍性。jellyfishjellyfish可以实现k-mer频率统计的第一步。特点是使用Hash表存储数据,可以多线程运行;它具有速度快和内存消耗低的特点。KMCKMC可以实现k-mer频率统计的第一步。GenomeScopeGenomeScope可以实现第二步,利用k-mer频率统计结果进行基因组特征评估。1.0版用于二倍体物种,2.0版用于多倍体物种。KAT(Thek-merAnalysisToolkit)KAT(Thek-merAnalysisToolkit)可以实现k-mer频率统计和基因组特征评估两个步骤。包含多个工具帮助用户利用k-mer对测序数据进行简单的分析,如装配完整性、测序错误、污染等。GCEGCE可以分别实现k-mer频率统计和基因组特征评估两步。KmerGenieKmerGenie可以同时实现k-mer频率统计和基因组特征评估两个步骤。最大的优势是可以实现多个预设k-mers下的自动分析。除了常规的k-mer频率统计外,它还可以根据不同的k-mers自动计算基因组大小,评估最佳基因组组装。组装k-mer值作为替代。3.2.使用软件的一些经验总结【推荐】用Smudgeplot评估物种倍性后,使用jellyfish+GenomeScope1.0组合进行二倍体物种基因组调查,使用KMC+GenomeScope2.0组合进行多倍体基因组调查物种。k-mer的长度通常为17/21。软件KmerGenie、GCE和jellyfish得到的频率分布表可以用于软件genomescope和GCE的第二步分析。由于老版本GCE第一步支持的最大k-mer频率为255,大于255的数据合并(新版本GCE好像更新了,目前没有这个限制);而jellyfish统计10000行,估计的结果会更准确。GenomeScope的高度重复序列的基因组统计数据的基因组大小将很小。建议将maxkmercoverage设置的大些,大于等于10000。有些软件(如GCE)还有一个参数需要注意和设置,无论是单倍体模式还是杂合子模式,你都可以分析这两种模式,看看有什么不同。实践经验发现,k-mer值设置的越高,估计的基因组大小就越大;另外,jellyfishhisto统计水母中的频率分布时,使用参数-h10000增加统计上限,在GenomeScope阶段,Maxkmercoverage设置的越大(即统计的kmer越多),估计的基因组大小会稍大。4.参考资料wiki_k-mer:https://en.wikipedia.org/wiki...genomesurvey:https://xuzhougeng.top/archiv...xuzhougeng的博客:https://www.jianshu.com/p/85d...k-mer与基因组组装:https://cloud.tencent.com/dev...k-mer分析与原理:https://www.bbsmax.com/A/lk5a...jellyfishpaper:https://academic.oup.com/bioi...jellyfishgithub:https://github.com/gmarcais/J...GenomeScope1.0github:https://github.com/schatzlab/...KATgithub:https://github.com/TGAC/KATGCEgithub:https://github.com/fanagislab...
