Entering edit mode
I am working on my portfolio in R, my dataset is Agilent Microarray data and i have issues with the gene ontology analysis. Below are my lines of code.
library(limma)
list.files(path = "C:/Users/User/Documents/Extracted_R")
alzh.port <- list.files(path = "C:/Users/User/Documents/Extracted_R", full.names = T)
adf <- read.maimages(files = alzh.port, source = "agilent.median", green.only = TRUE, other.columns = 'gIsWellAboveBG')
##########
library("GEOquery")
z <- getGEO("GPL16699")
alzhdf <- Table(z)
#Match the Elist to the columns on thegene annotation file alzhdf
adf$genes$Gene.ID <- (alzhdf$GENE_SYMBOL)
adf$genes$RefSeq <- (alzhdf$REFSEQ)
# Correct background and normalize between arrays
adf_bgc<- backgroundCorrect(adf, method = "normexp")
adf_nba <- normalizeBetweenArrays(adf_bgc, method = 'quantile')
#Gene filtering
Control <- adf_nba$genes$ControlType==1L
NoGene.ID <- is.na(adf_nba$genes$Gene.ID)
IsExpr <- rowSums(adf_nba$other$gIsWellAboveBG > 0) >= 4
adf_nbafilt <- adf_nba[!Control & !NoGene.ID & IsExpr, ]
dim(adf_nbafilt)
#Removing the unwanted columns
adf_nbafilt <- adf_nbafilt$genes[,c("ProbeName","Gene.ID",'RefSeq')]
head(adf_nbafilt)
#Differential Expression
Condition <- read.csv('C:/Users/User/Documents/BIS7017-B_R_Workshops/Portfolio/Input/Condition.csv')
Condition <- factor(c('Control', 'Control', 'Control', 'Control', 'Control', 'Control', 'Control', 'Control', 'Control', 'Control',
'Alzhemier','Alzhemier', 'Alzhemier','Alzhemier', 'Alzhemier', 'Alzhemier', 'Alzhemier', 'Alzhemier','Alzhemier' ),
levels = c ('Control', 'Alzhemier'))
table(Condition)
design <-model.matrix(~Condition, data=adf_nbafilt)
fit <- lmFit(data=adf_nbafilt,design)
fit <- eBayes(fit, proportion=0.01, stdev.coef.lim=c(0.1,4),
trend=FALSE, robust=FALSE, winsor.tail.p=c(0.05,0.1))
TableAlhz <- topTable(fit,coef=1,number=100)
summary(decideTests(fit))
AlhzGO <- goana(fit, coef=1, species="Hs", geneid="Gene.ID")
Thank you for the quick response. i ran the codes again and removed all the other arguments for eBayes except for fit. I ran the line of code for toptable with coef=4 as written in the limma guide but i encountered an error., which was the reason i used coef = 1, previously.
fit <- eBayes(fit)
topTable(fit,coef=4,n=20)
Error in fit$coefficients[, coef] : subscript out of bounds
Can you pls help identify where i went wrong. I look forward to your response. Thank you
Gordon said to let coef to its default value: