GEO数据挖掘—处理Affymetrix的表达量芯片的cel文件2

#01_start_GEO.R
rm(list = ls())
options(timeout = 100000) 
options(scipen = 20)
#传统下载方式
library(GEOquery)
eSet = getGEO("GSE6410", destdir = '.', getGPL = F)
#2.试试geoChina,只能下载2019年前的表达芯片数据
#library(AnnoProbe)
#eSet = geoChina("GSE6410") #选择性代替第8行
#研究一下这个eSet
class(eSet)
length(eSet)
eSet = eSet[[1]] 
class(eSet)
#(1)提取表达矩阵exp
exp <- exprs(eSet)
dim(exp)
range(exp)#看数据范围决定是否需要log，是否有负值，异常值，如有负值，结合箱线图进一步判断
#
#exp = log2(exp+1) #需要log才log，不需要log要注释掉这一句
#⭐第三个要检查的地方👇
boxplot(exp,las = 2) #看是否有异常样本

#step0 cel.R
rm(list = ls())

library(oligo)
celFiles<-list.celfiles('GSE6410_RAW/',listGzipped=T,
                          full.name=TRUE)
celFiles
options(BioC_mirror="https://mirrors.westlake.edu.cn/bioconductor")
exon_data<-oligo::read.celfiles(celFiles)
dim(exprs(exon_data))
exprs(exon_data)[1:4,1:4]
###标准化(一步完成背景校正、分位数校正标准化和RMA(robustmultichipaverage)表达整合)
exon_data_rma<-oligo::rma(exon_data)
exp_probe<-Biobase::exprs(exon_data_rma)
exp_probe[1:4,1:4]

dat=exp_probe
dim(dat)#看一下dat这个矩阵的维度
dat[1:4,1:4]#查看dat这个矩阵的1至4行和1至4列，逗号前为行，逗号后为列

exp <- dat
range(exp)
boxplot(exp,las = 2)

summary(as.vector(exp))
save(dat,file = "step0.cel.Rdata")

load("step0.cel.Rdata")
colnames(dat)
exp <- dat