RNA-seq差异表达基因分析之TopHat篇

TopHat是基于Bowtie的将RNA-Seq数据mapping到参考基因组上,从而鉴定可变剪切(exon-exon splice junctions)。

安装

最简单的安装方法,注意版本

  • 下载Bowtie、TopHat、Cufflinks的二进制发布包,解压到相同的目录
  • 下载samtools,make,将生成的可执行samtools程序也cp到同一个目录
  • 增加该目录到PATH

参数与使用

Usage: tophat [options]* <index_base> <reads1_1[,…,readsN_1]> [reads1_2,…readsN_2]

  • -o  输出目录,默认值为 “./tophat_out”。
  • –solexa-quals/solexa1.3-quals 质量编码,关于质量编码格式请参考《Fastq格式详解
  • -p 线程数,默认值为单线程1.,可以使用多线程
  • -G/–GTFSupply TopHat with a set of gene model annotations and/or known transcripts, as a GTF 2.2 or GFF3 formatted file.指定已有转录本信息
  • –no-novel-juncs 不查找新的可变剪切
  • -r 比对时两成对引物间的距离中值。比如说,如果你的插入片段有300bp,而每个引物有50bp,那么r值就应该是200=(300+50*2)/2。没有默认值,如果是末端配对比对时这个值是必须的。
  • –mate-std-dev 末端配对时中间插入片段的长度的标准差,默认值为20bp

paired-end数据应该如何做 Continue reading

blast+与blast的差异

很早就打算写一下blast+的差异,因为用的少了,所以草稿了很久,无意看到一篇介绍这方面的文章,写的挺好,就转载一下,作为补充了。

BLAST已经成为序列比对软件代名词,且其词性也已经开始变化,诸如BLASTing之类的词汇在各种论文中已是屡见不鲜,可见其影响之深,使用之广,如同分子生物学领域中的PCR。

自从1997年释出现有的BLAST版本后,这十多年来,BLAST经历了多次的升级,功能、性能一版比一版好,相应的其Source code也被修改的凌乱不堪,难于维护,极大的限制了对BLAST进一步 的修改、功能提升。再加上NCBI C++ Toolkit项目的开展,促使BLAST的维护者们决定从头开始,重新编写 BLAST代码。

2009年7月,NCBI发布了BLAST升级版——BLAST+,BLAST+使用了BLAST的核心算法,延 续了BLAST的优势功能,发展并增强了如BLAST的fastacmd程序,新增了如update_blastdb.pl等 程序。下面简单列举此次修改的主要内容: Continue reading

使用速铂Aspera下载NGS数据

关于速铂Aspera

速铂Aspera是一套商业的高速文件传输解决方案,随着高通量数据的大量产生,从而对于大文件快速传输的需求,开始应用到生物领域,目前NCBI、EBI的SRA库都提供这样的服务。

传统的FTP、HTTP等数据传输协议都是基于TCP的,TCP在远距离数据传输中存在一些先天的不足,文件越大、距离越远,其丢包、延时等问题对于传输速度的影响就越大。速铂Aspera通过应用了一个名为fasp™ 的底层技术,替换了传统的TCP传输协议。它彻底克服了TCP固有瓶颈,实现了在各种共享和私有网络环境中传输速度的最大化。这种技术可以获得完美的传输效率,不为网络延迟和丢包所限制。并且,用户享有对传输速度以及不同传输流之间带宽共享的无以伦比的控制。不管网络距离和动态性能如何,即便是在最困难的网络条件下(例如卫星,无线和洲际远程链接),文件传输时间仍然可以得到保障。FASP具有内置的,完整的安全性,包括连接节点安全验证,传输中数据加密以及数据完整性验证。与FTP传输相比快了3-184倍。它可以灵活地部署在C/S 或B/S构架的应用上,并利用普通的IP网络最大限度地利用带宽进行高效传输。同时,它也有着极好的跨平台性,支持几乎所有的主要操作系统。该软件同时也 包含一种文件接力技术,使得在传输大量极小文件时,其效率与传输单个大尺寸文件有着相同的效率与速度。

一句话,远距离,大文件,Aspera优势巨大。

客户端的下载与安装

即便Aspera是商业软件,但是作为客户应用方(相对于NCBI),我们使用其客户对进行数据的上传与下载是不用支付费用的。

  • Aspera Connect下载, 下载地址:http://www.asperasoft.com/downloads,根据不同的操作系统,下载相应的版本。注意下载的是Aspera Connect。Aspera Connect
  • 安装,windows下直接双击,下一步,安装,注意安装目录有别与常规软件,安装目录为C:\Users\[usename]\AppData\Local\Programs\Aspera\Aspera Connect。linux下
sh aspera-connect-xx-linux-64.sh

浏览器下使用Aspera下载SRA数据(win 7) Continue reading

RNA-seq拼接结果数据提交NCBI

RNA-seq的拼接结果也可以向NCBI提交,第一次提交,还是费了不少事,这里简单总结一下。RNA-seq的拼接结果应该提交到TSA库,TSA全称Transcriptome Shotgun Assembly Sequence Database,TSA is an archive of computationally assembled sequences from primary data such as ESTs, traces and Next Generation Sequencing Technologies.

对于注释信息的要求

TSA数据提交前,首先需要将原始的序列提交到SRA数据库,与提交普通核酸、EST类似,TSA还需要提供DBlink关于BioProject、SRA接收号、BioSample;提供拼接信息以及对于拼接过程的描述。

对于序列的要求

  • 必须是原始的测序结果的拼接数据
  • 需要去除载体或者测序引物
  • 序列长度不能少于200bp
  • 序列不能包括太多的N,少于10%或者小于14个N Continue reading

Fastq格式详解

FASTQ是基于文本的,保存生物序列(通常是核酸序列)和其测序质量信息的标准格式。其序列以及质量信息都是使用一个ASCII字符标示,最初由Sanger开发,目的是将FASTA序列与质量数据放到一起,目前已经成为高通量测序结果的事实标准。

格式说明

FASTQ文件中每个序列通常有四行:

  1. 序列标识以及相关的描述信息,以‘@’开头;
  2. 第二行是序列
  3. 第三行以‘+’开头,后面是序列标示符、描述信息,或者什么也不加
  4. 第四行,是质量信息,和第二行的序列相对应,每一个序列都有一个质量评分,根据评分体系的不同,每个字符的含义表示的数字也不相同。
例如:
@SEQ_ID
GATTTGGGGTTCAAAGCAGTATCGATCAAATAGTAAATCCATTTGTTCAACTCACAGTTT
+
!''*((((***+))%%%++)(%%%%).1***-+*''))**55CCF>>>>>>CCCCCCC65

Continue reading

paired-end reads的拼接

Velvet中paired-end reads的拼接

文件格式

要将两头测序(paired-end)的reads放到同一个文件当中,fastq格式,必须成对的依次放置reads [interleaved],velvet是成对读取的,另外Velvet假设来自两头read是反向互补的,如果不是,需要用反向互补序列来代替第一个read。Fastq格式中paired-end reads的编号相同,但是其有/1或者/2的后缀,通过这种方式来标示paired-end reads。

如果两端测序的reads放在不同的两个文件中,可以使用Velvet提供的perl脚本shuffleSequences fasta.pl进行转换合并,命令格式如下:

> ./shuffleSequences_fasta.pl forward_reads.fa reverse_reads.fa output.fa

低质序列过滤

在拼接前,首要要进行去除低质序列、接头等预处理,比如使用FASTX-Toolkit中的fastq_quality_filter去除低质序列:

fastq_quality_filter  -q 20 –p 100 -i s_1_1_sequence.txt -o s_1_1_sequence.txt_filtered_q20_p100.fastq
fastq_quality_filter  -q 20 –p 100 -i s_1_2_sequence.txt -o s_1_2_sequence.txt_filtered_q20_p100.fastq

这样势必带来一个问题,有些paired-end的前面序列被剔除,有些后面的序列被剔除,paired-end序列无法成对的错落出现,下面需要做的就是必须将单独的reads挑出来,方法有很多,下面是其中一个: Continue reading

文章来源数据拼接的一下实践

进来发现许多文章随之发表的都是其原始的二代测序的结果,很少将拼接好的序列一并发布,当然做可能有许多原因,比如编辑没有要求、拼接结果是此生数据、或者增加了工作量等等,所以需要用数据,还得拼一把,当然也是好事,可以用新的方法和文章中的处理方法进行比较,可以对于结果有一个验证。许多时候,我们自己的测序的拼接结果都是公司一起做的,我们拿到的就是拼接好的结果。而又或者对于自己的测序数据各方面了解的都非常清楚,完全下载的数据,如何进行拼接,需要注意什么样的问题,如何进行结果的比较,这里进行一些总结。以一篇以测序数据拼接与数据分析为主题的文章为例(Illumina RNA-seq测序),从NCBI SRA下载数据对其进行拼接,使用的拼接软件是velvet。

Continue reading

也谈统计分析

最近看到一篇博客《最新研究质疑RNA测序数据的统计分析》,Nature上一篇关于印记基因研究的文章的数据分析方法以及结论遭到质疑,其统计分析及其整体的数据分析不够仔细,导致数据的假阳性排除方面不够严格。激起我很强的共鸣,或者说触及到了我们这些生物方面做数据分析的软肋。同样的研究对象,不同的人,不同的处理方法,会得到不一样的数据;同样的数据,不同的人,不同的分析流程,不同的平台软件,会得到不同的结果。而之所以少有这样的结论,或许因为证明别人的结论不正确要远远的比自己得出一个结论难,况且还有面子和气,以及项目资助等等。

而国内更为深有感触,而目前的普遍的形式,公司流水线的程序,给出一套结果,然后是下面交给一个学生负责,包括数据分析、思路拟定,文章撰写。而预想的结果就是能有些数据,凑成文章,发表就可以了,而至于结果的正确性,只能听天意了,如果碰巧测序样本制备的很好,测序的质量很好,数据又没有特殊的偏好性,也符合目前的归一化的分析流程,结果OK,但是却无法证明这个结果OK,也没有精力或者能力去证明。

再反思,自己从前到后学活的统计学,分析过程中接触的统计学,研究结论中使用的统计学,也从来都是一种形式,学过的统计学有多少的统计方法,用这样的一组名词解释另一组的名词。而分析过程而完全依赖于软件,哦,这个值应该大于多少,统计学上才有意义,结论分析也同样,现成的别人的套路,这样的分析,就得到这样的结果,这样的结果就说明了这样的结论。没有压根搞明白,为什么?

回过头,想再补补课,深奥的统计学,什么时候可以真正变成自己的武器,而不是一个模子,变成一种思想,而不是一种形式。当然也希望这方面的高手多奉献些为什么的帖子。

SAM格式定义

SAM是一种序列比对格式标准, 由sanger制定,是以TAB为分割符的文本格式。主要应用于测序序列mapping到基因组上的结果表示,当然也可以表示任意的多重比对结果。

不同的软件,不同的时期,不同的研究方向,都会创建一种或者多种格式标准,当然根据当时的需要,创建符合需求的标准,也是最容易的事情,而反过来想要真正的理解标准,也必须理解为什么要创建这样的标准,解决什么样的需要。我前面的有篇文章已经对于现有的多重比对的格式进行总结,但其更多的站在比较基因组学的角度。当我们去了解sam标准格式是什么的时候,就要思考既然以及有了这么多得标准,为什么还要定义SAM标准,当然拿所有的格式进行比较也并非易事,但是简单的对比,就可以了解其中一二,比如aln格式,是比对视图化的展示,存储的信息不够结构化,无法方便的作为另外程序的输入;表示信息的有限性,如果100个多重比对序列放到一个文件中,查看维护就会非常困难;还有些格式标准挺强大,但是太繁琐,同时不够灵活。那么反过来就是SAM格式的优点,那么SAM如何做到这一点的呢? Continue reading

生物信息之心态分析

每一个分析之旅,尽管分析千差外别,而那份感觉或者心情却有几份的相像,轻轻的,却能感觉到。人类基因组测序已经十多年了,而近几年随着测序技术的发展,测序的成本呈指数下降,而测序产生的数据量呈指数上升,数据分析的需求开始激增,这就是我们所处的这个行业这个时代。

常常看到关于2000年那轮互联网泡沫的描述,经历了一场泡沫,互联网还成就了许多著名的公司,真真切切的融入与影响着我们的生活。而生物领域,同样,我们更为真切的感觉到了分子生物学泡沫对我们一代或者几代人的影响,想当年分子生物学、生物工程专业是多么的热门,大学了纷纷增设相关专业或者现有的专业改名与之挂靠,而以至于现在生物学方面的硕士博士的泛滥,同样相对的待遇也是同样的低。“二十一世纪是生物学的时代”,多么响亮的口号,而目前就公司于行业产值,还丝毫没有看到这方面的影子。

基因、分子克隆,转基因,生物工程到如今基因组、蛋白质组、转录组等各种组学、以及二代测序、个性化医疗,更多的应用还属于科研领域。科研、课题、博导、文章、圈子以及种种的外围切割,就是那一滩水。对于病毒、癌症、生物发酵、育种等等,关于生命本身的认识,并没有因为数据的增加,而有什么本质的改变。而数据仅仅只是给了我们对于生命的另外一个层面的描述,就如同对于绿叶、红花、肌肉的描述一样。数据堆砌着,文章累计着,有人从时间推测着,百年了,该有所质的突破了。。。

铛铛…..,又有新数据了,从测序公司拿到或者网上下载到,期盼已久的数据,怀着无比的希望打开,开始迫不及待的进入分析之旅。使用各种工具,加载各种算法,机器轰鸣着,一番运行,似乎经过运算,我们就能读懂数据中的一切,而几天下来,给我们的确实相反的沮丧,心里或许开始抱怨糟糕的设计、糟糕的数据质量、糟糕的算法。。。和想想中的差之千里。而生活还得继续,带着沮丧之情,肩负着艰巨的使命,继续挖掘。

而几轮的沮丧之后,开始渐渐明白

一个混沌进入另一个混沌,而或许正因为太广,太容易目空一切

测序的物种多了,测序的数量多了,我们却没有看到突飞猛进,反而平添了几分忧愁。

信息本身的复杂性,看看研究最多的人类基因组,看看身披无数文章的模式生物

观察数据,也就是测序样本、测序数据本身的有限性

我们需要在两者中间取得一个平衡,得到一个可信的结论。

即便是突破也是层级式的,一层层,或许看到的还是那份混沌,最少你应该明白,那就是你所期望的。

告别沮丧,怀着敬畏,尝试着,不断寻找与尝试着挖掘

如果,测序如水,不在制约,有人说,分析将成为制约,制约则意味着价值,不知道价值变成现金的日子,还会不会太远。