Velvet 1.2.06 发布

这一版本版本最大的改进就是可以不需要再将PE的两端数据合并到一起了,使用-separate参数可以使用两端分离的PE数据,其他没啥变化(难道这还不够你尖叫的吗?!),使用方法如下:

velveth Assem 31 -shortPaired -fasta -separate left.fa right.fa
Posted in 生物信息学. Tags: , . 4 Comments »7,091 views

nucleotide blast(blastn)中空位罚分gap costs的计算方法

为什么要把这个看上去很简单的问题写下来呢,因为里面有些地方很容易引起误解,写下来给大家一个参考。

首先,blastn中默认的linear的罚分是如何计算的?所谓linear线性罚分,就是gap中每个碱基都扣相同的分数,那么到底是扣几分呢?我没有找到官方的说法,但是通过我的试验,我得出的结论是每个gap中的每个碱基扣2.5分,得到一个总分后还要四舍五入一下,于是1个碱基空位扣3分,两个扣5分,三个扣8分,以此类推。

其次,选择其他gap costs选项时,比如“existence 5 extension 2”,存在一个空位扣5分,每个延伸扣2分。一开始我想当然地认为是空位扣分后每个延伸是从第二个空位碱基开始算的,那么单个碱基的空位应该是扣5分,但是实际情况却是空位延伸罚分是从第一个碱基开始算的,单个碱基的空位就要扣5+2分,两个碱基的空位是5+2*2,n个碱基的空位就要扣5+n*2了。

最后,没有了~

Posted in 生物信息学. Tags: . No Comments »3,901 views

序列拼接软件velvet 1.1.01发布

2011-3-29日发布了velvet最新版本1.1.01,在这一版本中终于支持多线程运行了!!
速度可以提高4-10倍,还有一些另外的修改,附上作者的原话:

Dear Velvet users,

it is my pleasure to announce the release of Velvet 1.1.01!

Sylvain Forêt from ANU did a excellent job multithreading Velvet, 
meaning that it can use multiple CPUs on the same computer, thus freeing 
up the computer that much faster. 4 to 10 fold accelerations have been 
reported. He also reduced the memory footprint, and improved the use of 
mate-pairs, filtering out paired-end contamination.

Also in this new version Yasubumi Sakakibara and colleagues shared their 
MetaVelvet code, which is designed to analyze metagenomic data.

Regards,

Daniel
Posted in 生物信息学. Tags: , . 10 Comments »8,730 views

用Amos软件包里面的minimus2合并454和Illumina/Solexa拼接得到的contig

用罗氏454测序得到的序列用newbler拼接的效果最好,而用短序列拼接软件velvet拼接效果很差,所以不能将454的原始reads和Illumina产生的reads合到一起后用velvet进行拼接。在用newbler和velvet分别拼接454和Illumina的reads得到contigs之后,我们就需要将两者的contig再合并起来,得到更好的拼接结果。这里就介绍一个简单易用的软件minimus2。

minimus2是amos拼接软件包里面的一个组件,它的功能就是将两组contig进行合并,延伸contig的长度,减少contig的数量。Amos是A Modular, Open-Source whole genome assembler的缩写,致力于打造成一个拼接软件的基础软件系统。minimus2用的是基于nucmer overlap检测的算法,速度上比Smith-Waterman hash-overlap的算法要快,下面就介绍一下用法。

首先当然是下载amos软件包进行安装,下载地址为:http://sourceforge.net/projects/amos/files/

安装啥的就不说了,根据说明来就行。安装完成之后,minimus2软件位于amos安装文件夹下的bin里面。在运行minimus2之前首先要准备好文件,比如现在有s1.fa和s2.fa两组包含contig的文件,首先要知道里面包含的contig数目,针对fasta格式,用

grep -c "^>" s1.fa s2.fa  命令得到,比如分别为100和200个contig。

然后用cat命令合并到一个文件:

cat s1.fa s2.fa >s1_s2.fa

再用amos里面的另一个软件toAmos转换成Amos格式,这个软件也位于bin文件夹下面

./toAmos -s s1_s2.fa -o s1_s2.afg
这里的-s是指输入的为fasta格式。

然后就可以运行minimus2了

minimus2的运行参数为:

 minimus2 prefix  \
   -D REFCOUNT=n  \  # Number of sequences is the first set
   -D OVERLAP=n   \  # Minimum overlap (Default 40bp)
   -D CONSERR=f   \  # Maximum consensus error (0..1) (Def 0.06)
   -D MINID=n     \  # Minimum overlap %id for align. (Def 94)
   -D MAXTRIM=n      # Maximum sequence trimming length (Def 20bp)
最简单的命令为:
./minimus2 s1_s2 -D REFCOUNT=100

这里只要告诉文件名(不要后缀)和作为参考序列的第一组contig的数目就可以了。会生成一堆以s1_s2开头的文件,其中s1_s2.fasta就是合并之后得到的contig文件。

大功告成!

应用于第二代测序技术的生物信息学工具[zz]

刚才闲逛发现了这篇被多家博客转载的文章,来自于SEQanswer,总结地比较全面,不过已经比较老了,随便看看~

Integrated solutions
* CLCbio Genomics Workbench – de novo and reference assembly of Sanger, Roche FLX, Illumina, Helicos, and SOLiD data. Commercial next-gen-seq software that extends the CLCbio Main Workbench software. Includes SNP detection, CHiP-seq, browser and other features. Commercial. Windows, Mac OS X and Linux.
* Galaxy – Galaxy = interactive and reproducible genomics. A job webportal.
* Genomatix – Integrated Solutions for Next Generation Sequencing data analysis.
* JMP Genomics – Next gen visualization and statistics tool from SAS. They are working with NCGR to refine this tool and produce others.
* NextGENe – de novo and reference assembly of Illumina, SOLiD and Roche FLX data. Uses a novel Condensation Assembly Tool approach where reads are joined via “anchors” into mini-contigs before assembly. Includes SNP detection, CHiP-seq, browser and other features. Commercial. Win or MacOS.
* SeqMan Genome Analyser – Software for Next Generation sequence assembly of Illumina, Roche FLX and Sanger data integrating with Lasergene Sequence Analysis software for additional analysis and visualization capabilities. Can use a hybrid templated/de novo approach. Commercial. Win or Mac OS X.
* SHORE – SHORE, for Short Read, is a mapping and analysis pipeline for short DNA sequences produced on a Illumina Genome Analyzer. A suite created by the 1001 Genomes project. Source for POSIX.
* SlimSearch – Fledgling commercial product. Continue reading “应用于第二代测序技术的生物信息学工具[zz]” »

Posted in 生物信息学. Tags: , , . 7 Comments »83,222 views

得分矩阵PAM与BLOSUM的比较与区别

  对于蛋白质序列,计分矩阵主要用于记录在做序列比对时两个相对应的残基的相似度,一旦这个矩阵定义好了以后,比对程式就可以利用这个矩阵,尽量将相似的残基排在一起,以达到最好的比对。
  得分矩阵主要有两种,第一种就是PAM(Point Accepted Multation),另一种就是BLOSUM。
1、PAM矩阵(Point Accepted Mutation)
   基于进化的点突变模型,如果两种氨基酸替换频繁,说明自然界接受这种替换,那么这对氨基酸替换得分就高。一个PAM就是一个进化的变异单位, 即1%的氨基酸改变,但这并不意味100次PAM后,每个氨基酸都发生变化,因为其中一些位置可能会经过多次突变,甚至可能会变回到原来的氨基酸。
PAM矩阵的制作步骤:
  构建序列相似(大于85%)的比对
  计算氨基酸 j 的相对突变率mj(j被其它氨基酸替换的次数)
  针对每个氨基酸对 i 和 j , 计算 j 被 i 替换次数
  替换次数除以相对突变率(mj)
  利用每个氨基酸出现的频度对j 进行标准化
  取常用对数,得到PAM-1(i, j)
  将PAM-1自乘N次,可以得到PAM-N。

  这种矩阵的缺点是一旦PAM1的矩阵有效地误 差,那么自乘250后得到的PAM250矩阵的误差就会变得很大。如,PAM120矩阵用于比较相距120个PAM单位的序列。
一个PAM-N矩阵元素(i,j)的值:
反应两个相距N个PAM单位的序列中第i种氨基酸替换第j种氨基酸的频率。
针对不同的进化距离采用PAM 矩阵
序列相似度 = 40%       50%       60% 
                   |              |              |
打分矩阵 = PAM120 PAM80 PAM60
PAM250 → 14% – 27% 
2、BLOSUM 矩阵
此矩阵与PAM矩阵的不同之处在于:
(1)用于产生矩阵的蛋白质家族及多肽链数目,BLOSUM比PAM大约多20倍。  
(2)PAM:家族内成员相比,然后把所有家族中对某种氨基酸的比较结果加和在一起,产生“取代”数据(PAM-1 );PAM-1自乘n次,得PAM-n。
BLOSUM:首先寻找氨基酸模式,即有意义的一段氨基酸片断(如一个结构域及其相邻的两小段氨基酸序列) 
,分别比较相同的氨基酸模式之间氨基酸的保守性(某种氨基酸对另一种氨基酸的取代数据),然后,以所有 60%保守性的氨基酸模式之间的比较数据为根据,产生BLOSUM60;以所有80%保守性的氨基酸模式之间的比 较数据为根据,产生BLOSUM80。 
(3)PAM-n中,n 越小,表示氨基酸变异的可能性越小;相似的序列之间比较应该选用n值小的矩阵,不太相似 的序列之间比较应该选用n值大的矩阵。PAM-250用于约20%相同序列之间的比较。BLOSUM-n中,n越小,表示氨基酸相似的可能性越小;相似的序列之间比较应该选用 n 值大的矩阵,不太相似的序列之间比较应该选 用n值小的矩阵。BLOSUM-62用来比较62%相似度的序列,BLOSUM-80用来比较80%左右的序列。

Posted in 生物信息学. 2 Comments »8,744 views

序列拼接软件velvet 1.0.14发布

这个版本主要是修正了几个bug,加了一些小功能:

1、重新加入了每个标准输出信息的时间戳;

2、在帮助增添了对于k-mer长度的说明;

3、捕获写文件时发生的错误;

4、修正文件访问权限;

5、消除编译时的警告和错误;

6、修正一个死循环的bug;

7、修正一个内存溢出问题;

8、加入-clean和-very_clean参数,-clean选项会清除Graph和Graph2之外的文件,-very_clean选项会清除所有文件。

Posted in 生物信息学. Tags: , , . 6 Comments »6,794 views

Illumina/Solexa测序Reads质量过滤软件包

刚刚在群组中看到SolexaQA软件包的发布,发表在BMC Bioinformatics上的文章,三大组件用途:

* SolexaQA -- 从Fastq文件统计reads质量并可视化显示。
* DynamicTrim -- 根据用户设定的测序质量阈值过滤低质量碱基。
* LengthSort -- 根据用户设定的Reads长度过滤质量过滤后低于此长度的Reads。 
软件下载地址:http://solexaqa.sourceforge.net/
文章地址: http://www.biomedcentral.com/1471-2105/11/485
Posted in 生物信息学. Tags: , . 2 Comments »6,488 views

Protein Blast中Blastp/PSI-BLAST/PHI-BLAST区别与选择

说完blastn,接着说blastp~blsatp中也有三个不同的算法可以选择,如下:

blastp (protein-protein BLAST)就是简单地进行蛋白与蛋白的比对,寻找蛋白质相似序列;

PSI-BLAST (Position-Specific Iterated BLAST)叫做位点特异性迭代比对,它在蛋白质数据库中循环搜索查询蛋白质,所有前一次被psi-blast发现的统计显著蛋白质序列将整合成新记分矩阵,通过多次迭代比对,直到不再发现统计显著的新蛋白质;

PHI-BLAST (Pattern Hit Initiated BLAST)可以在搜索的时候限定蛋白质的模式(pattern),只给出包含此模式的比对结果。

Posted in 生物信息学. Tags: . 1 Comment »10,198 views

nucleotide blast中MegaBlast/Discontiguous MegaBlast/BlastN的区别与选择

今天在准备学院培训班的试讲,复习了这部分内容,于是写下来吧~

从blastn页面上的简单帮助可以看到Highly similar sequences (megablast)多用于比较相似性比较高(相似性在95%以上)的序列,速度快;More dissimilar sequences (discontiguous megablast)用于相似性稍低于megablast的比对,但是灵敏度和精确度更高,多用于不同物种间的同源比对;而Somewhat similar sequences (blastn)用于比对相似性较差的序列,可以比对最短7个碱基的长度,所以比对精确度最高,比对结果最多,速度最慢。

所以,在选择的时候根据你提交的序列和搜索的目的进行选择,如果是想看这段序列在数据库当中是否有收录,可以用megablast,如果想用其他物种的基因注释信息来注释一个未注释物种的序列,可以选择discontiguous megablast,如果想得到更多更全面的结果,可以选择blastn。

更详细的解释可以查看:http://www.ncbi.nlm.nih.gov/blast/producttable.shtml#tab31

Posted in 生物信息学. Tags: , . 1 Comment »6,051 views