Entering edit mode
Cécile Laurent
▴
30
@cecile-laurent-2766
Last seen 10.3 years ago
Hi,
I would like to know why when I run SAM analysis with SamR and
Siggenes
packages, I don't have the same results for the small pvalues.
I test here on golub data, and I don't have the same differential
genes
number when I observe the delta tables (delta.table for samr and
sam.out at mat.fdr for siggenes) for a FDR 0% (250 are differential
with
samr package and 9 with siggenes package).
On my personnal data, I observe the same difference in differential
gene
number at FDR 5%.
I also observe the same s0, but not the same pi0, in the 2 sam
objects.
I think, in my code, the parameters are the same... So, I don't know
which package used...
Is it something wrong???
Can you explain me these differences in the results..
Thanks,
C?cile
Here is my sessionInfo()
R version 2.6.0 (2007-10-03)
i486-pc-linux-gnu
locale:
LC_CTYPE=fr_FR.UTF-8;LC_NUMERIC=C;LC_TIME=fr_FR.UTF-8;LC_COLLATE=fr_FR
.UTF-8;LC_MONETARY=fr_FR.UTF-8;LC_MESSAGES=fr_FR.UTF-8;LC_PAPER=fr_FR.
UTF-8;LC_NAME=C;LC_ADDRESS=C;LC_TELEPHONE=C;LC_MEASUREMENT=fr_FR.UTF-8
;LC_IDENTIFICATION=C
attached base packages:
[1] splines tools stats graphics grDevices utils
datasets
[8] methods base
other attached packages:
[1] siggenes_1.12.0 samr_1.25 impute_1.0-5 multtest_1.18.0
[5] survival_2.33 Biobase_1.16.0
loaded via a namespace (and not attached):
[1] rcompgen_0.1-15
### My code :
library(multtest)
library(samr)
library(siggenes)
data(golub)
golubCl=golub.cl
golubCl[whichgolub.cl==1)]=2
golubCl[whichgolub.cl==0)]=1
## SamR package
samrData=list(x=golub, y=golubCl, geneid=golub.gnames[,3],
genesnames=golub.gnames[,3], logged2=T)
samr.obj=samr(samrData, resp.type="Two class unpaired",
testStatistic="standard", nperms=1000, random.seed=123)
delta.table=samr.compute.delta.table(samr.obj, dels=seq(0.1,5,0.05))
pv.samr=samr.pvalues.from.perms(samr.obj$tt, samr.obj$ttstar)
## Siggenes package
sam.out=sam(golub, golub.cl, rand=123,B=1000,
gene.names=golub.gnames[,3], method="d.stat", var.equal=T, s0=NA,
include.zero=F, delta=seq(0.1,5,0.05) )
## plot to observe the values (the version with -log() to see the
difference in the small pvalues)
#plot(sam.out at d[order(sam.out at
d)],samr.obj$tt[order(samr.obj$tt)])
#plot(-log(sam.out at d[order(sam.out at
d)]),-log(samr.obj$tt[order(samr.obj$tt)]))
plot(sam.out at p.value[order(sam.out at
p.value)],pv.samr[order(pv.samr)])
plot(-log(sam.out at p.value[order(sam.out at
p.value)]),-log(pv.samr[order(pv.samr)]))