Dear users. Hello? I have a question about ComBat results using SVA and bladderbatch package. The followings are the codes and the results I got. library(sva) library(pamr) library(limma) library(bladderbatch) data(bladderdata) pheno = pData(bladderEset) edata = exprs(bladderEset) mod = model.matrix(~as.factor(cancer), data=pheno) mod0 = model.matrix(~1,data=pheno) pValues = f.pvalue(edata,mod,mod0) qValues = p.adjust(pValues,method="BH") sum(qValues<=0.05) [1] 15193 sum(qValues<=0.05)/length(qValues) [1] 0.6818202 Note that nearly 70% of the genes are strongly differentially expressed at an FDR of less than 5% between groups. This number seems artificially high, even for a strong phenotype like cancer. This is the point that the authors emphasized in the bladderbatchTutorial. For ComBat, batch = pheno$batch mod = model.matrix(~as.factor(cancer), data=pheno) combat_edata = ComBat(data=edata, batch=batch, mod=mod, numCovs=NULL, par.prior=TRUE, prior.plots=TRUE) pValuesComBat = f.pvalue(combat_edata,mod,mod0) qValuesComBat = p.adjust(pValuesComBat,method="BH") sum(qValuesComBat<=0.05) [1] 18300 sum(qValuesComBat<=0.05)/length(qValuesComBat) [1] 0.8212539 After ComBat adjustment, 80% of the genes are differentially expressed at an FDR of less than 5% between groups. (The authors did not provide their results.) Is this result correct? -- output of sessionInfo(): R version 2.14.2 (2012-02-29) Platform: i386-pc-mingw32/i386 (32-bit) locale: [1] LC_COLLATE=Korean_Korea.949 LC_CTYPE=Korean_Korea.949 [3] LC_MONETARY=Korean_Korea.949 LC_NUMERIC=C [5] LC_TIME=Korean_Korea.949 attached base packages: [1] splines stats graphics grDevices utils datasets methods [8] base other attached packages: [1] bladderbatch_1.0.2 Biobase_2.14.0 limma_3.10.3 [4] pamr_1.54 survival_2.37-4 cluster_1.14.3 [7] sva_3.0.2 mgcv_1.7-22 corpcor_1.6.4 loaded via a namespace (and not attached): [1] grid_2.14.2 lattice_0.20-10 Matrix_1.0-5 nlme_3.1-108 -- Sent via the guest posting facility at bioconductor.org.