#三个包同时加载
library(AnnoProbe)
library(GEOmirror)
library(GEOquery)
#下载获取GSE13507数据
gset=AnnoProbe::geoChina('GSE13255')
gset
#检查表达式集
eSet=gset[[1]] # 提取表达矩阵
probes_expr <- exprs(eSet)
dim(probes_expr)
#如果样本已经进行了log2处理，这里就不需要再进行归一化
probes_expr=log2(probes_expr+1) #表达矩阵行log2+1归一化处理
## 对表达芯片的探针进行基因注释
(gpl=eSet@annotation)
checkGPL(gpl)
printGPLInfo(gpl)
probe2gene=idmap(gpl) #获得探针注释
genes_expr <- filterEM(probes_expr,probe2gene )
class(genes_expr)
write.table(genes_expr,file="genes_expr.txt",quote=F,sep = "\t")
保存后，对数据进行整理，只保留想要分析的样本数据。
如果样本已经进行了log2处理，这里就不需要再进行归一化
报错
没有的话就自己整理
下载平台注释文件
整理成如下的文件，然后自己替换探针名
自行替换，然后再分组做差异分析
anno=read.table("anno.txt",sep="\t",header=T,check.names=F)
probes_expr<-as.data.frame(probes_expr)
ID<-as.data.frame(rownames(probes_expr))
probes_expr<-cbind(ID,probes_expr)
colnames(probes_expr)[1]<-'ID'
common<-merge(anno,probes_expr,by='ID')
common<-common[,-1]
write.table(common,file="genes_expr.txt",quote=F,sep = "\t")
差异分析
## 分组
group_list=factor(c(rep('Normal',91),rep('Cancer',137))) #自己定义分组和数量
table(group_list)
#limma包差异分析
library(limma)
library(data.table)
rt<-fread("genes_expr.txt",sep="\t")
rt=as.matrix(rt)
rownames(rt)=rt[,1]
exp=rt[,2:ncol(rt)]
dimnames=list(rownames(exp),colnames(exp))
data=matrix(as.numeric(as.matrix(exp)),nrow=nrow(exp),dimnames=dimnames)
data=avereps(data)
data<-as.data.frame(data)
nrow(data)
data=data[rowMeans(data)>0.5,]
design=model.matrix(~factor(group_list))
design
fit=lmFit(data,design)
fit=eBayes(fit)
DEG=topTable(fit,coef=2,n=Inf)
head(DEG)
write.table(DEG,file="DEG.txt",quote=F,sep = "\t")
绘图
##############绘制火山图
need_deg=data.frame(symbols=rownames(DEG), logFC=DEG$logFC, p=DEG$P.Value)
deg_volcano(need_deg,1) #第一种图
deg_volcano(need_deg,1,logFC_thred = 0) #设置logFC范围
deg_volcano(need_deg,2) # 第二种图
#################热图
deg_heatmap(DEG,data,group_list)
deg_heatmap(DEG,data,group_list,5) #显示前5对，数字可以自己定义
###########箱图
check_diff_genes('TLR7',data,group_list)
check_diff_genes('TLR7',data,group_list)+
stat_compare_means(method = "t.test") #加个P值，手到擒来
临床数据下载
临床数据从sangerbox里面进行下载
