How to get gene level of differential expression from 2 color channel array with limma?

Question

[MatthewP]

Hello, I can get probe level differential expression from 2 color channel array with _limma_ use code below. But I want to get gene level data, how can I get that? Thanks.

RG <- read.maimages(files = fileList, source = "agilent", path = rawDir, names = fileName)
RGb <- backgroundCorrect(RG, method = "normexp")
MA <- normalizeWithinArrays(RGb, method = "loess")
MA_Q <- normalizeBetweenArrays(MA, method = "Aquantile")
notControl <- MA_Q$genes$ControlType == 0
# remove Control probe
MA_F <- MA_Q[notControl, ]
MA_A <- avereps(MA_F, ID = MA_F$genes$ProbeName)
MA_A$targets$Cy3 <- rep_len(c("Resistant", "Sensitive"), length.out = 12)
MA_A$targets$Cy5 <- rep_len(c("Sensitive", "Resistant"), length.out = 12)
desigN <- modelMatrix(MA_A$targets, ref = "Resistant") %>% as.data.frame()
# add dye effect coef
desigN$DyeEffect <- rep_len(1, 12)

fit1 <- lmFit(MA_A, design = desigN)
fit2 <- eBayes(fit1)
DEGs <- topTable(fit2, coef = "Sensitive", number = Inf, genelist = MA_A$genes)
DEGs2 <- tibble::as_tibble(DEGs, rownames = "ProbeID")

Answer 1

[Gordon Smyth]

It is not possible to consolidate microarray probes into genes. If you insist on having only one probe for each gene then you have to choose one probe to keep for each gene (e.g., the most highly expressed) and discard the others. For most purposes (e.g., pathway analysis) a few duplicate gene symbols cause no problems so it is generally best to simply keep all the probes in the DE analysis.

PS. You don't need to specify genelist=MA_A$genes as that will be done automatically.

[MatthewP]

Hello Gordon Smyth. Do you know how to detect dye bias? I added dye effect in design model.

> desigN
           Sensitive DyeEffect
GSM4518466         1         1
GSM4518467        -1         1
GSM4518468         1         1
GSM4518469        -1         1
GSM4518470         1         1
GSM4518471        -1         1
GSM4518472         1         1
GSM4518473        -1         1
GSM4518474         1         1
GSM4518475        -1         1
GSM4518476         1         1
GSM4518477        -1         1

So I can estimate genes DEGs because diffferent condition or dye bias.

fit1 <- lmFit(MA_A, design = desigN)
fit2 <- eBayes(fit1)

DEGs <- topTable(fit2, coef = "Sensitive", number = Inf, genelist = MA_A$genes)
dyeDEGs <- topTable(fit2, coef = "DyeEffect", number = Inf)

But here DEGs has only 12 probes that adj.P.Val is less than 0.05 while dyeDEGs has 14368. So I want to check dye bias of this data, whether limma extracts correctly.

[Gordon Smyth]

You have detected dye bias. Dye effect = dye bias. What you see is what is to be expected. Agilent arrays in particular show a strong dye bias.

[MatthewP]

Yes I know dye effect = dye bias. Just this dataset condition "sensitive" has too little number of probes that expression differentially and "DyeEffect" has too many, makes me confusing.

[Gordon Smyth]

What you see is normal. Many genes are expected to show dye bias, because it is a technical rather than a biological effect, but for the biological comparison you should get far fewer DE genes.