基因组从头测序(de novo)相关软件

2 人赞了该文章

基因组从头测序(de novo)

回顾一下: 在使用aspera从SRA数据库得到二代测序的原始数据之后,我们使用ncbi的sratoolkit做过数据格式转化,然后使用fastqc做了质量分析控制, 之后使用trimmomatic对illumina的接头进行了剪切,得到了更加有效的数据。 为了便于分析,我们可以使用seqtk sample从原始数据中提取部分数据,作为测试使用。


下面理解一下,de novo测序中所需要的软件。 在得到大量的reads之后,我们需要对基因组进行从头测序,使用软件SPAdes,它是一个python脚本。 下载地址:

wget http://cab.spbu.ru/files/release3.12.0/SPAdes-3.12.0-Linux.tar.gz

具体看一下它的参数:

它的使用格式是:

Usage: ./spades.py [options] -o <output_dir>

它可以接收许多格式的原始reads数据,可以看帮助文档来获得具体信息,演示一下:

 ~/SPAdes-3.12.0-Linux/bin/spades.py   -1 ../test_7942raw_1.fq.gz  -2 ../test_7942raw_2.fq.gz -o OutFile~/SPAdes-3.12.0-Linux/bin/spades.py   -1 ../test_7942raw_1.fq.gz  -2 ../test_7942raw_2.fq.gz -o OutFile

 -1 -2参数分别表示的是PE测序的forward序列和reverse序列, -o是必须的,代表输出的结果

 当然,还有更多参数,比如老师ppt里面写的
 --careful表示的是 tries to reduce number of mismatches and short indels,减少错配和空位的次数;

然后脚本开始运行,过程如图所示: 运行完成之后,会得到contigs.fasta和scaffolds.fasta等一系列文件,其中:

contigs:拼接软件基于reads之间的overlap区,拼接获得的序列称为Contig(重叠群)。
Scaffold:基因组de novo测序,通过reads拼接获得Contigs后,往往还需要构建454 Paired-end库或Illumina Mate-pair库,以获得一定大小片段(如3Kb、6Kb、10Kb、20Kb)两端的序列。基于这些序列,可以确定一些Contig之间的顺序关系,这些先后顺序已知的Contigs组成Scaffold。

scaffold是最后的结果,fasta格式。 scaffold N50表示的是,将所有的scaffold长度相加,当达到scaffold总长度的0.5时的scaffold的长度就是scaffoldN50的值,它用来表示scaffold的好坏,contig N50的意思是一样的。

scaffold和染色体的区别

全基因组鸟枪测序的期望结果是把整个基因组整合到一个scaffold上面,然而这几乎是不可能的!

所以呢,即使是全基因组的某一个染色体也通常需要很多个scaffold来表示,但也有例外 Human chromosome 19就是一个单独的scaffold。

scaffold的长度不等,但是有一个Scaffold N50的指标,来判断该组装结果的好坏, 根据这个图其实就能很明确这些名词的意义了

简化来讲:如果基因组大小是100G,那么我可以选择几种打断方式. 如果打断片段是200bp,因为我们只能测序每个片段的首位和末尾的100bp,所以这个片段可以测通,那些只能根据各个read的overlap来拼接成contig,这个时候左右read之间没有gap。 如果打断片段是2Kb,而我们还是只能测到2KB长的片段的首位和末尾的100bp,我们依然是根据各个read的overlap来拼接成contig,但这时候左右read直接有着1800bp的gap,根据这左右read的关系,如果这两个read属于不同的contig,可以把这两个左右read各自属于的contig连接起来,这样变成了scaffold的,如此继续,直接没有更多的左右read来连接,最终的scaffold就是我们组装的scaffold 这样如果打断的片段更大,比如是15KB,那样它能把更远的read连接起来,这样可以把更多的read连接成大的scaffold,所以提高了scaffold N50. 但是打断的片段不可能无限大,所以scaffold是有极限的, 所以组装的测序结果总是不准确的,这些结果对于通常意义上的各种基因蛋白的研究已经足够了,(摘自新浪博客)

深度与覆盖度 测序深度:是指测序得到的总碱基数与待测基因组的大小的比值。假设一个基因大小是2M,测序深度为10×,那么获得的总数据量为20M 覆盖度:指测序获得的序列占整个基因组的比例。由于基因组中的高GC,重复序列等复杂结构的存在,最终获得的序列无法覆盖所有的区域,这部分没有获得的区域就称为Gap。如,覆盖度为98%,那么就还有2%的区域没有获得。 序列拼接评价软件quast 基本参数:

使用举例:

biocoder@iZ4owdgapta160Z:~/Seqs/test/OutFile$ ~/Biosofts/quast-5.0.0/quast.py  contigs.fasta -o quastReslut

最后就会形成一个质量报告,有html格式的,也有txt格式的,打开看一下:


在里面有具体的各个分段的contings的数量,以及NO50的大小。


评论

暂无评论