thr3ads.net - R devel - [Rd] C function is wrong under Windows 7 [Oct 2011]

If this information is useful, please help other people find it:
Share via:

Evarist Planet

2011-Oct-24 13:04 UTC

[Rd] C function is wrong under Windows 7

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> 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

As you see es.c is correct. I checked it doing the same computation with R.
It also 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> 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

es.c is not correct under Windows. It also gives a different result when i
rerun the same function.

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;
}

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@irbbarcelona.org
http://www.irbbarcelona.org/bioinformatics

	[[alternative HTML version deleted]]

Thomas Lumley

2011-Oct-24 16:36 UTC

head link

[Rd] C function is wrong under Windows 7

On Mon, Oct 24, 2011 at 6: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.
The fact that you extract pointers to the contents of fchr and sign
before coercing them to REALSXP is deeply suspicious to me, though I
admit I don't see how it causes the problem, since you don't use the
new versions, you don't overwrite them, and the old versions should
still be protected as arguments to the function.
> 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;
> }
>
> 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
>
> ? ? ? ?[[alternative HTML version deleted]]
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>


-- 
Thomas Lumley
Professor of Biostatistics
University of Auckland

Martin Morgan

2011-Oct-24 16:38 UTC

head link

[Rd] C function is wrong under Windows 7

On 10/24/2011 06:04 AM, Evarist Planet wrote:> Dear mailing list,
>
> I have a C function that gives me a wrong result when I run it under
Windows
Hi Evarist --

It seems like this can be written reasonably efficiently in R?

getEs <-function(fchr, sign) {
     nfchr <- length(fchr)
     nsign <- length(sign)

     nr <- sum(abs(fchr[sign]))

     phit <- numeric(nfchr)
     phit[sign] <- abs(fchr[sign]) / nr
     phit <- cumsum(phit)

     pmiss <- numeric(nfchr)
     pmiss[-sign] <- 1 / (nfchr - nsign)
     pmiss <- cumsum(pmiss)

     phit - pmiss
}

es.c <- .Call('getEs',score,s,PACKAGE='phenoTest')
all.equal(es.c, getEs(score, s))

(for your C code, it would help to have a simple reproducible example 
that doesn't rely on phenoTest).

Martin
> 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
>> 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
>
> As you see es.c is correct. I checked it doing the same computation with R.
> It also 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
>> 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
>
> es.c is not correct under Windows. It also gives a different result when i
> rerun the same function.
>
> 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;
> }
>
> 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
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel

-- 
Computational Biology
Fred Hutchinson Cancer Research Center
1100 Fairview Ave. N. PO Box 19024 Seattle, WA 98109

Location: M1-B861
Telephone: 206 667-2793

Seemingly Similar Threads

Search for more seemingly similar threads

R devel - Oct 2011 - C function is wrong under Windows 7

[Rd] C function is wrong under Windows 7

[Rd] C function is wrong under Windows 7

[Rd] C function is wrong under Windows 7

Seemingly Similar Threads