[Bioc-devel] C function is wrong under Windows 7.
I think I'd pay more attention to the architecture (32 vs. 64 bit). These days I find myself doing all of my C(++) <--> R bridge-work via Rcpp, so I'm a bit rusty w/ the "normal" conventions here, but out of curiosity: for the vars that essentially never "work with R", like pmiss for example, why do you go through the whole PROTECT/UNPROTECT allocVector rigmarole with it instead of just creating a double array of the appropriate length right off of the bat (then freeing it, obviously)? -steve On Mon, Oct 24, 2011 at 9:40 AM, Evarist Planet
<evarist.planet at irbbarcelona.org> wrote:
I run the same code on Windows 32-bit R version 2.13.x and got a similar result. I do not remember which exact version under 2.13 I used. I'll check it and try to run the code on the same version of R. Thanks for the quick reply, Evarist On Mon, Oct 24, 2011 at 3:20 PM, Steve Lianoglou <mailinglist.honeypot at gmail.com> wrote:
Hi, I didn't look at your code, but: Session info for you linux machine shows that it is 64-bit and R version 2.13.0 You window box is running in 32-bit and is R version 2.14.0-alpha Maybe you can can try to control one or both of these factors and test again? -steve On Mon, Oct 24, 2011 at 8:33 AM, Evarist Planet <evarist.planet at irbbarcelona.org> wrote:
Ok. This is the result of sessionInfo() under linux:
sessionInfo()
R version 2.13.0 (2011-04-13) Platform: x86_64-redhat-linux-gnu (64-bit) locale: ?[1] LC_CTYPE=en_US.UTF-8 ? ? ? LC_NUMERIC=C ?[3] LC_TIME=en_US.UTF-8 ? ? ? ?LC_COLLATE=en_US.UTF-8 ?[5] LC_MONETARY=C ? ? ? ? ? ? ?LC_MESSAGES=en_US.UTF-8 ?[7] LC_PAPER=en_US.UTF-8 ? ? ? LC_NAME=C ?[9] LC_ADDRESS=C ? ? ? ? ? ? ? LC_TELEPHONE=C [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] stats ? ? graphics ?grDevices utils ? ? datasets ?methods ? base other attached packages: [1] phenoTest_1.1.1 ? ? ?RSQLite_0.8-4 ? ? ? ?DBI_0.2-5 [4] Heatplus_1.22.0 ? ? ?annotate_1.30.1 ? ? ?AnnotationDbi_1.14.1 [7] Biobase_2.12.2 loaded via a namespace (and not attached): ?[1] affyio_1.20.0 ? ? ? biomaRt_2.8.1 ? ? ? Biostrings_2.21.6 ?[4] Category_2.18.0 ? ? cluster_1.13.3 ? ? ?gdata_2.7.1 ?[7] genefilter_1.34.0 ? gplots_2.10.1 ? ? ? graph_1.30.0 [10] grid_2.13.0 ? ? ? ? GSEABase_1.14.0 ? ? gtools_2.6.2 [13] hgu133a.db_2.5.0 ? ?Hmisc_3.8-3 ? ? ? ? hopach_2.12.0 [16] IRanges_1.11.11 ? ? lattice_0.19-23 ? ? limma_3.8.3 [19] Matrix_0.9996875-3 ?mgcv_1.7-5 ? ? ? ? ?nlme_3.1-100 [22] oligoClasses_1.14.0 RBGL_1.22.0 ? ? ? ? RCurl_1.6-10 [25] SNPchip_1.16.0 ? ? ?splines_2.13.0 ? ? ?survival_2.36-5 [28] tools_2.13.0 ? ? ? ?XML_3.4-0 ? ? ? ? ? xtable_1.5-6 And under Windows 7:
sessionInfo()
R version 2.14.0 alpha (2011-10-13 r57240)
Platform: i386-pc-mingw32/i386 (32-bit)
locale:
[1] LC_COLLATE=Spanish_Spain.1252 ?LC_CTYPE=Spanish_Spain.1252
[3] LC_MONETARY=Spanish_Spain.1252 LC_NUMERIC=C
[5] LC_TIME=Spanish_Spain.1252
attached base packages:
[1] stats ? ? graphics ?grDevices utils ? ? datasets ?methods ? base
other attached packages:
[1] phenoTest_1.1.1 ? ? ? RSQLite_0.10.0 ? ? ? ?DBI_0.2-5
[4] Heatplus_1.99.0 ? ? ? annotate_1.31.1 ? ? ? AnnotationDbi_1.15.36
[7] Biobase_2.13.10
loaded via a namespace (and not attached):
?[1] affyio_1.21.2 ? ? ? ?biomaRt_2.9.3 ? ? ? ?Biostrings_2.21.11
?[4] Category_2.19.1 ? ? ?cluster_1.14.0 ? ? ? gdata_2.8.2
?[7] genefilter_1.35.0 ? ?gplots_2.10.1 ? ? ? ?graph_1.31.2
[10] grid_2.14.0 ? ? ? ? ?GSEABase_1.15.0 ? ? ?gtools_2.6.2
[13] hgu133a.db_2.6.3 ? ? Hmisc_3.8-3 ? ? ? ? ?hopach_2.13.1
[16] IRanges_1.11.31 ? ? ?lattice_0.19-33 ? ? ?limma_3.9.21
[19] Matrix_1.0-0 ? ? ? ? mgcv_1.7-8 ? ? ? ? ? nlme_3.1-102
[22] oligoClasses_1.15.56 RBGL_1.29.0 ? ? ? ? ?RCurl_1.6-10.1
[25] SNPchip_1.17.0 ? ? ? splines_2.14.0 ? ? ? survival_2.36-10
[28] tools_2.14.0 ? ? ? ? XML_3.4-2.2 ? ? ? ? ?xtable_1.6-0
[31] zlibbioc_0.1.8
This is the C code:
#include "getEs.h"
#include <R.h>
#include <Rinternals.h>
double absolute(double x)
{
?if (x<0)
?return (-x);
?else
?return (x);
}
void cumsum(double *x, int len)
{
?int i;
?for (i = 1; i < len; ++i) {
? ?*(x + i) = *(x + i) + *(x + i -1);
?}
}
double getNr(double *fchr, int *sign, int signLen)
{
?int i;
?double nr;
?nr = 0.0;
?for (i = 0; i < signLen; ++i) {
? ?nr = absolute(fchr[sign[i] -1]) + nr;
? ?}
?return nr;
}
void getPhit(double *fchr, int *sign, int signLen, double nr, double
*phit)
{
?int i;
?for (i = 0; i < signLen; ++i) {
? ?*(phit + sign[i]-1) = absolute(*(fchr + sign[i]-1)) / nr;
? ?}
}
void getPmiss(int *sign, int fchrLen, int signLen, double *pmiss)
{
?int i;
?double tmp = 1.0 / (fchrLen-signLen);
?for (i = 0; i < fchrLen; ++i) {
? ?*(pmiss + i) = tmp;
? ?}
?for (i = 0; i < signLen; ++i) {
? ?*(pmiss + sign[i]-1) = 0;
? ?}
}
SEXP getEs(SEXP fchr, SEXP sign)
{
?int i, nfchr, nsign;
?double *rfchr = REAL(fchr), *res;
?int *rsign = INTEGER(sign);
?PROTECT(fchr = coerceVector(fchr, REALSXP));
?PROTECT(sign = coerceVector(sign, REALSXP));
?nfchr = length(fchr);
?nsign = length(sign);
?SEXP es;
?PROTECT(es = allocVector(REALSXP, nfchr));
?res = REAL(es);
?double nr = getNr(rfchr, rsign, nsign);
?SEXP phit;
?PROTECT(phit = allocVector(REALSXP, nfchr));
?double *rphit;
?rphit = REAL(phit);
?for(i = 0; i < nfchr; i++) rphit[i] = 0.0;
?getPhit(rfchr, rsign, nsign, nr, rphit);
?cumsum(rphit, nfchr);
?SEXP pmiss;
?PROTECT(pmiss = allocVector(REALSXP, nfchr));
?double *rpmiss;
?rpmiss = REAL(pmiss);
?for(i = 0; i < nfchr; i++) rpmiss[i] = 0.0;
?getPmiss(rsign, nfchr, nsign, rpmiss);
?cumsum(rpmiss, nfchr);
?for (i = 0; i < nfchr; ++i) {
? ?res[i] = rphit[i] - rpmiss[i];
? ?}
?UNPROTECT(5);
?return es;
}
Many thanks,
Evarist
On Mon, Oct 24, 2011 at 2:22 PM, Vincent Carey
<stvjc at channing.harvard.edu>wrote:
your attachment is automatically scrubbed so we cannot see it. ?this is not really a bioconductor question; you might consider r-devel. ?you would also do well to include the result of sessionInfo(). On Mon, Oct 24, 2011 at 8:04 AM, Evarist Planet <evarist.planet at irbbarcelona.org> wrote:
Dear mailing list, I have a C function that gives me a wrong result when I run it under
Windows
7. This is the code under Linux (RHEL5):
library(phenoTest) data(epheno) sign <- sample(featureNames(epheno))[1:20] score <- getFc(epheno)[,1] head(score)
1007_s_at ? 1053_at ? ?117_at ? ?121_at 1255_g_at ? 1294_at -1.183019 ?1.113544 ?1.186186 -1.034779 -1.044456 -1.023471
s <- which(names(score) %in% sign)
es.c <- .Call('getEs',score,s,PACKAGE='phenoTest')
head(es.c)
[1] -0.001020408 -0.002040816 -0.003061224 -0.004081633 -0.005102041 [6] -0.006122449
es.c <- .Call('getEs',score,s,PACKAGE='phenoTest')
head(es.c)
[1] -0.001020408 -0.002040816 -0.003061224 -0.004081633 -0.005102041 [6] -0.006122449 As you see es.c is correct. I checked it doing the same computation with R. It runs without problems under Mac. I run valgrind on the same piece
of
code and got no errors. This is the same piece of code under Windows 7:
library(phenoTest) data(epheno) sign <- sample(featureNames(epheno))[1:20] score <- getFc(epheno)[,1] head(score)
1007_s_at ? 1053_at ? ?117_at ? ?121_at 1255_g_at ? 1294_at -1.183019 ?1.113544 ?1.186186 -1.034779 -1.044456 -1.023471
s <- which(names(score) %in% sign)
es.c <- .Call('getEs',score,s,PACKAGE='phenoTest')
head(es.c)
[1] 1.447208e+215 1.447208e+215 1.447208e+215 1.447208e+215 1.447208e+215 1.447208e+215
es.c <- .Call('getEs',score,s,PACKAGE='phenoTest')
head(es.c)
[1] 3.176615e+170 3.176615e+170 3.176615e+170 3.176615e+170 3.176615e+170 3.176615e+170 es.c is not correct under Windows. It also gives a different result when
i
rerun the same function. I attached the C code. Could you please help me to find out what I am doing wrong? Many thanks in advance, -- Evarist Planet Research officer, Bioinformatics and Biostatistics unit IRB Barcelona Tel (+34) 93 402 0553 Fax (+34) 93 402 0257 evarist.planet at irbbarcelona.org http://www.irbbarcelona.org/bioinformatics
_______________________________________________ Bioc-devel at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/bioc-devel
-- Evarist Planet Research officer, Bioinformatics and Biostatistics unit IRB Barcelona Tel (+34) 93 402 0553 Fax (+34) 93 402 0257 evarist.planet at irbbarcelona.org http://www.irbbarcelona.org/bioinformatics ? ? ? ?[[alternative HTML version deleted]]
_______________________________________________ Bioc-devel at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/bioc-devel
-- Steve Lianoglou Graduate Student: Computational Systems Biology ?| Memorial Sloan-Kettering Cancer Center ?| Weill Medical College of Cornell University Contact Info: http://cbio.mskcc.org/~lianos/contact
-- Evarist Planet Research officer, Bioinformatics and Biostatistics unit IRB Barcelona Tel (+34) 93 402 0553 Fax (+34) 93 402 0257 evarist.planet at irbbarcelona.org http://www.irbbarcelona.org/bioinformatics
Steve Lianoglou Graduate Student: Computational Systems Biology ?| Memorial Sloan-Kettering Cancer Center ?| Weill Medical College of Cornell University Contact Info: http://cbio.mskcc.org/~lianos/contact