GEOquery_gene_expression.r

library(Biobase)
library(GEOquery)


# load series and platform data from GEO

gset <- getGEO("GSE34412", GSEMatrix =TRUE)
gset<- gset[[1]]

# make proper column names to match toptable 
fvarLabels(gset) <- make.names(fvarLabels(gset))
featureNames(gset)[1:10]
sampleNames(gset)[1:10]

# log2 transform
ex <- exprs(gset)
qx <- as.numeric(quantile(ex, c(0., 0.25, 0.5, 0.75, 0.99, 1.0), na.rm=T))
LogC <- (qx[5] > 100) ||
  (qx[6]-qx[1] > 50 && qx[2] > 0) ||
  (qx[2] > 0 && qx[2] < 1 && qx[4] > 1 && qx[4] < 2)
if (LogC) { ex[which(ex <= 0)] <- NaN
            exprs(gset) <- log2(ex) }

annotation(gset)
gpl <- getGEO('GPL8269')
#get annotation tabole
head(gpl)
gpl_ann<-Table(gpl)
names(gpl_ann)
head(gpl_ann)
gpl_ann[gpl_ann$"Blast Gene Symbol"=="TRIM29",]

exprs(gset)[1:4,1:4]

rownames(exprs(gset))<- gpl_ann$"Blast Gene Symbol"

index<- rownames(exprs(gset))=="TRIM29"

exprs(gset)[index,][1:4,1:4]

#only the first 25 xenografts data
Trim29_exprs<- exprs(gset)[index,][,1:25]
Trim29_exprs
#set the margin 
par(mar=c(8,8,4,2))
boxplot(Trim29_exprs, ylab="Trim29 log2 transformed expression level", las=2) #make the lable vertical by specifying las=2
#############################################
names(pData(gset))

pData(gset)[,c(1,2,28)]

names(fData(gset))

fData(gset)$ID[1:20]
fData(gset)$GB_ACC[1:20]
fData(gset)$Blast.Gene.Symbol[1:20]
"TRIM29" %in% fData(gset)$Blast.Gene.Symbol
which (fData(gset)$Blast.Gene.Symbol=="TRIM29")