每月一生信流程之rnaseqGene

每月一生信流程栏目灵感来自于《铁汉1991》博客的《每日一生信》,他那个时候介绍的主要是生信基础知识,包括数据结构,数据格式,数据库资源,计算机基础等等,所以每天都可以进步,每天都有成果。这些基础知识已经被分享的七七八八了,所以我这里推陈出新,来一个每月一生信流程,陪生信技能树的粉丝们一起进步!
前面我们推荐两个值得用一个月时间来学习的流程,包括
但是RNA-seq的分析肯定远不止那些啦,拿到基于基因的表达矩阵固然可以根据转录组经典表达量矩阵下游分析大全 里面的R包和代码进行统计可视化,但是表达矩阵并不是凭空产生,上游分析也需要我们有一定的认知,本次我们介绍的流程就会涵盖这些知识点。(很多朋友会下意识的认为RNA-seq数据的上游分析必然是基于Linux,其实也是可以使用bioconductor的全部R包来完成的哦!)
今天学习rnaseqGene
流程的原标题是:RNA-seq workflow: gene-level exploratory analysis and differential expression
  • http://www.bioconductor.org/packages/release/workflows/vignettes/rnaseqGene/inst/doc/rnaseqGene.html
这个同样的也没有中文版,讲解了把测序的FASTQ文件比对到参考基因组,或者使用salmon进行定量。介绍了2种得到表达矩阵的方法,涵盖的R包或者软件如下:
在R里面安装这个bioconductor流程
if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install("rnaseqGene")

因为原文写的实在是太详细了,我这里就不拷贝粘贴了,大家直接去阅读即可
全部目录如下;
  • 1 Introduction
  • 1.1 Experimental data
  • 2 Preparing quantification input to DESeq2
  • 2.1 Transcript quantification and tximport / tximeta
  • 2.2 Quantifying with Salmon
  • 2.3 Reading in data with tximeta
  • 2.4 DESeq2 import functions
  • 2.5 SummarizedExperiment
  • 2.6 Branching point
  • 3 The DESeqDataSet object, sample information and the design formula
  • 3.1 Starting from SummarizedExperiment
  • 3.2 Starting from count matrices
  • 4 Exploratory analysis and visualization
  • 4.1 Pre-filtering the dataset
  • 4.2 The variance stabilizing transformation and the rlog
  • 4.3 Sample distances
  • 4.4 PCA plot
  • 4.5 PCA plot using Generalized PCA
  • 4.6 MDS plot
  • 5 Differential expression analysis
  • 5.1 Running the differential expression pipeline
  • 5.2 Building the results table
  • 5.3 Other comparisons
  • 5.4 Multiple testing
  • 6 Plotting results
  • 6.1 Counts plot
  • 6.2 MA-plot
  • 6.3 Gene clustering
  • 6.4 Independent filtering
  • 6.5 Independent Hypothesis Weighting
  • 7 Annotating and exporting results
  • 7.1 Exporting results
  • 7.2 Plotting fold changes in genomic space
  • 8 Removing hidden batch effects(去除批次效应哦!)
  • 8.1 Using SVA with DESeq2
  • 8.2 Using RUV with DESeq2
  • 9 Time course experiments (时间序列转录组数据)
流程代码
我这里其实更加推崇我自己GitHub的airway数据集的分析代码:https://github.com/jmzeng1314/GEO/tree/master/airway_RNAseq 节选如下:
library(airway)

## 表达矩阵来自于R包:  airway
# 如果当前工作目录不存在文件:'airway_exprSet.Rdata' 
# 就运行下面 if 代码块内容,载入R包airway及其数据集airway后拿到表达矩阵和表型信息
if(!file.exists('airway_exprSet.Rdata')){
  library(airway)
  data(airway)
  exprSet=assay(airway)
  group_list=colData(airway)[,3]
  save(exprSet,group_list,file = 'airway_exprSet.Rdata')
}
# 大家务必注意这样的代码技巧,多次存储Rdata文件。
# 存储后的Rdata文件,很容易就加载进入R语言工作环境。
load(file = 'airway_exprSet.Rdata')
table(group_list)
# 下面代码是为了检查
if(T){
  colnames(exprSet)
  pheatmap::pheatmap(cor(exprSet))
  group_list
  tmp=data.frame(g=group_list)
  rownames(tmp)=colnames(exprSet)
  # 组内的样本的相似性理论上应该是要高于组间的
  # 但是如果使用全部的基因的表达矩阵来计算样本之间的相关性
  # 是不能看到组内样本很好的聚集在一起。
  pheatmap::pheatmap(cor(exprSet),annotation_col = tmp)
  dim(exprSet)
  # 所以我这里初步过滤低表达量基因。
  exprSet=exprSet[apply(exprSet,1, function(x) sum(x>1) > 5),]
  dim(exprSet)

exprSet=log(edgeR::cpm(exprSet)+1)
  dim(exprSet)
  # 再挑选top500的MAD值基因
  exprSet=exprSet[names(sort(apply(exprSet, 1,mad),decreasing = T)[1:500]),]
  dim(exprSet)
  M=cor(log2(exprSet+1))
  tmp=data.frame(g=group_list)
  rownames(tmp)=colnames(M)
  pheatmap::pheatmap(M,annotation_col = tmp)
  # 现在就可以看到,组内样本很好的聚集在一起
  # 组内的样本的相似性是要高于组间
  pheatmap::pheatmap(M,annotation_col = tmp,filename = 'cor.png')

library(pheatmap)
  pheatmap(scale(cor(log2(exprSet+1))))

}

# 以上代码,就是为了检查R包airway及其数据集airway里面的表达矩阵和表型信息

有一个质控的讨论,值得大家看看哦!
不同数据变换公式的差异
学习这样的流程是需要一定背景知识的

首先是LINUX学习

我在《生信分析人员如何系统入门Linux(2019更新版)》把Linux的学习过程分成6个阶段 ,提到过每个阶段都需要至少一天以上的学习:
  • 第1阶段:把linux系统玩得跟Windows或者MacOS那样的桌面操作系统一样顺畅,主要目的就是去可视化,熟悉黑白命令行界面,可以仅仅以键盘交互模式完成常规文件夹及文件管理工作。
  • 第2阶段:做到文本文件的表格化处理,类似于以键盘交互模式完成Excel表格的排序、计数、筛选、去冗余,查找,切割,替换,合并,补齐,熟练掌握awk,sed,grep这文本处理的三驾马车。
  • 第3阶段:元字符,通配符及shell中的各种扩展,从此linux操作不在神秘!
  • 第4阶段:高级目录管理:软硬链接,绝对路径和相对路径,环境变量
  • 第5阶段:任务提交及批处理,脚本编写解放你的双手
  • 第6阶段:软件安装及conda管理,让linux系统实用性放飞自我

然后是R学习

我在在生信分析人员如何系统入门R(2019更新版) 里面给初学者的知识点路线图如下:
  • 了解常量和变量概念
  • 加减乘除等运算(计算器)
  • 多种数据类型(数值,字符,逻辑,因子)
  • 多种数据结构(向量,矩阵,数组,数据框,列表)
  • 文件读取和写出
  • 简单统计可视化
  • 无限量函数学习

必备书籍及视频

书籍贪多不烂,下面2本必买,读5遍以上
视频必须强推生信技能树近30万学习量的基础合辑:
生信技能树关于RNA-seq上下游数据分析的教程的确不少了
因为做目录确实很浪费时间,差不多就下面这些,大家先学习吧:
(0)

相关推荐

  • 【1115.】PCA、PLS-DA、OPLS-DA到底啥关系?

    导读 代谢组学是一门十分火热的研究领域,在代谢组学的数据分析中,下图你一定不陌生.        图中的同一种颜色所覆盖的区域代表同一组的样本数据,如果同组的样本都聚在一起,不同组的数据分布在不同的颜 ...

  • RNA seq汇总篇,一文掌握RNA seq

    RNA测序(RNA-seq)在过往十年里逐渐成为全转录组水平分析差异基因表达和研究mRNA差异剪接必不可少的工具.RNA-seq帮助大家对RNA生物学的理解会越来越全面:从转录本在何时何地转录到RNA ...

  • AAAI 2020 | 首个使用 NAS 设计的 GCN,达到动作识别SOTA,代码将开源

    本文作者为52CV群友彭伟,现就读奥卢大学(芬兰) 博士二年级. Github:https://github.com/xiaoiker 知乎:https://www.zhihu.com/people/ ...

  • 手把手教你下载TCGA数据(代码+视频+答疑+服务)

    现在TCGA数据下载的代码满天飞,例如以使用TCGAbiolinks下载为例: if (!requireNamespace("BiocManager", quietly = TRU ...

  • 双重差分法(DID)安慰剂检验的做法:随机抽取500次?

    "安慰剂"(placebo)一词来自医学上的随机实验,比如要检验某种新药的疗效.此时,可将参加实验的人群随机分为两组,其中一组为实验组,服用真药:而另一组为控制组,服用安慰剂(比如 ...

  • 每月一生信流程之RNAseq123

    目前bioconductor社区有27个流程,早在2015/2016年我组织生信菜鸟团小伙伴建设bioconductor中文社区的时候就想系统性的学习和分享,一晃四五年过去了, 我们的biocondu ...

  • 每月一生信流程之rnaseqDTU(差异转录本)

    每月一生信流程栏目灵感来自于<铁汉1991>博客的<每日一生信>,他那个时候介绍的主要是生信基础知识,包括数据结构,数据格式,数据库资源,计算机基础等等,所以每天都可以进步,每 ...

  • 人一生要流多少眼泪,仿佛是苍天注定的更是无法逃避的

    ​2021年4月27日,晚21点左右到医院探望病人,路过急症室听到一女人一边嚎啕大哭,嘴里一边说着伤感的话,不知为什么我真想去安慰她,又怕引起她更加伤心,不经意间想到人一生要流多少眼泪,仿佛是苍天注定 ...

  • 诗词丨醉漾轻舟,信流引到花深处

    点绛唇·桃源 宋·秦观 醉漾轻舟,信流引到花深处. 尘缘相误,无计花间住. 烟水茫茫,千里斜阳暮. 山无数,乱红如雨.不记来时路. 译 文 醉酒后荡着小船,任流水引着轻舟飘向花草深处.现实世界的名利缠 ...

  • 我好恨!!!为什么没能早点看到这封信[流...

    我好恨!!!为什么没能早点看到这封信[流...

  • 西贝莜面村,以客户为中心的流程之美

    今天是周末,去西贝,用餐流程又简化了.一直有一个梦想...有时间的时候,我要整理开发一个餐厅的流程故事.扯远了,先说说今天感受到的西贝流程之美吧. 一.沙漏,背后是端到端流程绩效管理能力 说到餐厅行业 ...

  • 检测bam文件的完整度-流程之殇

    本来以为有bai文件就说明是流程运行是完整的,事实上我还是太年轻,最近处理一个 600个病人的肿瘤WES数据,走流程过程发现卡在CNVKIT,部分样本出现了:   File "pysam/l ...

  • 可变剪切流程之suppa的diffsplice太慢了

    我们发布了转录组产品线的一些服务,连接生信技能树粉丝群体的数据分析工程师和有数据分析的科研人员: 明码标价之转录组常规测序服务(仅需799每个样品) 明码标价之普通转录组上游分析 明码标价之转录组下游 ...

  • 如信流鼻血要仰头止血民间几大偏方会害了宝宝!

    多多今年一岁半,平时由奶奶照顾.不久前刚学会走路,多多揉着奶奶的时间去卫生间,跑到饮水机前玩.不小心碰到热水开关时,多多的左臂被热水烫伤.李奶奶立即找到牙膏抹在多多的手上.半个多小时后,她发现多多的手 ...