1、输入代码:括号内为“作者名/R包名称”
devtools::install_github("xjsun1221/tinyarray")
2、本地安装:从github官网上R包界面下载到本地,并放到当前工作目录下,使用如下命令:
devtools::install_local("tinyarray-master.zip",upgrade = F,dependencies = T)
注意:本地安装时要写全文件名称而并不只是R包名称!
3、安装历史版本的R包:下载后本地安装。
1、热图:输入数据是数值型矩阵/数据框;
2、散点图和箱线图
箱线图:输入数据是一个连续性向量和一个有重复值的离散型向量;可用来展示单个基因在两组之间的表达量差异
3、火山图:芯片差异分析的起点是一个取过log的表达矩阵
横坐标:Foldchange(FC):处理组平均值/对照组平均值
logFoldchange(logFC):Foldchange取log2
纵坐标:-log10(P.Value)
4、主成分分析:PCA样本聚类图,用于“预实验”,简单查看组间是否有差别
图上的点代表样本(中心除外),点与点之间的相对距离代表样本差异
dim1,dim2后的数据值表示主成分1和主成分2各能解释数据变化方向
理想实验设计:同一分组聚成一簇(组内重复好)、中心点之间有距离(组件差异大)
原理是用探针的表达量来代表基因的表达量,因此我们得到的表达矩阵的行名是探针名,需要转化为gene symol(常说的基因名)
需要说明的是symol与EnterzID并非一一对应,增加或损失部分属于正常。
本次为标准流程,要求:1)表达芯片数据;2)二分组数据。
1、准备工作:使用循环以及判断语句安装R包:
options("repos"="https://mirrors.ustc.edu.cn/CRAN/")
if(!require("BiocManager")) install.packages("BiocManager",update = F,ask = F)
options(BioC_mirror="https://mirrors.ustc.edu.cn/bioc/")
cran_packages <- c('tidyr',
'tibble',
'dplyr',
'stringr',
'ggplot2',
'ggpubr',
'factoextra',
'FactoMineR',
'devtools',
'cowplot',
'patchwork',
'basetheme',
'paletteer',
'AnnoProbe',
'ggthemes',
'VennDiagram',
'tinyarray')
Biocductor_packages <- c('GEOquery',
'hgu133plus2.db',
'ggnewscale',
"limma",
"impute",
"GSEABase",
"GSVA",
"clusterProfiler",
"org.Hs.eg.db",
"preprocessCore",
"enrichplot")
for (pkg in cran_packages){
if (! require(pkg,character.only=T) ) {
install.packages(pkg,ask = F,update = F)
require(pkg,character.only=T)
}
}
for (pkg in Biocductor_packages){
if (! require(pkg,character.only=T) ) {
BiocManager::install(pkg,ask = F,update = F)
require(pkg,character.only=T)
}
}
2、以GSE56649为例,如何判断它是哪种数据类型?
GEO官网——GEO accession内输入GSE编号——“Experiment type”显示为“Expression Profiling by array”
其中Platforms为GPL平台编号
3、表达矩阵的获取可以直接从网页Download family处下载Series Matrix File(s)(内存大小通常以M为单位),然后放在工作目录下;或者使用以下代码:
library(GEOquery)
gse_number = "GSE56649"
#换数据集时只需要修改引号内的内容即可
eSet <- getGEO(gse_number, destdir = '.', getGPL = F)
#destdir = '.'表示下载数据到工作目录下
class(eSet)
length(eSet)
eSet = eSet[[1]]
Tips:R语言里面广义的对象:向量、矩阵、数据框及列表等
R语言里面狭义的对象:1)由R包作者定义的以某种模式组织的数据;表示出自“Biobase”包的“ExpressionSet”对象;2)有帮助文档:?ExpressionSet;3)简单对象可以用@或者$取子集,复杂对象则需要函数提取子集。
4、提取表达矩阵:
exp <- exprs(eSet)
dim(exp)
exp[1:4,1:4]
exp = log2(exp+1)
boxplot(exp)
如何判断表达矩阵可用?
1)取值范围:是否取过log?;
2)箱线图的平直程度:上下四分位数线、中位数线;
异常数据的处理:①舍弃;②
exp=limma::normalizeBetweenArrays(exp)
3)有无负值:取过log后可以有负值,但没取过log则为错误数据,如果有一半负值表示作者上传数据时已经做过标准化,可用来画图,但不可以做差异分析。
5、使用pData函数提取临床信息,并使临床信息表格pd行名与表达矩阵列名完全一致
pd <- pData(eSet)
p = identical(rownames(pd),colnames(exp));p
if(!p) exp = exp[,match(rownames(pd),colnames(exp))]
tips:match函数的使用:谁在外边,谁在后边。
x[match(y,x)]
#表示以y为模版给c调顺序
6、提取芯片平台编号,并将以上输入数据储存为Rdata文件
gpl_number <- eSet@annotation;gpl_number
save(gse_number,pd,exp,gpl_number,file = "step1output.Rdata")
——————————————————图及代码均引自“生信技能树”小洁老师————————————————————
原创声明:本文系作者授权腾讯云开发者社区发表,未经许可,不得转载。
如有侵权,请联系 cloudcommunity@tencent.com 删除。
原创声明:本文系作者授权腾讯云开发者社区发表,未经许可,不得转载。
如有侵权,请联系 cloudcommunity@tencent.com 删除。