基因芯片分析流程

source("https://bioconductor.org/biocLite.R")
biocLite("affy")
biocLite("hgu133plus2.db")#注释文件
library(affy)
library("hgu133plus2.db")

#读取CEL文件
celFiles <- list.celfiles(path = "/yourdataPATH/", full.names=TRUE)
data.affy <- ReadAffy(filenames = celFiles)
data.affy

#Normalize data with RMA
data.rma <- rma(data.affy)

expr.rma <- exprs(data.rma) # format as table

# Convert gene names
Annot <- data.frame(REFSEQ=sapply(contents(hgu133plus2REFSEQ), paste, collapse=", "),
                    SYMBOL=sapply(contents(hgu133plus2SYMBOL), paste, collapse=", "),
                    DESC = sapply(contents(hgu133plus2GENENAME), paste, collapse=", ")  )
# Merge data frames together (like a database table join)
all <- merge(Annot, expr.rma, by.x=0, by.y=0, all = TRUE)
#remove probe ID,gene REFSEQ and DESC and NA 
all<-all[,c(-1:-2,-4)]
all<-all[which(all[,1] != "NA"),]
#use mean of probe ID as gene expression
all<-aggregate(.~SYMBOL,all,mean)
#change transcript ID to gene symbol in rownames
rownames(all)<-all$SYMBOL

source("https://bioconductor.org/biocLite.R")
biocLite("oligo")
library("oligo")

library("oligo")
#文件位置
data.dir <- "../../../test/GSE81580_RAW-2/"
##CEL文件读取
celfiles <- list.files(data.dir, "\\CEL.gz$")
data.raw <- read.celfiles(filenames = file.path(data.dir, celfiles))
#表达量计算
data.eset <- oligo::rma(data.raw)  #含背景处理、归一化和表达量计算
data.exprs <- exprs(data.eset)  #提取表达量矩阵
write.csv(data.exprs,"../../../test/GSE81580expr.csv")

biocLite("ALL")
library(ALL)
show(data.raw)