When I output the logCPM values for my genes, they are identical across contrasts. However, when I output the logFold values these seem to be correctly calculated for each contrast. Is this a bug or am I not understanding what the logCPM is that is returned by glmLRT?
Code:
library(edgeR) # Read in the quantification data data <- read.table('File.txt', header = TRUE, row.names = 1) keep <- rowSums(cpm(data) > 2) >= min(summary(3)) data <- data[keep,] y <- DGEList(counts = data) y <- calcNormFactors(y) rownames(design) <- colnames(y) y <- estimateGLMCommonDisp(y, design, verbose = TRUE) y <- estimateGLMTrendedDisp(y, design) y <- estimateGLMTagwiseDisp(y, design) # Fit the linear model to the design matrix fit <- glmFit(y, design) # Make Contrast Matrix my.contrasts <- makeContrasts(....... # Different contrasts levels = design) # Loop to calculate the contrasts for (i in 1:dim(my.contrasts)[2]) { # Create the log ratios lrt <- glmLRT(y, fit, contrast = my.contrasts[,i]) tempLR <- lrt$table[2] if(exists("LogCPM")){ LogCPM< <- cbind(LogCPM, tempLR) } else{ LogCPM <- tempLR } } # rows of LogCPM are identical?!?
Comment