Entering edit mode
6.6 years ago
lech.kaczmarczyk
▴
50
Hi All, I tried to use EdgeR GLM fit function on my data, but I only get some clearly incorrect results. I spent hours trying to spot my error. Limma voom lmFit worked fine on the same data. Here is the script. I figured that something is wrong with how the comparison are set, but don't know how and what to fix. You are my last hope.
> y_TU_TagGLM <- DGEList(txi_TPMs$counts[,3:8], samples = SampleTable[3:8,1:3], genes = rownames(txi_TPMs))
> colnames(y_TU_TagGLM) <- SampleTable[3:8,1]
> designGLM <- model.matrix(~0 + ~Type, y_TU_TagGLM$samples)
> colnames(designGLM) <- gsub("Type","", colnames(designGLM))
> names(y_TU_TagGLM)
[1] "counts" "samples"
> designGLM
Affinity_Purified_vGlut2 Unbound
vGluT2_TU1 1 0
vGluT2_TU2 1 0
Unbound_Gad2_TU1 0 1
Unbound_Gad2_TU2 0 1
Unbound_vGluT2_TU1 0 1
Unbound_vGluT2_TU2 0 1
attr(,"assign")
[1] 1 1
attr(,"contrasts")
attr(,"contrasts")$Type
[1] "contr.treatment"
> y_TU_TagGLM <- calcNormFactors(y_TU_TagGLM)
> y_TU_TagGLM <- estimateDisp(y_TU_TagGLM, designGLM)
>
> fitTU <- glmFit(y_TU_TagGLM, designGLM)
> lrtTU <- glmLRT(fitTU, coef=2)
> ttTU <- topTags(lrtTU, n=nrow(y_TU_TagGLM), p.value=0.1, adjust.method="BH", sort.by="PValue")
> tt10 <- as.data.frame(topTags(lrtTU, n=nrow(y_TU_TagGLM)))
> tt10
logFC logCPM LR PValue FDR
ENSMUSG00000100131 -25.788249 7.6141957137 3246.33983 0.000000e+00 0.000000e+00
ENSMUSG00000020140 -17.153947 3.5624865899 1522.37536 0.000000e+00 0.000000e+00
ENSMUSG00000008333 -16.922668 3.9244238753 1526.39683 0.000000e+00 0.000000e+00
ENSMUSG00000051920 -16.903855 3.8515603204 1679.09047 0.000000e+00 0.000000e+00
ENSMUSG00000025940 -16.880893 3.9362217669 1504.04685 0.000000e+00 0.000000e+00
ENSMUSG00000061559 -16.780391 3.7090760077 1660.99975 0.000000e+00 0.000000e+00
ENSMUSG00000026027 -16.769318 4.0481712621 1557.80775 0.000000e+00 0.000000e+00
ENSMUSG00000043999 -16.732672 3.5288708577 1584.29262 0.000000e+00 0.000000e+00
ENSMUSG00000005802 -16.727090 3.8779942891 1560.42342 0.000000e+00 0.000000e+00
> lrtTU
An object of class "DGELRT"
$coefficients
Affinity_Purified_vGlut2 Unbound
ENSMUSG00000000001 -11.35335 -11.54869
ENSMUSG00000000003 -18.50960 -18.50960
ENSMUSG00000000028 -13.05009 -13.57928
ENSMUSG00000000037 -14.73369 -14.83902
ENSMUSG00000000049 -15.92172 -16.37713
34834 more rows ...
$fitted.values
vGluT2_TU1 vGluT2_TU2 Unbound_Gad2_TU1 Unbound_Gad2_TU2 Unbound_vGluT2_TU1 Unbound_vGluT2_TU2
ENSMUSG00000000001 340.280203 312.299640 126.5293480 59.4623969 27.2048770 40.2990287
ENSMUSG00000000003 0.000000 0.000000 0.0000000 0.0000000 0.0000000 0.0000000
ENSMUSG00000000028 62.149383 57.038963 16.5101174 7.7589205 3.5498145 5.2583982
ENSMUSG00000000037 11.325269 10.394015 4.5840940 2.1542924 0.9856189 1.4600133
ENSMUSG00000000049 3.267397 2.998726 0.8712465 0.4094418 0.1873253 0.2774881
34834 more rows ...
$deviance
ENSMUSG00000000001 ENSMUSG00000000003 ENSMUSG00000000028 ENSMUSG00000000037 ENSMUSG00000000049
3.094813 0.000000 4.449726 6.071378 4.895752
34834 more elements ...
$method
[1] "oneway"
$unshrunk.coefficients
Affinity_Purified_vGlut2 Unbound
ENSMUSG00000000001 -1.135413e+01 -1.154966e+01
ENSMUSG00000000003 -1.000000e+08 -1.000000e+08
ENSMUSG00000000028 -1.305436e+01 -1.358616e+01
ENSMUSG00000000037 -1.475686e+01 -1.486754e+01
ENSMUSG00000000049 -1.599990e+01 -1.652796e+01
34834 more rows ...
$df.residual
[1] 4 4 4 4 4
34834 more elements ...
$design
Affinity_Purified_vGlut2 Unbound
vGluT2_TU1 1 0
vGluT2_TU2 1 0
Unbound_Gad2_TU1 0 1
Unbound_Gad2_TU2 0 1
Unbound_vGluT2_TU1 0 1
Unbound_vGluT2_TU2 0 1
attr(,"assign")
[1] 1 1
attr(,"contrasts")
attr(,"contrasts")$Type
[1] "contr.treatment"
$offset
[,1] [,2] [,3] [,4] [,5] [,6]
x 17.1839 17.09809 16.39013 15.635 14.85305 15.24598
attr(,"class")
[1] "compressedMatrix"
attr(,"repeat.row")
[1] TRUE
attr(,"repeat.col")
[1] FALSE
$dispersion
[1] 0.08905732 1.05908953 0.35835305 1.31137167 0.80502574
34834 more elements ...
$prior.count
[1] 0.125
$samples
group lib.size norm.factors Name Type Mouse
vGluT2_TU1 vGluT2_TU1 27893527 1.0408040 vGluT2_TU1 Affinity_Purified_vGlut2 vGluT2
vGluT2_TU2 vGluT2_TU2 24908681 1.0696863 vGluT2_TU2 Affinity_Purified_vGlut2 vGluT2
Unbound_Gad2_TU1 Unbound_Gad2_TU1 13566535 0.9675523 Unbound_Gad2_TU1 Unbound Gad2
Unbound_Gad2_TU2 Unbound_Gad2_TU2 6416577 0.9613711 Unbound_Gad2_TU2 Unbound Gad2
Unbound_vGluT2_TU1 Unbound_vGluT2_TU1 2907967 0.9705309 Unbound_vGluT2_TU1 Unbound vGluT2
Unbound_vGluT2_TU2 Unbound_vGluT2_TU2 4201912 0.9949466 Unbound_vGluT2_TU2 Unbound vGluT2
$prior.df
[1] 5.966951
$AveLogCPM
[1] 3.4006651 -2.7721353 0.8862272 -1.0586601 -2.0704172
34834 more elements ...
$table
logFC logCPM LR PValue
ENSMUSG00000000001 -16.66124 3.4006651 938.6024 3.983953e-206
ENSMUSG00000000003 -26.70370 -2.7721353 0.0000 1.000000e+00
ENSMUSG00000000028 -19.59076 0.8862272 286.6896 2.617393e-64
ENSMUSG00000000037 -21.40817 -1.0586601 135.8326 2.170791e-31
ENSMUSG00000000049 -23.62721 -2.0704172 275.9815 5.640562e-62
34834 more rows ...
$comparison
[1] "Unbound"
$df.test
[1] 1 1 1 1 1
34834 more elements ...
