c语言处理rna序列,RNAseq 完整操作流程以及后续例子操作

2.软件以及流程以及代码： 本实验使用常规的Tophat2 比对

(1)样本制备、建库：总RNA —> polyA富集mRNA —> 打断 —> 随机引物反转录成cDNA —> 末端修复、加A、加接头

RNA提取：方法-使用trizol-based方法或者试剂盒方法；

富集mRNA，除去rRNA等RNA:依据为-在测mRNA过程中，首先要去除rRNA。以人为例，在抽提的总RNA中，95%的RNA是rRNA，2%的RNA是mRNA，剩下的则是lncRNA、microRNA、siRNA等。rRNA整个人类当中是非常保守的，在各个组织器官中也是非常稳定的，因此这些测序结果对我们的研究是没有用处的。mRNA则是RNA中比较重要的部分。

具体：以人的为例：总RNA跟带有Poly(T)探针的磁珠结合-洗脱结合的mRNA-用Mg离子溶液打碎mRNA-随机引物反转第一条链的cDNA，之后再合成第二条cDNA,获得双链cDNA-对双链cDNA末端修复，加A加接头-片段选择，PCR扩增、纯化(如果样本中存在污染物，则需要结合试剂盒进一步纯化)。

###一般在会测序前对总RNA进行一次质检，根据电泳质检结果中的18S和28S(rRNA)两个峰的高度以及峰的尖度来判断RNA的质量，峰越高越尖(RIN > 8.0)表示RNA的完整度越好。当然，浓度以及A260/A280比值也是需要的。

1、真核生物种常规的去除rRNA的方法是通过oligo(dT)富集带有polyA尾的mRNA来实现的，2、不含有polyA尾的转录本序列以及存在部分降解的总RNA样本，所以这种方法针对福尔马林(Formalin-Fixed)样本和FFPE(Paraffin-Embedded)石蜡包埋样本是不适用的，否则对获得样本中最全面的转录本信息会产生显著影响，一般采用需结合RiboZero、RiboMinus等是结合来开展去除。针对FFPE样本还有结合双链特异性核酸酶构建文库来降低后续测序数据中的rRNA序列比例的。

建库：去除rRNA之后获得的mRNA进行构建文库，先对mRNA打断再进行反转录的文库构建方法。之后反转的cDNA再末端修复到平末端，加上ployA和接头。

###当然，里面涉及的蛋白相关知识，例如蛋白变性失活、键的破坏、复性和盐析等涉及的蛋白四级结构、各种化学键和相互作用的内容。

(2)测序：SE测序-单端测序；PE测序-双端测序；一般现在用的是PE测序。

##里面有测序的方法、测序原理、不同品牌的区别

.sra-> .fastq  代码：fastq-dump  –gzip  –split-3 –O ../fastq/ -A ../xx.sra

(3)质量控制:当然用IGV也可以做一些质控。

Fastqc  自己命名.fastq.gz   -o保存文件夹/

(4)数据与处理：

比对质量过滤、修剪：trimmomatic  PE输入文件.1.fastq.gz 输入文件.2.fastq.gz  paired1.fq.gz unpaired.1.fq.gz  paired2.fq.gz unpaired.2.fq.gz AVGQUAL:20 MINLEN:50(删去质量小于20，且删除读段小于50的片段)

除去那些测不出来或者说未被识别的序列 标记为N，很多的话配对时需要删除这样的序列：prinseq-lite.pl –fastq read1.fastq

–fastq2 read2.fastq –ns_max_n 2 –out_god nfiltered –out_bad null –no_qual_header

– log –verbose (删除每条序列上含有2个N的序列)

去接头：trimmomatic  PE输入文件.1.fastq.gz 输入文件.2.fastq.gz  paired1.fq.gz unpaired.1.fq.gz  paired2.fq.gz unpaired.2.fq.gzILLUMINACLIP:TruSeq2 –PE .fa:2:30:10:1:true (删除TruSeq2接头，允许有2个不匹配，回文剪接阈值是30，简单剪接阈值10，回文模式检测到的最低接头长度是1，反向读段被保留-默认它被删除)

重复：prinseq-lite.pl–fastq read1.fastq –fastq2 read2.fastq –derep_min 101

–out_god nfiltered –out_bad

null –no_qual_header – log –verbose

(5)比对or 的de novo：

tophat2 –o 保存文件夹/ –transcriptome-index  转录本索引bt2/ -p 8 基因组索引bt2/  输入文件

(6)比对结果注释，量化RSeQc

1.比对统计项：bam_stat.py –i 比对文件

2.比对到基因组各个原件上的情况：

Read_distribution.py –r 基因组的gtf对应的bed文件–I 比对后的文件bam

3.转录本覆盖度是不是有偏差

geneBody-coverage.py –r 基因组的gtf对应的bed文件–I 比对后的文件bam

4.测序深度上表达丰度检测

RPKM_saturation.py  –r 基因组的gtf对应的bed文件–I 比对后的文件bam

5. 剪接点

Junction_annotation.py–r 基因组的gtf对应的bed文件–I 比对后的文件bam

6.剪接点饱和状态检测

Junction_saturation.py –r 基因组的gtf对应的bed文件–I 比对后的文件bam

量化：

1.   基因水平上-reads落在那些基因上  HTSeq软件

比对结果按照名字排序：samtools sort –n 比对文件 排序后的文件

量化：htseq-count –f bam –stranded=no 排过序的文件bam基因组索引文件gtf > counts.txt

2.   转录水平量化：cufflinks –G 基因组文件gtf  -b基因组–u –p 8 比对后为排序的文件bam  -o保存文件夹

3.   外显子水平  DEXSeq软件

扁平化(即，先将基因组注释文件扁平化，拉开距离，形成不重叠的外显子区域，进而将比对数据进行比对)

Python2 dexseq_prepare_annotation.py 基因组索引文件gtf  扁平化文件.gtf

量化  python2 dexseq_count.py –p yes–s no –r name扁平化的文件.gtf 按照名字排序后的sam文件 输出文件.txt

(7)组装：

Cufflinks -P 8 –O保存文件夹/ 比对后的文件bam/

De novo:

Trinity.pl –seqType fq –JM 10G –left 1.fq–-right 2.fq –CPU 4

(8)差异表达分析：R语言作图

参照：http://www.360doc.com/content/18/0309/18/33459258_735717104.shtml        http://www.360doc.com/content/17/0911/22/46164085_686360709.shtml   http://www.biotrainee.com/thread-1084-1-1.html

后续理解上：

步骤分为如下过程，为了好理解：

1. 数据下载：geo数据库的测序数据；UCSC 站基因组数据(chromFa.tar.gz)；gencode 站基因组注释文件(gtf)；hisat2 站的index文件；RSeQc软件 站作覆盖度的文件(.bed)

2.数据比对：hisat2  只是用了index文件； samtools的格式转换–排序等

3.比对结果质检：RSeQc的质检；

4.read计数，read归类于哪个基因(HTSeq)、转录(cufflinks)、外显子区域(DEXSeq)；

5. 差异表达分析：归类后要进行比较了；DESeq2包

6. 富集分析：因为某一个基因不能仅凭借表达量多少就判断多少了，如果是低表达的突然高表达一点呢，高表达但是比正常状态却少了呢，这时候就需要看看富集到一起的结果是咋样的了。Y叔的clusterfiler包。GO KEGG

7. 其他分析：聚类分析图；主成分分析；(R语言中高级技能包括四类：广义线性模型、聚类分析、时间序列、主成分分析)；别的嘛，需要什么就怎么操作吧。

文章知识点与官方知识档案匹配，可进一步学习相关知识C技能树首页概览113320 人正在系统学习中 相关资源：丝柏人像美肤处理软件CPAC Imaging Pro 3绿化汉化破解版