velvet序列的组装.docx

每日一生信-velvet序列的组装（终结版） (2014-03-02 10:48:06)Velvet用于拼接短序列组装安装1，下载地址：http:/www.ebi.ac.uk/zerbino/velvet/2，首先需要安装连个关联软件包zlib1g-dev和pdflatexsudo apt-get installzlib1g-dev（我的系统自带这个了）sudo apt-get install texlive（这个得下）3，解压下载的velvet软件包，tar xvzfvelvet_version.tgzcd velvet_versionmake colorCATEGORIES=2MAXKMERLENGTH=127BIGASSEMBLY=1LONGSEQUENCES=1OPENMP=1（参数怎么修改见下面）make clean 就是删除，然后再重新安装其他一些设置：make color设置不同的颜色？make MAXKMERLENGTH=127#默认的最大kmer为31，这个需要修改make CATEGORIES=57#处理的数据的管道数，不同的管道可以用来处理来自不同插入文库或者不同样品（这里面有一个问题：宏基因组算是多个样品吗？），而这个值越大，对内存消耗越大，而我的样品，默认的为2，我怕内存不够，暂且设为2吧make BIGASSEMBLY=1#read的ID保存在一个32字节的字符串中，如果你有大于22亿（2.2G）的reads，你就需要更多的内存，设置为1，是调用更多的内存。make LONGSEQUENCES=1#read的长度保存于16字节的字符串中，如果长度大于32kb,就意味着要更多的字节来存放字长度，所以可以直接设置为1，调用更多的内存，但是跟上面一样也会消耗很多的内存。make OPENMP=1#可以利用多个cpu，但是速度不会有指数的变化，因为只有一部分可以利用多线程make BUNDLEDZLIB=1#如果不能使用系统默认的zlib，也可以造常运行。也可以直接这样，一致性到位：make colorCATEGORIES=57MAXKMERLENGTH=127BIGASSEMBLY=1LONGSEQUENCES=1OPENMP=1环境变量：export PATH=$PATH:/sam/velvet/bin#注意我把安装好后的velveth,velvetg两个程序罗到了新建文件夹/bin里面，这样调用这两个程序的时候，就直接可以在终端的任何地方输入velveth,而不用cd 到程序目录，然后./velveth罗。也可以这样来改环境变量：sudo cp velvetg /usr/bin/sudo cp velveth /usr/bin/运行环境：64bit，gcc（GNU Compiler Collection，GNU编译器套装，是自由的类Unix及苹果计算机Mac OS X操作系统的标准编译器），12G内存主要工作原理图：1、Bruijn算法给定的数据总是存在着冗余信息，因此，算法的使用很重要。当前大多数是用的overlap-consensus-layout算法，它的每一个read都是单独的实体。而Bruijn算法，将分析基于observed words（或者k-mers）。给定k-mers，不管发现多少次，都是唯一的node。此外，他可以无差别的容纳long和short read的混合物。构建Bruijn图关键在于对所有的reads进行hash。Hash值通过设置参数给定，或者是默认的21。这对所有的序列进行成对比对时会减少时间，一旦所有的reads已经被hash，由k-mer node可以追溯他们的每一个路径。增值覆盖率和创建arcs就属于这种方法。2、错误的移除错误移除算法：Tour Bus，对没有错误连接的graph进行移除。这要确保基因组中低覆盖率的唯一的点没有被破坏。因此，这一步在画图之后进行。这里有两种类型的错误：tips和bubbles。Tips：起源于低质量的read的结尾部分没有任何的重叠。Bubbles：他产生于长read的中间或两个错误的read结尾部分重叠。首先去除tips，由长度和覆盖率来定义。然后去除bubbles，运用的是和Dijkstras算法相似的算法，最后低覆盖率的点会被淘汰。使用：合并两个fastq文件，如果是fasta,则需将下面的代码中fastq改为fasta/sam/velvet/contrib/shuffleSequences_fasta/shuffleSequences_fastq.pl output_forward_paired.fq output_reverse_paired.fq paired.fq第一步：建立哈希索引直接对fastq格式的原始文件进行处理，首先是用velveth命令建立hash表子集输入./velveth会出来使用帮助：./velveth directory hash_length -file_format-read_type-separate|-interleaved filename1 filename2 . . optionsdirectory: directory name for output fileshash_length: EITHER an odd integer（奇数）(if even, it will be decremented) = 127 (if above, will be reduced): OR: m,M,s where m and M are odd（奇数）integers (if not, they will be decremented) with m M = 127 (if above, will be reduced)and s is a step (even number). Velvet will then hash from k=m to k=M with a step of sfilename: path to sequence file or - for standard inputFile format options:-fasta-fastq-raw-fasta.gz-fastq.gz-raw.gz-sam -bam -fmtAuto(Note: -fmtAuto will detect fasta or fastq, and will try the following programs for decompression : gunzip, pbunzip2, bunzip2（默认的是fasta格式）File layout options for paired reads (only for fasta and fastq formats):-interleaved: File contains paired reads interleaved in the one file (default)-separate: Read 2 separate files for paired readsRead type options:-short-shortPaired-short2独立的插入大小的文库-shortPaired2同上-long -longPaired针对sanger和454-referenceOptions:-strand_specific: for strand specific transcriptome sequencing data (default: off)原来这个可以针对转录组数据啊-reuse_Sequences: reuse Sequences file (or link) already in directory (no need to provide original filenames in this case (default: off)-reuse_binary: reuse binary sequences file (or link) already in directory (no need to provide original filenames in this case (default: off)-noHash: simply prepare Sequences file, do not hash reads or prepare Roadmaps file (default: off)-create_binary: create binary CnyUnifiedSeq file (default: off)Synopsis:（大纲）- Short single end reads:velveth Assem 29 -short -fastq s_1_sequence.txt- Paired-end short reads (PE1和PE2是合并的remember to interleave paired reads):velvethAssem31-shortPaired-fastainterleaved.fna#Assem为生成文件所在文件夹，31为设置的kmer，一般为奇数，默认的是31，最开始默认最大的也为31，我已经改为127了，-shortPaired双末端模式，interleaved.fna？难道是是它指定的排列顺序？默认的输入文件是fasta格式- Paired-end short reads (PE1和PE2是分开的)velveth Assem 31 -shortPaired -fasta -separate left.fa right.fa- Two channels and some long reads:velveth Assem 43 -short -fastq unmapped.fna -longPaired -fasta SangerReads.fasta- Three channels:velveth Assem 35 -shortPaired -fasta pe_lib1.fasta -shortPaired2 pe_lib2.fasta -short3 se_lib1.faOutput:directory/Roadmapsdirectory/SequencesBoth files are picked up by graph, so please leave them there第二步：根据建立的索引来组装在终端输入./velvetg查看说明Usage:./velvetg directory optionsdirectory: working directory nameStandard options:-cov_cutoff : removal of low coverage nodes AFTER tour bus or allow the system to infer it(default: no removal)-ins_length : expected distance between two paired end reads (default: no read pairing)-read_trkg : tracking of short read positions in assembly (default: no tracking)追踪段read的位置（默认的是关）-min_contig_lgth : minimum contig length exported to contigs.fa file (default: hash length * 2)-amos_file : export assembly to AMOS file (default: no export)-exp_cov : expected coverage of unique regions or allow the system to infer it(default: no long or paired-end read resolution)-long_cov_cutoff : removal of nodes with low long-read coverage AFTER tour bus(default: no removal)Advanced options:-ins_length* : expected distance between two paired-end reads in the respective short-read dataset (default: no read pairing)-ins_length_long : expected distance between two long paired-end reads (default: no read pairing)-ins_length*_sd : est. standard deviation of respective dataset (default: 10% of corresponding length)replace * by nothing, 2 or _long as necessary-scaffolding : scaffolding of contigs used paired end information (default: on)-max_branch_length : maximum length in base pair of bubble (default: 100)-max_divergence : maximum divergence rate between two branches in a bubble (default: 0.2)-max_gap_count : maximum number of gaps allowed in the alignment of the two branches of a bubble (default: 3)-min_pair_count : minimum number of paired end connections to justify the scaffolding of two long contigs (default: 5)-max_coverage : removal of high coverage nodes AFTER tour bus (default: no removal)-coverage_mask : minimum coverage required for confident regions of contigs (default: 1)-long_mult_cutoff : minimum number of long reads required to merge contigs (default: 2)-unused_reads : export unused reads in UnusedReads.fa file (default: no)-alignments : export a summary of contig alignment to the reference sequences (default: no)-exportFiltered : export the long nodes which were eliminated by the coverage filters (default: no)-clean : remove all the intermediary files which are useless for recalculation (default : no)-very_clean : remove all the intermediary files (no recalculation possible) (default: no)-paired_exp_fraction : remove all the paired end connections which less than the specified fraction of the expected count (default: 0.1)-shortMatePaired* : for mate-pair libraries, indicate that the library might be contaminated with paired-end reads (default no)-conserveLong : preserve sequences with long reads in them (default no)Output:directory/contigs.fa: fasta file of contigs longer than twice hash lengthdirectory/stats.txt: stats file (tab-spaced) useful for determining appropriate coverage cutoffdirectory/LastGraph: special formatted file with all the information on the final graphdirectory/velvet_asm.afg: (if requested) AMOS compatible assembly filevelvetg Assem -cov_cutoff 30 -ins_length 350 -ins_length_sd 100 -exp_cov auto -min_contig_lgth 500#-ins_length为Insert size加上reads的长度，我前者为160，后者为91+91，所以后面的参数为350；默认的#-ins_length_sd 100为什么设置这两个参数就有点不明白了，Assem注意和第一步建立的文件夹的对应velvet拼接过程中主要的几个需要调试的参数主要有三个：Kmer值，exp_cov和cov_cutoff三个。本文分别设置不同的参数值，进行比较。Kmer值是基于De Bruijn算法的高通量读段拼接软件中最重要的一个参数，Kmer必须为奇数，在执行velveth时设置。在velvet的使用说明中，有一段关于Kmer设置的一般要求，如下：Ck= C(L-K+1)/LL表示读段长度，K表示Kmer值，C表示碱基的覆盖深度，Ck表示kmer的覆盖深度。L值已知，C值通过基因组规模和测序量来估计，当K值设定之后，可算得Ck值。根据经验，Ck值应当大于10，才能较好地完成拼接，而当Ck值大于20时，就浪费了测序深度。Kmer值越大Ck值就越小。当测序深度加大的时候，Kmer值也可以设的更大，有利于拼接。这个公式只能确定Kmer的取值范围，在此范围内选择采用哪个Kmer值，则需要进行Kmer参数值的调试。在执行velvetg时，可设定exp_cov和cov_cutoff两个选项。选项cov_cutoff用来过滤低覆盖深度的contigs，选项exp_cov为预期的覆盖深度。这两个参数可设定具体的值，默认值为“auto”。当设为auto时，软件在拼接时自动取值，如果数据的覆盖比较均匀，设为auto，拼接结果会比较好；exp_cov设为auto时，cov_cutoff也将强制为auto，auto表示覆盖深度的阈值为exp_cov的一半。如何使用多个Kmer值进行拼接If you wish for Velvet to e?ciently test various k-mer lengths without doi