scheike / timereg

timereg R package, survival analysis methods
29 stars 5 forks source link

Error msg in haplo.surv function #4

Closed shirleydaidai closed 5 years ago

shirleydaidai commented 6 years ago

Hello, I'm trying to use your R package "HaploSurvival" and "timereg" to perform haplotype-based survival analysis. When trying to repeat code from "haplo.surv.html" and performing haplo.surv function, errors occurred. Not sure if here's the right place report this, but I cannot find vignettes or support for the "HaploSurvival". Thanks Shirley

> library(HaploSurvival)
> library(timereg)
X<-matrix(rbinom(n,1,0.5),n,1); X<-cbind(1,X);
Zh<- apply((g[,c(1,3)]==c(0,0)),1,prod)+
> ### Simulating data
             x1=X[,1],x2=X[,2],z1=Z[,1],z2=Z[,2],zh=Zh)

### specifies Cox-Aalen model
> n<-400;
> haplo.types<-rbind(c(0,0),c(1,0),c(0,1),c(1,1))
> ### randomly selects haplotypes
> pairs<-matrix(sample(1:4,n*2,replace=TRUE),n,2)
>
> ### observed genotype, haplotype-pair 1 = g[c(1,3),]
> g<-matrix(haplo.types[pairs,],n,4)
>
> gs<-geno.setup(g);
>
> Z<-matrix(rnorm(n),n,1); Z<-cbind(Z,Z*0.5+1+rnorm(n))
> X<-matrix(rbinom(n,1,0.5),n,1); X<-cbind(1,X);
> Zh<- apply((g[,c(1,3)]==c(0,0)),1,prod)+
+      apply((g[,c(2,4)]==c(0,0)),1,prod)  ### counts 0,0
>
> ### simulates model with baselines depending on X and
> ### proportional in Z and Zh
> ### Zh observed only through genotype
>
> time<-exp(Z %*% c(-0.3,0.3) + Zh * 0.5 ) * rexp(n)/ X %*% c(0.2,0.5)
> status<-(time<3); time[status==0]<-3
> simdata=data.frame(time=time,status=status,
+              x1=X[,1],x2=X[,2],z1=Z[,1],z2=Z[,2],zh=Zh)
>
> ### specifies Cox-Aalen model
> designX<-function(x,z,h) { ### design for baseline, does not depend on haplotype
+ return(x)
+ }
> designZ<-function(x,z,h) { ### design for proportional part of model
+ h<-round(h);
+ zh<-(h[1]==0)+(h[2]==0) ### counts number of haplotype 0 = "0,0" see gs
+ y<-c(z,zh)
+ return(y)
+ }
> designZ(c(1,3),c(1,1),c(0,0)); designZ(c(1,3),c(1,1),c(1,0))
[1] 1 1 2
[1] 1 1 1
>
> hapfit<-haplo.freqs(g,geno.setup=gs)
> out<-haplo.surv(Surv(time,status)~x2+prop(z1)+prop(z2),data=simdata,
+ designfuncX=designX,designfuncZ=designZ,geno.setup=gs,
+ two.stage=1,alpha=hapfit$haplo.alpha)
Warning messages:
1: In haplo.surv(Surv(time, status) ~ x2 + prop(z1) + prop(z2), data = simdata,  :
  passing an object of type 'language' to .C (arg 56) is deprecated
2: In haplo.surv(Surv(time, status) ~ x2 + prop(z1) + prop(z2), data = simdata,  :
  passing an object of type 'language' to .C (arg 57) is deprecated
> ### MLE starting value alpha
> out
Error in if (cox.aalen.object$prop.odds == 0) p.o <- FALSE else p.o <- TRUE :
  argument is of length zero

> sessionInfo()
R version 3.4.3 (2017-11-30)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 16.04.3 LTS

Matrix products: default
BLAS: /media/storage1/software/install/R-3.4.3/lib/R/lib/libRblas.so
LAPACK: /media/storage1/software/install/R-3.4.3/lib/R/lib/libRlapack.so

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=en_US.UTF-8    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] HaploSurvival_1.0-4 timereg_1.9.2       survival_2.41-3

loaded via a namespace (and not attached):
[1] compiler_3.4.3    Matrix_1.2-12     tools_3.4.3       splines_3.4.3
[5] grid_3.4.3        lava_1.6.1        numDeriv_2016.8-1 lattice_0.20-35
scheike commented 6 years ago

In newer versions of R it not possible to use the design-functions in the .C call when going to C for the computations. So it should run on an older version of R.

Updating to newer R versions is quite a bit of work, but I will keep this in mind.

best


From: shirleydaidai [notifications@github.com] Sent: Monday, July 30, 2018 9:29 AM To: scheike/timereg Cc: Subscribed Subject: [scheike/timereg] Error msg in haplo.surv function (#4)

Hello, I'm trying to use your R package "HaploSurvival" and "timereg" to perform haplotype-based survival analysis. When trying to repeat code from "haplo.surv.html" and performing haplo.surv function, errors occurred. Not sure if here's the right place report this, but I cannot find vignettes or support for the "HaploSurvival". Thanks Shirley

library(HaploSurvival) library(timereg) X<-matrix(rbinom(n,1,0.5),n,1); X<-cbind(1,X); Zh<- apply((g[,c(1,3)]==c(0,0)),1,prod)+

Simulating data

x1=X[,1],x2=X[,2],z1=Z[,1],z2=Z[,2],zh=Zh)

specifies Cox-Aalen model

n<-400; haplo.types<-rbind(c(0,0),c(1,0),c(0,1),c(1,1))

randomly selects haplotypes

pairs<-matrix(sample(1:4,n*2,replace=TRUE),n,2)

observed genotype, haplotype-pair 1 = g[c(1,3),]

g<-matrix(haplo.types[pairs,],n,4)

gs<-geno.setup(g);

Z<-matrix(rnorm(n),n,1); Z<-cbind(Z,Z*0.5+1+rnorm(n)) X<-matrix(rbinom(n,1,0.5),n,1); X<-cbind(1,X); Zh<- apply((g[,c(1,3)]==c(0,0)),1,prod)+

  • apply((g[,c(2,4)]==c(0,0)),1,prod) ### counts 0,0

simulates model with baselines depending on X and

proportional in Z and Zh

Zh observed only through genotype

time<-exp(Z %% c(-0.3,0.3) + Zh 0.5 ) rexp(n)/ X %% c(0.2,0.5) status<-(time<3); time[status==0]<-3 simdata=data.frame(time=time,status=status,

  • x1=X[,1],x2=X[,2],z1=Z[,1],z2=Z[,2],zh=Zh)

specifies Cox-Aalen model

designX<-function(x,z,h) { ### design for baseline, does not depend on haplotype

  • return(x)
  • } designZ<-function(x,z,h) { ### design for proportional part of model
  • h<-round(h);
  • zh<-(h[1]==0)+(h[2]==0) ### counts number of haplotype 0 = "0,0" see gs
  • y<-c(z,zh)
  • return(y)
  • } designZ(c(1,3),c(1,1),c(0,0)); designZ(c(1,3),c(1,1),c(1,0)) [1] 1 1 2 [1] 1 1 1

hapfit<-haplo.freqs(g,geno.setup=gs) out<-haplo.surv(Surv(time,status)~x2+prop(z1)+prop(z2),data=simdata,

  • designfuncX=designX,designfuncZ=designZ,geno.setup=gs,
  • two.stage=1,alpha=hapfit$haplo.alpha) Warning messages: 1: In haplo.surv(Surv(time, status) ~ x2 + prop(z1) + prop(z2), data = simdata, : passing an object of type 'language' to .C (arg 56) is deprecated 2: In haplo.surv(Surv(time, status) ~ x2 + prop(z1) + prop(z2), data = simdata, : passing an object of type 'language' to .C (arg 57) is deprecated

    MLE starting value alpha

    out Error in if (cox.aalen.object$prop.odds == 0) p.o <- FALSE else p.o <- TRUE : argument is of length zero

sessionInfo() R version 3.4.3 (2017-11-30) Platform: x86_64-pc-linux-gnu (64-bit) Running under: Ubuntu 16.04.3 LTS

Matrix products: default BLAS: /media/storage1/software/install/R-3.4.3/lib/R/lib/libRblas.so LAPACK: /media/storage1/software/install/R-3.4.3/lib/R/lib/libRlapack.so

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=en_US.UTF-8 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] HaploSurvival_1.0-4 timereg_1.9.2 survival_2.41-3

loaded via a namespace (and not attached): [1] compiler_3.4.3 Matrix_1.2-12 tools_3.4.3 splines_3.4.3 [5] grid_3.4.3 lava_1.6.1 numDeriv_2016.8-1 lattice_0.20-35

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHubhttps://github.com/scheike/timereg/issues/4, or mute the threadhttps://github.com/notifications/unsubscribe-auth/AJvzwGC_kcgM-TuliazL9A1OjMZxFpWRks5uLrX1gaJpZM4Vl6Yx.