----- Forwarded Message -----
From: Jichun Chan <jichun@y7mail.com>
To: David Winsemius <dwinsemius@comcast.net>
Sent: Friday, 14 December 2012 5:30 PM
Subject: Re: [R] simulate time data
Hi David,
Thank you for replying my email. I have to apology that I haven't been very
clear with my question.
I am actually trying to simulate a set of data with two correlated variables
(e.g. x1 and x2) , one of which (says x1) is associated with the outcome with a
predefined OR. Then, I want to calculate the HR of the other variable in
relation to the outcome using Cox's regression. But I don't know how to
simulate the "time" data as it is not normally distributed.
I tried to treat it like normally distributed and put it the covariance matrix,
it somehow produced the HR that is closely resembled that of real data. Am I
doing the right thing. The following is the Rcommond I used:
> set.seed(12345)
> for (i in 1:1000) {
> sim <- rmvnorm( n = 1000, mean = c(0.81,26.5,88.07), sigma =
matrix(c(0.019, 0.41*sqrt(0.48), 1.27,
0.41*sqrt(0.48), 25.48, -23.82, 1.27, -23.82, 5297.64),3,3))
> x1 <- sim[,1]
> x2 <- sim[,2]
> time<- sim[,3]
> L <- log(1.89) * (x1/0.14) #0.14 is the SD and 1.89 is the preset OR
> event <- ifelse(runif(1000) < plogis(L), 1, 0)
> sdata<- data.frame(x1, x2, event)
> m1<- coxph(Surv(time, event) ~ x2, data=sdata)
I also tried to simulate the time using Weibull distribution from the R help as
follows, but I don't quite understand what they are doing. And I don't
know how to combine it to the dataset so that the time variable is related to
the outcome. I hope you know what I mean.
> n = 1000
> beta1 =0.18 ; beta2 = -0.65
> lambdaT = .007 # baseline hazard
> lambdaC = .003 # hazard of censoring
> T = rweibull(n, shape=1, scale=lambdaT*exp(-beta1*bmi-beta2*fnbmd))
> C = rweibull(n, shape=1, scale=lambdaC) #censoring time
> time = pmin(T,C) #observed time is minimum of censored and true
> event = time==T # set to 1 if event is observed
Would you please kindly help me with this and forgive me if I ask silly question
as I am still new to statistics.
Regards
Scomet
From: David Winsemius <dwinsemius@comcast.net>
To: Jichun Chan <jichun@y7mail.com>
Cc: "r-help@r-project.org" <r-help@r-project.org>
Sent: Thursday, 13 December 2012 3:10 PM
Subject: Re: [R] simulate time data
On Dec 12, 2012, at 6:34 PM, Jichun Chan wrote:
> Hi,
> Does anyone know how to write a command to generate time-to-event data for
Cox's regression?
>
There are examples in the rms package for the cph function. I also see example
in the help pages for the survival package.
--David Winsemius, MD
Alameda, CA, USA
[[alternative HTML version deleted]]