library(affy)
library(tcltk)
filters<-matrix(c("CEL file", ".[Cc][Ee][Ll]", "All",".*"), ncol = 2,byrow = T)
cel.files <-tk_choose.files(caption = "Select CELs", multi =TRUE,filters= filters, index = 1)
data.raw<-ReadAffy(filenames = cel.files)
n.cel <-length(cel.files)
sampleNames(data.raw)
==============================================
sampleNames(data.raw)<-paste("S",1:n.cel,sep='')
pData(data.raw)$treatment <-rep(c("0h", "1h","24h", "7d"),each=2)
#生成0h,1h,24h,7d四个值依次重复两次所组成的数列,数列命名为treatment
pData(data.raw)
#指针pData函数读取文件
==============================================
1、计算基因表达量
eset.rma <-rma(data.raw)
eset.mas5<-mas5(data.raw)
%%注意rma的eset结果是经过对数变换的,而mas5的eset结果是原始信号强度。虽然表达量是用对数变换的信号值表示的,但是有些计算过程要用到未经变换的原始值,应该把它们都计算出来:
emat.rma.log2<-exprs(eset.rma)
emat.mas5.nologs <-exprs(eset.mas5)
emat.rma.nologs<-2^emat.rma.log2
emat.mas5.log2 <-log2(emat.mas5.nologs)
results.rma<-data.frame((emat.rma.log2[,c(1,3,5,7)] + emat.rma.log2[,c(2,4,6,8)])/2)
#计算平均值,并转换为数据框格式#计算表达量差异倍数
results.rma$fc.1h<- results.rma[,2]-results.rma[,1]
results.rma$fc.24h <- results.rma[,3]-results.rma[,1]
results.rma$fc.7d <- results.rma[,4]-results.rma[,1]
subset.logic <-results.rma$fc.7d>0
subset.data <- results.rma[subset.logic,]
#要注意的是逻辑向量的长度要和相应维度的数据长度一致,逻辑向量中为TRUE的就保留,FALSE的就丢弃。
hea