On Nov 25, 2010, at 7:27 AM, Ben Rhelp wrote:
I manage to achieve similar results with a Cox model as follows but
I don't
really understand why we have to take the inverse of the linear
prediction with
the Cox model
Different parameterization. You can find expanded answer(s) in the
archives and in the documentation of survreg.distributions.
and why we do not need to divide by the number of days in the year
anymore?
Here I'm guessing (since you don't offer enough evidence to confirm)
that the difference is in the time scales used in your Aidsp$survtime
versus some other example to which you are comparing .
> data(Aidsp)
Warning message:
In data(Aidsp) : data set 'Aidsp' not found
Am I getting a similar result out of pure luck?
thanks in advance,
Ben
# MASS example with the proportional hazard model
par(mfrow = c(1, 2));
(aids.ps <- survreg(Surv(survtime + 0.9, status) ~ state + T.categ +
pspline(age, df=6), data = Aidsp))
zz <- predict(aids.ps, data.frame(state = factor(rep("NSW", 83),
levels =
levels(Aidsp$state)),
T.categ = factor(rep("hs", 83), levels = levels(Aidsp$T.categ)),
age =
0:82), se = T, type = "linear")
plot(0:82, exp(zz$fit)/365.25, type = "l", ylim = c(0, 2), xlab =
"age", ylab =
"expected lifetime (years)")
lines(0:82, exp(zz$fit+1.96*zz$se.fit)/365.25, lty = 3, col = 2)
lines(0:82, exp(zz$fit-1.96*zz$se.fit)/365.25, lty = 3, col = 2)
rug(Aidsp$age + runif(length(Aidsp$age), -0.5, 0.5), ticksize = 0.015)
# same example but with a Cox model instead
(aids.pscp <- coxph(Surv(survtime + 0.9, status) ~ state + T.categ +
pspline(age, df=6), data = Aidsp))
zzcp <- predict(aids.pscp, data.frame(state = factor(rep("NSW", 83),
levels =
levels(Aidsp$state)),
T.categ = factor(rep("hs", 83), levels = levels(Aidsp$T.categ)),
age =
0:82), se = T, type = "lp")
plot(0:82, 1/exp(zzcp$fit), type = "l", ylim = c(0, 2), xlab =
"age", ylab =
"expected lifetime (years)")
lines(0:82, 1/exp(zzcp$fit+1.96*zzcp$se.fit), lty = 3, col = 2)
lines(0:82, 1/exp(zzcp$fit-1.96*zzcp$se.fit), lty = 3, col = 2)
rug(Aidsp$age + runif(length(Aidsp$age), -0.5, 0.5), ticksize = 0.015)
David Winsemius, MD
West Hartford, CT
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
