This might be a quite basic question, but has been troubling me for a while. How is fitting a negative binomial glm in R different from edgeR DE analysis. I think edgeR allows one to model a different variance for all genes, which is not feasible in the standard glm.nb function in R. Other than this, do they both use the same model? I have drug response data for cases and controls and would like to estimate the effect of drug as well as the disease. With glm.nb I can analyze all the samples together using drug and disease as two factors (0/1). However, edgeR only allows pairwise comparisons, such as between cases with no drug vs controls with no drug OR cases with drug vs controls with drug to look at the effect of disease, and cases with no drug vs cases with drug OR controls with no drug vs controls with drug to look at the effect of drug. This tremendously reduces the power of the data I have. I just want to confirm if I am correct in understanding this. Any comments will be highly appreciated.