R help - Apr 2011 - help on pspline in coxph

Hi there,

I have a question on how to extract the linear term in the penalized 
spline in coxph. Here is a sample code:

n=100

set.seed(1)

x=runif(100)

f1 = cos(2*pi*x)

hazard = exp(f1)

T = 0

for (i in 1:100) {

       T[i] = rexp(1,hazard[i])

}

     C = runif(n)*4

     cen = T<=C

     y = T*(cen) + C*(1-cen)

     data.tr=cbind(y,cen,x)

    fit=coxph(Surv(data.tr[,1], data.tr[,2])~pspline(data.tr[,3]) )

If I use summary(fit), it will show the following results:

summary(fit)
Call:
coxph(formula = Surv(data.tr[, 1], data.tr[, 2]) ~ pspline(data.tr[,
     3]))

   n= 100
                           coef  se(coef) se2   Chisq DF   p
pspline(data.tr[, 3]), li 0.495 0.437    0.437  1.28 1.00 2.6e-01
pspline(data.tr[, 3]), no                      43.79 3.08 1.9e-09

                    exp(coef) exp(-coef) lower .95 upper .95
ps(data.tr[, 3])2     0.4404      2.270   0.08164     2.376
ps(data.tr[, 3])3     0.2065      4.842   0.01701     2.507
ps(data.tr[, 3])4     0.0951     10.512   0.00695     1.302
ps(data.tr[, 3])5     0.0493     20.274   0.00387     0.628
ps(data.tr[, 3])6     0.0280     35.741   0.00230     0.340
ps(data.tr[, 3])7     0.0192     52.068   0.00156     0.237
ps(data.tr[, 3])8     0.0219     45.605   0.00178     0.271
ps(data.tr[, 3])9     0.0473     21.156   0.00399     0.561
ps(data.tr[, 3])10    0.1432      6.983   0.01250     1.640
ps(data.tr[, 3])11    0.3936      2.541   0.03436     4.509
ps(data.tr[, 3])12    0.9449      1.058   0.06885    12.969
ps(data.tr[, 3])13    2.2406      0.446   0.07643    65.683

Iterations: 3 outer, 9 Newton-Raphson
      Theta= 0.697
Degrees of freedom for terms= 4.1
Rsquare= 0.385   (max possible= 0.994 )
Likelihood ratio test= 48.6  on 4.08 df,   p=7.74e-10
Wald test            = 45.1  on 4.08 df,   p=4.29e-09

My question is how to extract the linear coefficient (0.495) in 
pspline(data.tr[, 3]). I tried coef(fit) but it fails to display this 
term. Your help is greatly appreciated!

Lei Liu
Associate Professor
Division of Biostatistics and Epidemiology
Department of Public Health Sciences
University of Virginia School of Medicine

http://people.virginia.edu/~ll9f/

> Hi there,
> I have a question on how to extract the linear term in the penalized 
> spline in coxph. Here is a sample code:
Here is a slightly more general example: 
   fit <- coxph(Surv(time, status) ~ ph.ecog + pspline(age) + sex, lung)
(The lung data set is part of the survival package.)
There are 14 coefficients, of which 2-13 are the for the spline.  

The slope logic is subtle. Luckily I left comments on this in the source
code for pspline, because I remember it taking a long time for me to
figure it out.  
      slope = (X'VX)^{-1} X'V c
A weighted regression on the coefficients where c= spline coefficients
fit$coef[2:13] in this case, V= fit$var[2:13,2:13], and X is a two
column matrix with column 1=1 (intercept) and column 2 = centers of the
pspline basis functions.  These centers are saved away in the print
function:
	get("cbase", environment(fit$printfun[[1]]))
If there were multiple psline terms in the fit use [[2]], [[3]] etc as
appropriate.  No, you won't find this without reading the source code.

  Note that this is the linear slope from a fit with only a linear term,
i.e., what would have resulted from a fit with "age" instead of
pspline(age).  I thought it quite interesting that this quantity could
be derived at all.  I leave the details as a homework exercise for the
reader.

Terry Therneau