[R] numerical approximation of survival function
Carlo Giovanni Camarda
carlo-giovanni.camarda at ined.fr
Fri Nov 18 17:30:37 CET 2016
Dear R-users,
my question concerns numerical approximation and, somehow, survival
analysis.
Let’s assume we have a density function and we aim to numerically
compute the hazard, which is, in theory, the ratio between density and
survival.
In the following example, I take a pdf from a log-normal and proceed as
I would have known only this function. I numerically approximate the
survival function by splinefun() and integrate().Finally I compute the
hazard. It’s easy to see that (1) true hazard is different from the
numerically approximated one and (2) discrepancy depends upon the
approximation carried out to compute the survival. Is there any way to
overcome or attenuate this issue? Or do I just miss something obvious?
In my real study I do not deal with analytically defined pdf, so I would
like to be sure that approximation is done at its best.
Thanks in advance for your help.
Giancarlo Camarda
## x-values
delta <- 0.1
x <- seq(0, 200, by=delta)
m <- length(x)
## true pdf
fxT <- dlnorm(x, meanlog=3.5, sdlog=0.5)
fxT <- fxT/sum(fxT)
## true survival
SxT <- 1 - plnorm(x, meanlog=3.5, sdlog=0.5)
## true hazard
hxT <- fxT/SxT
## numerical approximation of the pdf
fx <- spline(x, fxT, n=m)$y
## numerical approximation of the survival
Sx <- numeric(m)
fxfun <- splinefun(x, fx/delta)
for(j in 1:m){
Sx[j] <- 1 - integrate(fxfun, x[1], x[j])$value
}
## hazard: pdf divided by survival
hx <- fx/Sx ## both numerical approx
hx1 <- fx/SxT ## only the survival numerically approx
hx2 <- fxT/Sx ## only the pdf numerically approx
## plotting hazards
plot(x, hxT)
lines(x, hx, col=2, lwd=2)
lines(x, hx1, col=3, lwd=2)
lines(x, hx2, col=4, lty=2, lwd=2)
## it seems that we have got something
## at the survival approx level:
plot(x, Sx/SxT)
[[alternative HTML version deleted]]
More information about the R-help
mailing list