Compare edgeR and DESeq2 RNA-seq results
Mark Ziemann
2020-10-05
library("edgeR")
library("DESeq2")
library("eulerr")x <- read.table("kevinwatt_YAP_gnames_filt.mx",row.names=1,header=TRUE)
xx <- x[which(rowMeans(x)>10),]
ss <- as.data.frame(colnames(xx))
ss$trt <- as.integer(grepl("YAP",ss[,1]))
rownames(ss) <- ss[,1]
ss[,1]=NULLdesign <- model.matrix(~ss$trt)
rownames(design) <- rownames(ss)
z <- DGEList(counts=xx)
z <- calcNormFactors(z)
z <- estimateDisp(z, design,robust=TRUE,prior.df=1)
fit <- glmFit(z, design)
lrt<-glmLRT(fit)
dge<-as.data.frame(topTags(lrt,n=Inf))
dge$dispersion<-lrt$dispersion
dge<-merge(dge,lrt$fitted.values,by='row.names')
rownames(dge)=dge$Row.names
dge$Row.names=NULL
dge<-dge[order(dge$PValue),]
head(dge,20)## logFC logCPM LR PValue
## ENSMUSG00000031880_Rrad 3.1972029 3.815439748 37.44393 9.407888e-10
## ENSMUSG00000053110_Yap1 -1.1042249 5.649211556 31.05114 2.513188e-08
## ENSMUSG00000097418_Mir155hg 2.1415047 5.482651770 25.21228 5.135404e-07
## ENSMUSG00000101567_Txn-ps1 0.8827607 1.797137631 22.07069 2.627924e-06
## ENSMUSG00000027861_Casq2 -0.4688354 5.953738012 19.46127 1.026601e-05
## ENSMUSG00000046993_Gm5637 0.9111893 1.663006881 19.20742 1.172568e-05
## ENSMUSG00000052837_Junb 0.7465992 4.061391526 17.65175 2.652711e-05
## ENSMUSG00000107265_Gm15469 -1.5043605 0.427842274 15.82002 6.966140e-05
## ENSMUSG00000035385_Ccl2 2.5369480 1.505463348 15.26810 9.327855e-05
## ENSMUSG00000041688_Amot -0.3896997 6.386455416 13.89128 1.936949e-04
## ENSMUSG00000030374_Strn4 0.3138776 4.658630776 13.75888 2.078364e-04
## ENSMUSG00000006221_Hspb7 0.7432434 9.015278247 13.72887 2.111834e-04
## ENSMUSG00000068614_Actc1 -1.1450147 8.177982870 13.51411 2.367760e-04
## ENSMUSG00000108322_Gm45123 -0.7047329 4.984158174 13.50436 2.380096e-04
## ENSMUSG00000091845_Gm4604 0.5071542 3.201228912 12.86429 3.349117e-04
## ENSMUSG00000035686_Thrsp -0.3222182 4.229477712 12.77015 3.521950e-04
## ENSMUSG00000044139_Prss53 -1.4032631 -0.008487417 12.76735 3.527225e-04
## ENSMUSG00000019066_Rab3d 0.6518751 2.756840015 12.51044 4.046848e-04
## ENSMUSG00000040950_Mgl2 0.5350511 2.482316308 11.13139 8.487904e-04
## ENSMUSG00000096544_Gm4617 1.0367471 0.650123912 11.04018 8.915808e-04
## FDR dispersion C1 C2
## ENSMUSG00000031880_Rrad 1.241371e-05 1.679345e-02 35.023453 32.603978
## ENSMUSG00000053110_Yap1 1.658076e-04 9.765625e-05 877.429785 816.815543
## ENSMUSG00000097418_Mir155hg 2.258722e-03 3.075308e-02 211.330529 196.731480
## ENSMUSG00000101567_Txn-ps1 8.668863e-03 9.372075e-01 29.867018 27.803757
## ENSMUSG00000027861_Casq2 2.578674e-02 7.757997e-03 922.176490 858.471074
## ENSMUSG00000046993_Gm5637 2.578674e-02 1.858088e-02 26.737477 24.890410
## ENSMUSG00000052837_Junb 5.000361e-02 2.672517e-02 158.580275 147.625298
## ENSMUSG00000107265_Gm15469 1.148978e-01 3.442563e-03 22.611343 21.049316
## ENSMUSG00000035385_Ccl2 1.367567e-01 2.391600e-02 9.947464 9.260278
## ENSMUSG00000041688_Amot 2.243240e-01 1.494407e-02 1216.769601 1132.713225
## ENSMUSG00000030374_Strn4 2.243240e-01 7.936348e-02 287.510199 267.648538
## ENSMUSG00000006221_Hspb7 2.243240e-01 1.064160e-01 4970.669785 4627.288025
## ENSMUSG00000068614_Actc1 2.243240e-01 9.073636e-03 5121.663333 4767.850699
## ENSMUSG00000108322_Gm45123 2.243240e-01 6.846747e-03 501.561507 466.912842
## ENSMUSG00000091845_Gm4604 2.737749e-01 2.502558e-02 95.811986 89.193142
## ENSMUSG00000035686_Thrsp 2.737749e-01 4.284278e-03 265.646919 247.295608
## ENSMUSG00000044139_Prss53 2.737749e-01 3.596123e-03 15.625829 14.546372
## ENSMUSG00000019066_Rab3d 2.966565e-01 1.088111e-03 65.888363 61.336690
## ENSMUSG00000040950_Mgl2 5.882204e-01 5.465013e-03 56.938685 53.005270
## ENSMUSG00000096544_Gm4617 5.882204e-01 5.762788e-03 11.798868 10.983784
## C3 YAP1 YAP2 YAP3
## ENSMUSG00000031880_Rrad 33.22027 299.008433 326.747348 314.563345
## ENSMUSG00000053110_Yap1 832.25538 378.599669 413.722236 398.295048
## ENSMUSG00000097418_Mir155hg 200.45019 865.508590 945.801539 910.533774
## ENSMUSG00000101567_Txn-ps1 28.32932 51.195975 55.945409 53.859274
## ENSMUSG00000027861_Casq2 874.69831 618.163719 675.510565 650.321615
## ENSMUSG00000046993_Gm5637 25.36090 46.756765 51.094374 49.189129
## ENSMUSG00000052837_Junb 150.41578 246.937609 269.845930 259.783710
## ENSMUSG00000107265_Gm15469 21.44720 7.316886 7.995671 7.697523
## ENSMUSG00000035385_Ccl2 9.43532 54.135962 59.158137 56.952204
## ENSMUSG00000041688_Amot 1154.12432 861.646395 941.581050 906.470661
## ENSMUSG00000030374_Strn4 272.70776 331.607541 362.370665 348.858312
## ENSMUSG00000006221_Hspb7 4714.75525 7719.770314 8435.930887 8121.365497
## ENSMUSG00000068614_Actc1 4857.97491 2148.623551 2347.950657 2260.398492
## ENSMUSG00000108322_Gm45123 475.73865 285.465933 311.948515 300.316342
## ENSMUSG00000091845_Gm4604 90.87911 126.388678 138.113714 132.963625
## ENSMUSG00000035686_Thrsp 251.97011 197.106472 215.391975 207.360276
## ENSMUSG00000044139_Prss53 14.82133 5.406397 5.907947 5.687647
## ENSMUSG00000019066_Rab3d 62.49611 96.116515 105.033212 101.116654
## ENSMUSG00000040950_Mgl2 54.00720 76.599490 83.705599 80.584321
## ENSMUSG00000096544_Gm4617 11.19140 22.584418 24.679567 23.759297dge_edger <- dge
sig <- subset(dge_edger,FDR<0.05)
dge_edger_up <- rownames(subset(sig,logFC>0))
dge_edger_dn <- rownames(subset(sig,logFC<0))
length(dge_edger_up)## [1] 4length(dge_edger_dn)## [1] 2dds <- DESeqDataSetFromMatrix(countData = xx , colData = ss, design = ~ trt )## the design formula contains one or more numeric variables with integer values,
## specifying a model with increasing fold change for higher values.
## did you mean for this to be a factor? if so, first convert
## this variable to a factor using the factor() functionres <- DESeq(dds)## estimating size factors## estimating dispersions## gene-wise dispersion estimates## mean-dispersion relationship## final dispersion estimates## fitting model and testingz<- results(res)
vsd <- vst(dds, blind=FALSE)
zz<-cbind(as.data.frame(z),assay(vsd))
dge<-as.data.frame(zz[order(zz$pvalue),])
head(dge,20)## baseMean log2FoldChange lfcSE stat
## ENSMUSG00000053110_Yap1 618.29299 -1.1010629 0.1671646 -6.586699
## ENSMUSG00000031880_Rrad 172.06141 3.2074253 0.5374085 5.968319
## ENSMUSG00000035385_Ccl2 32.82971 2.5541601 0.5072145 5.035661
## ENSMUSG00000006221_Hspb7 6396.77727 0.7467404 0.1534628 4.865939
## ENSMUSG00000097418_Mir155hg 550.84257 2.1441472 0.4506509 4.757890
## ENSMUSG00000031765_Mt1 305.41498 1.0540847 0.2540339 4.149386
## ENSMUSG00000019102_Aldh3a1 13.60618 3.9486760 0.9587922 4.118385
## ENSMUSG00000000031_H19 7488.07626 -0.6431462 0.1564750 -4.110218
## ENSMUSG00000108322_Gm45123 389.37924 -0.7022134 0.1748880 -4.015218
## ENSMUSG00000052837_Junb 204.46624 0.7506752 0.1910985 3.928211
## ENSMUSG00000027861_Casq2 764.33020 -0.4663011 0.1303572 -3.577102
## ENSMUSG00000001473_Tubb6 415.94170 0.6107095 0.1715146 3.560685
## ENSMUSG00000023067_Cdkn1a 530.33652 0.5867572 0.1718318 3.414718
## ENSMUSG00000019987_Arg1 15.85252 -2.5716010 0.7760504 -3.313703
## ENSMUSG00000068614_Actc1 3578.10497 -1.1428977 0.3457881 -3.305197
## ENSMUSG00000021268_Meg3 324.21314 -0.7263349 0.2210472 -3.285881
## ENSMUSG00000002985_Apoe 1498.99524 -0.5785025 0.1783652 -3.243359
## ENSMUSG00000035373_Ccl7 18.19394 2.1951511 0.6821945 3.217779
## ENSMUSG00000041688_Amot 1032.09822 -0.3860624 0.1238270 -3.117756
## ENSMUSG00000028011_Tdo2 14.79331 -2.8382828 0.9216300 -3.079634
## pvalue padj C1 C2
## ENSMUSG00000053110_Yap1 4.497128e-11 5.931262e-07 10.007683 9.912018
## ENSMUSG00000031880_Rrad 2.397104e-09 1.580771e-05 7.488156 7.096244
## ENSMUSG00000035385_Ccl2 4.762029e-07 2.093547e-03 6.981225 6.644697
## ENSMUSG00000006221_Hspb7 1.139146e-06 3.756050e-03 12.381124 12.225958
## ENSMUSG00000097418_Mir155hg 1.956274e-06 5.160259e-03 8.131147 8.546271
## ENSMUSG00000031765_Mt1 3.333689e-05 6.516788e-02 8.475911 8.571212
## ENSMUSG00000019102_Aldh3a1 3.815361e-05 6.516788e-02 6.589624 6.316077
## ENSMUSG00000000031_H19 3.952863e-05 6.516788e-02 13.438067 13.204287
## ENSMUSG00000108322_Gm45123 5.939082e-05 8.703394e-02 9.118699 9.392328
## ENSMUSG00000052837_Junb 8.557996e-05 1.128714e-01 8.224063 8.174546
## ENSMUSG00000027861_Casq2 3.474240e-04 4.065386e-01 9.924278 10.145364
## ENSMUSG00000001473_Tubb6 3.698887e-04 4.065386e-01 9.074874 8.808251
## ENSMUSG00000023067_Cdkn1a 6.384814e-04 6.477640e-01 9.183153 8.985862
## ENSMUSG00000019987_Arg1 9.206909e-04 8.345101e-01 6.882715 6.834009
## ENSMUSG00000068614_Actc1 9.490979e-04 8.345101e-01 12.218105 12.580638
## ENSMUSG00000021268_Meg3 1.016639e-03 8.380280e-01 9.484791 8.849026
## ENSMUSG00000002985_Apoe 1.181294e-03 9.164756e-01 10.877842 10.602191
## ENSMUSG00000035373_Ccl7 1.291873e-03 9.465843e-01 6.943790 6.480653
## ENSMUSG00000041688_Amot 1.822339e-03 9.999897e-01 10.463829 10.310809
## ENSMUSG00000028011_Tdo2 2.072552e-03 9.999897e-01 6.589624 6.683288
## C3 YAP1 YAP2 YAP3
## ENSMUSG00000053110_Yap1 9.959846 9.085615 8.781778 9.318499
## ENSMUSG00000031880_Rrad 7.053093 8.670771 8.483967 9.229548
## ENSMUSG00000035385_Ccl2 6.744658 7.159812 7.608564 7.641595
## ENSMUSG00000006221_Hspb7 12.192662 12.595351 13.034944 13.279479
## ENSMUSG00000097418_Mir155hg 8.559521 10.066668 10.493912 9.318499
## ENSMUSG00000031765_Mt1 8.131356 8.550862 9.404581 9.288702
## ENSMUSG00000019102_Aldh3a1 6.545764 6.989038 6.667314 7.476028
## ENSMUSG00000000031_H19 12.839366 12.742975 12.300236 12.578160
## ENSMUSG00000108322_Gm45123 9.381360 8.869400 8.913492 8.545178
## ENSMUSG00000052837_Junb 8.137006 8.643744 8.447239 8.815580
## ENSMUSG00000027861_Casq2 9.989051 9.542747 9.636404 9.709875
## ENSMUSG00000001473_Tubb6 8.772395 9.267580 9.202124 9.553867
## ENSMUSG00000023067_Cdkn1a 9.294325 9.520218 9.420877 9.882694
## ENSMUSG00000019987_Arg1 7.519124 6.600885 6.811590 6.476649
## ENSMUSG00000068614_Actc1 12.090555 10.732864 11.691244 11.071601
## ENSMUSG00000021268_Meg3 8.897645 8.381369 8.742614 8.641523
## ENSMUSG00000002985_Apoe 11.239795 10.338322 10.154812 10.676414
## ENSMUSG00000035373_Ccl7 6.597237 6.905800 7.361071 7.223916
## ENSMUSG00000041688_Amot 10.331778 10.087584 9.904777 10.100821
## ENSMUSG00000028011_Tdo2 7.621604 6.644767 6.667314 6.543043dge_deseq2 <- dge
sig <- subset(dge,padj<0.05)
dge_deseq2_up <- rownames(subset(sig,log2FoldChange>0))
dge_deseq2_dn <- rownames(subset(sig,log2FoldChange<0))
length(dge_deseq2_up)## [1] 4length(dge_deseq2_dn)## [1] 1v1 <- list("edgeR up"=dge_edger_up, "edgeR dn"=dge_edger_dn,
"DESeq2 up"=dge_deseq2_up,"DESeq2 dn"=dge_deseq2_dn)
plot(euler(v1),quantities = TRUE)
Let's look at some genes of interest.
# YAP1
# edgeR
t(dge_edger[grep("_Yap1",rownames(dge_edger)),])## ENSMUSG00000053110_Yap1
## logFC -1.104225e+00
## logCPM 5.649212e+00
## LR 3.105114e+01
## PValue 2.513188e-08
## FDR 1.658076e-04
## dispersion 9.765625e-05
## C1 8.774298e+02
## C2 8.168155e+02
## C3 8.322554e+02
## YAP1 3.785997e+02
## YAP2 4.137222e+02
## YAP3 3.982950e+02# DESeq2
t(dge_deseq2[grep("_Yap1",rownames(dge_deseq2)),])## ENSMUSG00000053110_Yap1
## baseMean 6.182930e+02
## log2FoldChange -1.101063e+00
## lfcSE 1.671646e-01
## stat -6.586699e+00
## pvalue 4.497128e-11
## padj 5.931262e-07
## C1 1.000768e+01
## C2 9.912018e+00
## C3 9.959846e+00
## YAP1 9.085615e+00
## YAP2 8.781778e+00
## YAP3 9.318499e+00# Amotl2
# edgeR
t(dge_edger[grep("Amotl2",rownames(dge_edger)),])## ENSMUSG00000032531_Amotl2
## logFC 0.14971276
## logCPM 3.90871010
## LR 1.24619392
## PValue 0.26428066
## FDR 1.00000000
## dispersion 0.07555595
## C1 181.02453045
## C2 168.51906848
## C3 171.70449700
## YAP1 186.32976108
## YAP2 203.61551221
## YAP3 196.02294255# DESeq2
t(dge_deseq2[grep("Amotl2",rownames(dge_deseq2)),])## ENSMUSG00000032531_Amotl2
## baseMean 183.6288612
## log2FoldChange 0.1538704
## lfcSE 0.1827634
## stat 0.8419107
## pvalue 0.3998380
## padj 0.9999897
## C1 8.3915279
## C2 8.2612137
## C3 8.2082033
## YAP1 8.4142979
## YAP2 8.4263773
## YAP3 8.3043599# Glul
# edgeR
t(dge_edger[grep("Glul",rownames(dge_edger)),])## ENSMUSG00000026473_Glul
## logFC 2.069183e-01
## logCPM 6.612465e+00
## LR 3.274722e+00
## PValue 7.035484e-02
## FDR 1.000000e+00
## dispersion 8.011969e-03
## C1 1.165194e+03
## C2 1.084701e+03
## C3 1.105204e+03
## YAP1 1.247785e+03
## YAP2 1.363542e+03
## YAP3 1.312697e+03# DESeq2
t(dge_deseq2[grep("Glul",rownames(dge_deseq2)),])## ENSMUSG00000026473_Glul
## baseMean 1.207812e+03
## log2FoldChange 2.097094e-01
## lfcSE 1.262018e-01
## stat 1.661699e+00
## pvalue 9.657309e-02
## padj 9.999897e-01
## C1 1.027081e+01
## C2 1.041576e+01
## C3 1.025589e+01
## YAP1 1.038405e+01
## YAP2 1.065161e+01
## YAP3 1.045704e+01# Dkk3
# edgeR
t(dge_edger[grep("Dkk3",rownames(dge_edger)),])## ENSMUSG00000030772_Dkk3
## logFC 0.12908672
## logCPM 3.19097498
## LR 0.06524524
## PValue 0.79838986
## FDR 1.00000000
## dispersion 0.01513808
## C1 110.15074339
## C2 102.54135515
## C3 104.47964118
## YAP1 111.77274506
## YAP2 122.14186614
## YAP3 117.58734759# DESeq2
t(dge_deseq2[grep("Dkk3",rownames(dge_deseq2)),])## ENSMUSG00000030772_Dkk3
## baseMean 110.9302950
## log2FoldChange 0.1323659
## lfcSE 0.3403270
## stat 0.3889373
## pvalue 0.6973225
## padj 0.9999897
## C1 8.1418227
## C2 7.5435385
## C3 7.9283002
## YAP1 8.2404625
## YAP2 7.5021224
## YAP3 8.0369196