Hello. I have a few groups of samples from RNA-Seq. I'd like to find out which group has the biggest change in expression comparing to others. I'm following this code.
myCPM = cpm(count)
thresh = myCPM > 0.5
keep <- rowSums(thresh) >= 20
counts.keep = countdata[keep,]
dgeObj = DGEList(counts.keep)
logcounts <- cpm(dgeObj,log=TRUE)
dgeObj <- calcNormFactors(dgeObj)
sampleinfo = samples_data
group = sampleinfo$type
dgeObj = estimateDisp(dgeObj)
design <- model.matrix(~group, data=dgeObj$samples)
fit = glmFit(dgeObj, design)
head(coef(fit))
results <- glmLRT(fit, coef=2)
topTags(results)
robust<- estimateGLMRobustDisp(dgeObj, design)
summary(robust$tagwise.dispersion)
Is there a way to compare all groups with each other automatically? I know that there's a glmQLFTest function but it compares only the first group to the rest of them and what I want to obtain is the overall comparison result.
And btw how to change the group that is the intercept? I changed the levels but it has no effect on the design matrix and the intercept is still the first group.
I read the EdgeR guide and searched the biostars but I'm still not sure how to solve these problems. Thank you in advance.