大家好，又见面了，我是你们的朋友全栈君。如果您正在找激活码,请点击查看最新教程,关注关注公众号 “全栈程序员社区” 获取激活教程,可能之前旧版本教程已经失效.最新Idea2022.1教程亲测有效,一键激活。

Jetbrains全系列IDE使用 1年只要46元售后保障童叟无欺

转录组分析流程(有参和无参de novo)

获得测序数据，Fastq格式，称之为Raw data。

质量检测

比对Mapping

Quantification|Quantitation

差异表达分析

补充：开始项目之前，先确立合理的文件目录结构。

【1】Raw Data 处理

理论知识

高通量测序之所以能够能够达到如此高的通量的原因就是他把原来几十M，几百M，甚至几个G的基因组通过物理或化学的方式打算成几百bp的短序列，然后同时测序。在测序过程中，机器会对每次读取的结果赋予一个值，用于表明它有多大把握结果是对的。从理论上都是前面质量好，后面质量差。并且在某些GC比例高的区域，测序质量会大幅度降低。因此，我们在正式的数据分析之前需要对分析结果进行质控。

Fastq 文件

测序给的“原始数据”，称之为Raw Data。

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

FASTQ文件中以四行最为一个基本单元，并对应一条序列的测序信息，各行记录信息如下：

第一行记录序列标识以及相关的描述信息，以‘@’开头，为了保证后续分析软件能够区分每条序列，单个序列的标识必须具有唯一性；
第二行为碱基序列；
第三行以‘+’开头，后面是序列标示符、描述信息，或者什么也不加；
第四行，是质量信息，长度和第二行的序列相对应，每一个序列都有一个质量评分，根据评分体系的不同，每个字符的含义表示的数字也不相同。
碱基质量得分与错误率的换算关系： Q = -10log10p（p表示测序的错误率，Q表示碱基质量分数）
ASCII值与碱基质量得分之间的关系：
Phred64 Q=ASCII转换后的数值-64
Phred33 Q=ASCII转换后的数值-33
目前illumina使用的碱基质量格式为phred+33, 和Sanger的质量基本一致（老数据建议查看清楚再进行后续处理）。

SRA 文件

Sequence Read Archive (SRA) makes biological sequence data available to the research community to enhance reproducibility and allow for new discoveries by comparing data sets. The SRA stores raw sequencing data and alignment information from high-throughput sequencing platforms, including Roche 454 GS System®, Illumina Genome Analyzer®, Applied Biosystems SOLiD System®, Helicos Heliscope®, Complete Genomics®, and Pacific Biosciences SMRT®.

具体参考NCBI关于SRA的介绍

SRA文件转换成fastq文件

我们从NCBI等数据库下载的原始数据很多为SRA格式，需要转换成fastq文件。常用工具为：NCBI SRA Toolkit
简单介绍使用方法（具体的坑踩过才能记住！！！）

单个文件转换(PE)

fastq-dump –gzip –split-3 -O outputdir -A file1.sra

多个文件批量转换

for I in seq 56 62
do
fastq-dump –gzip –split-3 -O ./fastq/ -A SRR35899${I}.sra
done

–split-3：如果是双端测序数据，则输出两个文件，如果不是则只输出一个文件
–gzip：输出格式为gzip的压缩文件（fastqc软件可以直接识别gzip压缩的文件）
-A：accession序列号，输入的文件
-O：outdir输出文件夹，指定输出路径

FastQC（测序质量分析）：多个文件批量进行

$ fastqc -q -t 4 -o ./fastqc_result/ *.fastq.gz &
-t：调用核心数目
-q：安静运行，运行过程中不会生成报告，在结束时将报告生成一个文件
-o ../path/to/file ：文件输出位置
*. fq.gz：输入文件，当前目录下所有名字中有“ .fq.gz ”的文件

查看QC结果

1、单个查看：鼠标双击打开html文件查看

2、批量查看：使用 moltiqc软件: moltiqc *fastqc.zip

Fastqc结果报告关注重点：
1）.basic statistics
2）.per base sequence quality
3）.per base sequcence content
4）.adaptor content
5）.sequence duplication levels
最主要的几个指标是GC含量，Q20和Q30的比例以及是否存在接头（adaptor）、index以及其他物种序列的污染等。

Per base sequence quality，各位置碱基质量，每个read各位置碱基的测序质量。横轴碱基的位置，纵轴是质量分数，Quality score=-10log10p（p代表错误率），所以当质量分数为40的时候，p就是0.0001，质量算高了。红色线代表中位数，蓝色代表平均数，黄色是25%-75%区间，触须是10%-90%区间。若任一位置的下四分位数低于10或者中位数低于25，出现“警告”；若任一位置的下四分位数低于5或者中位数低于20，出现“失败，Fail”。

Per base sequence content，碱基分布，对所有reads的每一个位置，统计ATCG四种碱基的分布，横轴为位置，纵轴为碱基含量，正常情况下每个位置每种碱基出现的概率是相近的，四条线应该平行且相近。当部分位置碱基的比例出现bias时，即四条线在某些位置纷乱交织，往往提示我们有overrepresented sequence的污染。本结果前10个位置，每种碱基频率有明显的差别，说明有污染。当任一位置的A/T比例与G/C比例相差超过10%，报‘WARN’；当任一位置的A/T比例与G/C比例相差超过20%，报‘FAIL’。

Adapter content，接头含量：在此处可查看数据中使用接头信息
具体接头查询地点：github-fastqc 或者：Download common Illumina adapters

TruSeq Universal Adapter：
AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT
Illumina Small RNA 3p Adapter：
1 ATCTCGTATGCCGTCTTCTGCTTG

另外，一般而言RNA-Seq数据在sequence deplication levels 未通过是比较正常的。毕竟一个基因会大量表达，会测到很多遍

Trimmomatic处理fastq数据

根据FastQC质控报告，利用Trimmomatic软件处理数据。trimmomatic，是java软件，所以直接Google找到其官网，然后下载二进制版本解压即可使用！这个软件设计就是为了illumina的测序数据的，因为它自带的adaptor文件有限！一般只去除TruSeq Universal Adapter 这个接头，运行的时候，不报错才算是成功的！
官网有例子，很简单的：http://www.usadellab.org/cms/?page=trimmomatic

Paired End:
java -jar trimmomatic-0.35.jar PE -phred33 input_forward.fq.gz input_reverse.fq.gz output_forward_paired.fq.gz output_forward_unpaired.fq.gz output_reverse_paired.fq.gz output_reverse_unpaired.fq.gz ILLUMINACLIP:TruSeq3-PE.fa:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:36 ## 所以只需要把参数放对位置即可！

This will perform the following:

Remove adapters (ILLUMINACLIP:TruSeq3-PE.fa:2:30:10)
Remove leading low quality or N bases (below quality 3) (LEADING:3)
Remove trailing low quality or N bases (below quality 3) (TRAILING:3)
Scan the read with a 4-base wide sliding window, cutting when the average quality per base drops below 15 (SLIDINGWINDOW:4:15)
Drop reads below the 36 bases long (MINLEN:36)

！！！记住指定接头文件一定要用全路径哦！！！
当然，也可以用cutadapt这个python软件来去除接头的，但是它有一个弊端，需要自己指定接头文件。

【2】下载参考基因组及基因注释

在 UCSC 下载 hg19 参考基因组；

从 gencode 数据库下载基因注释文件，并且用 IGV 去查看感兴趣的基因的结构，比如TP53,KRAS,EGFR 等等。

截图几个基因的 IGV 可视化结构

下载 ENSEMBL，NCBI 的 gtf，也导入 IGV 看看，截图基因结构

了解 IGV 常识
来源于生信技能树：http://www.biotrainee.com/forum.php?mod=viewthread&tid=1750#lastpost
Y大宽链接：https://www.jianshu.com/p/02a92e4ead4b

【3】比对Mapping

把整理好的数据和参考基因组序列进行比对，寻找每个reads的最佳匹配位置。可以使用HISAT2，tophat2，STAR等软件。

具体参照Y叔介绍序列比对：Hisat2 pipeline（全流程，不仅仅是Mapping）
Hisat2(比对,sam) ==> SAMtools (bam, sort, index) ==> htseq-count(reads 计数, 表达矩阵) ==> R(合并多个矩阵) ==> DEseq2(筛选差异表达基因并注释)

【4】表达量分析（软件：HTSeq，RSEM等）

RSEM (RNA-Seq by Expectation-Maximization)

RSEM^1,2 is an RNA-Seq transcript quantification program developed in 2009. You need a server with Linux/Mac OS. To run RSEM, your server should have C++, Perl and R installed. In addition, you need at least one aligner to align RNA-Seq reads for you. RSEM can call Bowtie, Bowtie 2 or STAR for you if you have them installed. Last but not least, you need to install the latest version of RSEM.
github链接

HTseq

htseq-count 是一款用于reads计数的软件，他能对位于基因组上的一些单位的reads数进行统计，这里所说的单位主要是指染色体上的一组位置区间（我们常见的就是gene exon）。
HTSeq follows install conventions of many Python packages. In the best case, it should install from PyPI like this:

pip install HTSeq

If this does not work, please open an issue on Github .

基本用法：
htseq-count [options] <alignment_file> <gff_file>

htseq-counts跟bedtools的区别

bedtools不管三七二十一，只要你的reads比对到基因组的坐标跟目的基因坐标有交叉，就算你一个reads，不需要管你是不是multiple mapping的。但是htseq就谨慎很多，而且还可以挑选model，一般来说，它会把multiple mapping的reads归类到 not unique aligned里面。

最后htseq-counts使用的时候有一些参数尤其需要注意：

软件官网说明书： https://htseq.readthedocs.io/en/release_0.11.1/

参考gtf文件可以是gencode或者是ensembl数据库的，但是尤其要注释chr的问题，而且版本问题，gtf/gff格式无所谓。比对后的文件一定要进行sort，推荐一定要sort -n，根据reads的name来sort
-f sam/bam: 如果对bam文件进行counts，必须保证服务器的python安装了正确的pysam模块
-r name/pos: 一般情况下我们的bam都是按照参考基因组的pos来sort的，但是这个软件默认却是reads的name，很坑，一般建议重新把bam文件sort一下，而不是选择 -r pos，因为-r pos实在是太消耗内存了。
-s yes/no/reverse, 这也是巨坑的参数，默认是yes，一般人拿到的数据都是no，所以千万要注意！！！
-t 选择gff/gtf文件的第3列，一般是exon，也可以是gene，transcript ，这个很少调整的。
-i 这个需要修改，不然默认是ensembl的基因ID，一般人看不懂，可以改为gene_name，前提是你的gff文件里面有gene_name这个属性。

【5】差异表达分析（edgeR, DEseq2, EBseq等）

DEseq2

安装（坑自己踩！！！）

if (!requireNamespace(“BiocManager”))
install.packages(“BiocManager”)
BiocManager::install(“DEseq2”)

需要准备2个table：一个是countData，一个是colData

DESeq流程

> dds <- DESeqDataSetFromMatrix(mycounts, colData, design= ~ condition)
> dds <- DESeq(dds)
> res = results(dds, contrast=c("condition", "control", "treat"))
> 或下面命令
> res= results(dds)
> res = res[order(res$pvalue),]

> head(res)
> summary(res)
> 所有结果先进行输出
> write.csv(res,file="All_results.csv")
> table(res$padj<0.05)

提取差异基因

> diff_gene_deseq2 <-subset(res, padj < 0.05 & abs(log2FoldChange) > 1)
    或
> diff_gene_deseq2 <-subset(res,padj < 0.05 & (log2FoldChange > 1 | log2FoldChange < -1))
> dim(diff_gene_deseq2)
> head(diff_gene_deseq2)
> write.csv(diff_gene_deseq2,file= "DEG_treat_vs_control.csv")

EBSeq

Bioconductor软件包安装：

     if (!requireNamespace("BiocManager", quietly = TRUE))
             install.packages("BiocManager")
     BiocManager::install("EBSeq", version = "3.8")

TPR relatively independent of sample size and presence of outliers.

Poor FDR control in most situations, relatively unaffected by outliers.

Medium computational time requirement, increases slightly with sample size.

标准化方法不同+检验不同=多种组合/软件，用之前需要结合自己的样本量来考虑，多参考有相似实验设计的文献，常用的方法都跑一下，自己评估下结果差异，再做定夺。（研究本来就是充满了不确定性，一切都只能用“可能性”来定义，所以，采用同样参数仍然无法完全重复出文献中的结果也是常见。即使你心有芥蒂…）by fatlady：https://www.jianshu.com/p/364650e8bd3e

其他软件包

DEGseq(http://www.bioconductor.org/packages/release/bioc/html/DEGseq.html)
NOISeq(http://www.bioconductor.org/packages/release/bioc/html/NOISeq.html)
Ballgown (https://github.com/alyssafrazee/ballgown)是R语言中基因差异表达分析的工具，能利用RNA-Seq实验的数据(StringTie, RSEM, Cufflinks)的结果预测基因、转录本的差异表达。然而Ballgown并没有不能很好地检测差异外显子，而 DEXseq、rMATS和MISO可以很好解决该问题。

【6】基因注释

正在纠结ID转换…..

【7】基因富集分析(gene set enrichment analysis, GSEA)

用clusterProfiler做富集分析

（1）Bioconductor安装软件包：http://www.bioconductor.org/install/

【7】数据可视化处理

根据不同的PIPELINE选择合适的方式(R)进行可视化。
比较好的教程：

https://www.jianshu.com/p/72c871663e62

【8】转录组分析流程举例

有参转录组分析

1.首推Y叔分享的教程：RNA-seq分析简洁版

文章很精辟，在此不赘述。

2. Hisat2-Stringtie-Ballgown

案例来源：Pertea, M., Kim, D., Pertea, G. M., Leek, J. T., & Salzberg, S. L. (2016). Transcript-level expression analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nature protocols, 11(9), 1650-67.doi: 10.1038/nprot.2016.095
数据和软件准备：

    $ sudo yum install axel        
    $ axel ftp://ftp.ccb.jhu.edu/pub/RNAseq_protocol/chrX_data.tar.gz
    $ tar –zxvf chrX_data.tar.gz
    #下载和安装miniconda
    $ wget https://repo.continuum.io/miniconda/Miniconda3-latest-Linux-x86_64.sh
    #下载完成后在终端中安装
    $ bash Miniconda-latest-Linux-x86_64.sh
    #按照提示安装，完成后
    $ source ~/.bashrc   #使以上的安装立即生效
    #输入以下命令检验miniconda是否安装成功
    $ conda list

利用conda install 软件名+版本号安装软件即可，我们需要安装hisat2、stringtie、samtools三个软件，安装的命令为：

    $ conda install hisat2
    $ conda install stringtie 
    $ conda install samtools

具体分析流程：
首先使用hisat2进行比对：

    $ for i in {188044,188104,188234,188245,188257,188273,188337,188383,188401,188428,188454,204916};do hisat2 -p 4 -x chrX_data/indexes/chrX_tran -1 chrX_data/samples/ERR${i}_chrX_1.fastq.gz -2 chrX_data/samples/ERR${i}_chrX_2.fastq.gz -S ERR${i}_chrX.sam
    done

脚本执行完可得到12个sam文件。
通过samtools将sam文件转换为bam文件，作为stringtie的输入文件，具体脚本为：

    $ for i in {188044,188104,188234,188245,188257,188273,188337,188383,188401,188428,188454,204916};do samtools sort -@ 4 -o ERR${i}_chrX.bam ERR${i}_chrX.sam done

此处sort默认输出的bam文件是按基因组位置排序，同样tophat的输出的bam文件也是按此顺序排序的，而sort -n 则是按reads的ID排序。

接下来利用stringtie对转录组进行组装，会针对每个bam文件生成一个gtf文件：

    $ for i in {188044,188104,188234,188245,188257,188273,188337,188383,188401,188428,188454,204916};do
    stringtie -p 8 -G ./genes/chrX.gtf -o ERR${i}_chrX.gtf -l ERR${i} ERR${i}_chrX.bam
    done

将输出的12个GTF文件的文件名录入到mergelist.txt文件中，然后利用软件stringtie将12个含有转录本信息的gtf文件合并成一个gtf。

$ stringtie --merge -p 8 -G ./genes/chrX.gtf -o stringtie_merged.gtf ./mergelist.txt

参数–merge 为转录本合并模式。在合并模式下，stringtie将所有样品的GTF/GFF文件列表作为输入，并将这些转录本合并/组装成非冗余的转录本集合。这种模式被用于新的差异分析流程中，用以生成一个跨多个RNA-Seq样品的全局的、统一的转录本。
接下来，重新组装转录本并估算基因表达丰度，并为ballgown创建读入文件。利用组装好的非冗余的转录本文件即stringtie_merged.gtf 和12个bam文件，执行下面的脚本:

 $ for i in {188044,188104,188234,188245,188257,188273,188337,188383,188401,188428,188454,204916};do
    stringtie -e -B -p 8 -G stringtie_merged.gtf -o ballgown/ERR${i}/ERR${i}_chrX.gtf ERR${i}_chrX.bam
  done

接下来要用到R语言分析（除dplyr，都为Bioconductor包）：

> library(RSkittleBrewer)
> library(ballgown)
> library(genefilter)
> library(dplyr)
> library(devtools)
#设置R语言的工作路径
> setwd("F:/data/R") 
#读取表型数据如右图所示
> read.csv("geuvadis_phenodata.csv")
> pheno_data<-read.csv("F:/data/R/geuvadis_phenodata.csv ")
#dataDir告知数据路径， samplePattern则依据样本的名字来，pheno_data    则指明了样本数据的关系，这个里面第一列样本名需要和ballgown下面的文件夹的样本名一样，不然会报错
> bg_chrX = ballgown(dataDir = “F:/data/R/ballgown", samplePattern = "ERR", pData=pheno_data) 
#滤掉低丰度的基因，这里选择过滤掉样本间差异少于一个转录本的数据
> bg_chrX_filt=subset(bg_chrX,"rowVars(texpr(bg_chrX))>1",genomesubset=TRUE)
#确认组间有差异的转录本，在这里我们比较male和famle之间的基因差异，指定的分析参数为“transcripts”，主变量是“sex”，修正变量是“population”，getFC可以指定输出结果显示组间表达量的foldchange。
> result_transcripts=stattest(bg_chrX_filt,feature = "transcript",covariate = "sex",adjustvars = c("population"), getFC=TRUE,meas="FPKM")
#查看有差异转录本的输出结果，如下图
 > result_transcripts

确定各组间有差异的基因

>result_genes=stattest(bg_chrX_filt,feature = "gene",covariate = "sex",adjustvars = c("population"), getFC=TRUE,meas="FPKM")
#为组间有差异的转录本添加基因名，如下图
> result_transcripts = data.frame(geneNames=ballgown::geneNames(bg_chrX_filt), geneIDs = ballgown::geneIDs(bg_chrX_filt), result_transcripts)
#按照p-value从小到大排序
> result_transcripts=arrange(result_transcripts,pval) 
> result_genes=arrange(result_genes,pval)
#把两个结果写入到csv文件中
> write.csv(result_transcripts, "chrX_transcript_results.csv",row.names=FALSE)
> write.csv(result_genes, "chrX_gene_results.csv",row.names=FALSE)
#从以上的输出中筛选出q值小于0.05的transcripts和genes，即性别间表达有差异的基因
> subset(result_transcripts,result_transcripts$qval<0.05)
> subset(result_genes,result_genes$qval<0.05)
#赋予调色板五个指定颜色
> tropical= c('darkorange', 'dodgerblue','hotpink', 'limegreen', 'yellow') 
> palette(tropical)
#以FPKM为参考值作图，以性别作为区分条件
> fpkm = texpr(bg_chrX,meas="FPKM")
#方便作图将其log转换，+1是为了避免出现log2(0)的情况
> fpkm = log2(fpkm+1) 
> boxplot(fpkm,col=as.numeric(pheno_data$sex),las=2,ylab='log2(FPKM+1)')

#查看单个转录本在样品中的分布，如图，并绘制箱线图
> ballgown::transcriptNames(bg_chrX)[12] 
> ballgown::geneNames(bg_chrX)[12]
>plot(fpkm[12,] ~ pheno_data$sex, border=c(1,2), main=paste(ballgown::geneNames(bg_chrX)[12],' : ', ballgown::transcriptNames(bg_chrX)[12]),pch=19, xlab="Sex", ylab='log2(FPKM+1)')
>points(fpkm[12,] ~ jitter(as.numeric(pheno_data$sex)), col=as.numeric(pheno_data$sex))

# plotTranscripts函数可以根据指定基因的id画出在特定区段的转录本，可以通过sample函数指定看在某个样本中的表达情况，这里选用id=1750, sample=ERR188234
> plotTranscripts(ballgown::geneIDs(bg_chrX)[1750], bg_chrX, main=c('Gene MSTRG.575 in sample ERR188234'), sample=c('ERR188234'))

3. 用tophat和cufflinks分析RNAseq数据

这个流程比较老，参考文章：

Trapnell, C., Roberts, A., Goff, L., Pertea, G., Kim, D., Kelley, D. R., Pimentel, H., Salzberg, S. L., Rinn, J. L., … Pachter, L. (2012). Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nature protocols, 7(3), 562-78. doi:10.1038/nprot.2012.016
https://www.nature.com/articles/nprot.2012.016

无参转录组分析

de novo组装

1. Trinity

Trinity是由 the Broad Institute 开发的转录组de novo组装软件，由三个独立的软件模块组成：Inchworm,Chrysalis和Butterfly。三个软件依次来处理大规模的RNA-seq的reads数据。Trinity的简要工作流程为：

Inchworm: 将RNA-seq的原始reads数据组装成Unique序列；
Chrysalis: 将上一步生成的contigs聚类，然后对每个类构建Bruijn图；
Butterfly: 处理这些Bruijn图，依据图中reads和成对的reads来寻找路径，从而得到具有可变剪接的全长转录子，同时将旁系同源基因的转录子分开。

目前最常用 Illumina RNA-Seq data de novo组装软件。案例：

Trinity download

Build Trinity by typing ‘make’ in the base installation directory.
Assemble RNA-Seq data like so:

     Trinity --seqType fq --left reads_1.fq --right reads_2.fq --CPU 6 --max_memory 20G

Find assembled transcripts as: ‘trinity_out_dir/Trinity.fasta’

    Trinity --seqType fq --max_memory 100G --CPU 50 --min_kmer_cov 3 --left   FCHK2FVCCXY_L3_WHDAVllgEAAARAAPEI-96_1.fq.gz,FCHK2FVCCXY_L3_WHDAVllgEAABRAAPEI-97_1.fq.gz,FCHK2FVCCXY_L3_WHDAVllgEAABRAAPEI-97_1.fq.gz  --right FCHK2FVCCXY_L3_WHDAVllgEAAARAAPEI-96_2.fq.gz,FCHK2FVCCXY_L3_WHDAVllgEAABRAAPEI-97_2.fq.gz,FCHK2FVCCXY_L3_WHDAVllgEAACRAAPEI-98_2.fq.gz --output gongtong_trinity_out  --group_pairs_distance 230 --no_version_check  --verbose --min_contig_length 250 --min_glue 3 --no_distributed_trinity_exec ~/bio/trinityrnaseq-Trinity-v2.4.0/trinity-plugins/parafly/bin/ParaFly -c recursive_trinity.cmds -CPU 50 -v

对于Trinity得到的转录本序列，Trinity官网推荐取每条基因中最长的转录本作为Unigene，并用这些Unigene去进行功能注释，但是在计算表达量的时候我们依然会用到所有的转录本。
后续流程参考：无参转录组GO、KEGG富集分析——diamond+idmapping+GOstats

注意事项

Trinity分步运行
当数据量比较大的时候，trinity运行的时间会很长，同时，内存不够等情况出现的时候有可能程序运行崩溃。最好是分步运行。下一步会接着前一步进行下去。

Stage 1: generate the kmer-catalog and run Inchworm: –no_run_chrysalis

Stage 2: Chrysalis clustering of inchworm contigs and mapping reads: –no_run_quantifygraph

Stage 3: Chrysalis deBruijn graph construction: –no_run_butterfly

Stage 4: Run butterfly, generate final Trinity.fasta file. (exclude –no_ options)

计算资源
Ideally, you will have access to a large-memory server, ideally having ~1G of RAM per 1M reads to be assembled (but often, much less memory may be required).

The assembly from start to finish can take anywhere from ~1/2 hour to 1 hour per million reads (your mileage may vary). 个人记录了一次，使用dell服务器，64GB RAM，24 threads : 53M 的reads，运行了16.5h（平均3.2M/h），内存使用峰值为43G.

2. est2Assembly 和gsassembler 软件

案例： 2014 巴西橡胶树的研究，是一个综合多组织样本的RNA库，ployT建库，454测序，用的是est2Assembly 和gsassembler 软件做组装，用 NCBI RefSeq, Plant Protein Database 做注释，因为没有分组，所以不必做差异分析，只需要找SNV和SSR标记即可，最后也是做GO/KEGG注释。

参考资料

发布者：全栈程序员-用户IM，转载请注明出处：https://javaforall.cn/167627.html原文链接：https://javaforall.cn

【正版授权，激活自己账号】： Jetbrains全家桶Ide使用，1年售后保障，每天仅需1毛

【官方授权正版激活】： 官方授权正版激活支持Jetbrains家族下所有IDE 使用个人JB账号...

转录组分析学习笔记（持续补充）

转录组分析流程(有参和无参de novo)

【1】Raw Data 处理

理论知识

Fastq 文件

SRA 文件

SRA文件转换成fastq文件

FastQC（测序质量分析）：多个文件批量进行

查看QC结果

Trimmomatic处理fastq数据

This will perform the following:

【2】下载参考基因组及基因注释

【3】比对Mapping

【4】表达量分析（软件：HTSeq，RSEM等）

RSEM (RNA-Seq by Expectation-Maximization)

HTseq

【5】差异表达分析（edgeR, DEseq2, EBseq等）

DEseq2

EBSeq

其他软件包

【6】基因注释

【7】基因富集分析(gene set enrichment analysis, GSEA)

用clusterProfiler做富集分析

【7】数据可视化处理

【8】转录组分析流程举例

有参转录组分析

1.首推Y叔分享的教程：RNA-seq分析简洁版

2. Hisat2-Stringtie-Ballgown

3. 用tophat和cufflinks分析RNAseq数据

无参转录组分析

de novo组装

1. Trinity

注意事项

2. est2Assembly 和gsassembler 软件

参考资料

发表回复

转录组分析学习笔记（持续补充）

转录组分析流程(有参和无参de novo)

【1】Raw Data 处理

理论知识

Fastq 文件

SRA 文件

SRA文件转换成fastq文件

FastQC（测序质量分析）：多个文件批量进行

查看QC结果

Trimmomatic处理fastq数据

This will perform the following:

【2】下载参考基因组及基因注释

【3】比对Mapping

【4】表达量分析（软件：HTSeq，RSEM等）

【5】差异表达分析（edgeR, DEseq2, EBseq等）

DEseq2

其他软件包

【6】基因注释

【7】基因富集分析(gene set enrichment analysis, GSEA)

用clusterProfiler做富集分析

【7】数据可视化处理

【8】转录组分析流程举例

有参转录组分析

1.首推Y叔分享的教程：RNA-seq分析简洁版

2. Hisat2-Stringtie-Ballgown

3. 用tophat和cufflinks分析RNAseq数据

无参转录组分析

de novo组装

1. Trinity

注意事项

2. est2Assembly 和gsassembler 软件

参考资料

相关推荐

发表回复