[R] how to check linearity in Cox regression

2008-11-26 Thread Terry Therneau
 On examining non-linearity of Cox coefficients with penalized splines -  I
 have not been able to dig up a completely clear description of the test
 performed in R or S-plus.

 One iron clad way to test is to fit a model that has the variable of 
interest 
x as a linear term, then a second model with splines, and do a likelihood 
ratio test with 2*(difference in log-likelihood) on (difference in df) degrees 
of freedom.  With a penalized model this test is conservative: the chi-square 
is 
not quite the right distribution, the true dist has the same mean but smaller 
variance.

 The pspline function uses an evenly spaced set of symmetric basis functions.  
A 
neat consequence of this is that the Wald test for linear vs 'more general' is 
a 
test that the coefficients of the spline terms fall in a linear series.  That 
is, a linear trend test on the coefficients.  This is what coxph does.  As with 
the LR test, the chi-square dist is conservative.  I have not worked at putting 
in the more correct distribution.  See Eilers and Marx, Statistical Science 
1986.
 
  And what is the null for the non-linear test?

The linear test is is a linear better than nothing, the non-linear one is a 
sequential test is the non-linear better than the linear.  The second test of 
course depends on the total number of df you allowed for the pspline fit.  As a 
silly example adding + pspline(x, df=200) would likely show that the 
nonlinear 
term was not a significant addition, i.e., not worth 199 more degrees of 
freedom.

Terry Therneau

__
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.


Re: [R] how to check linearity in Cox regression

2008-11-25 Thread Michael Margolis
On examining non-linearity of Cox coefficients with penalized splines -  I
have not been able to dig up a completely clear description of the test
performed in R or S-plus.

From the Therneau and Grambsch book (2000 - page 126) I gather that the test
reported for linear has as its null hypothesis that the spline coefficient
is the same at the center of basis. Thus, in the example quoted below (from
the thread) one might say treating age as linear is not too bad, to the
extent that failure to reject a null can be so interpreted.

Is that right? And what is the null for the nonlinear test?

--Michael Margolis


..QUOTE..
  Use pspline within a Cox model.  It includes a fairly general test for
nonlinearity, that is similar to GAM models.
Terry Therneau


 coxph(Surv(time, status) ~ ph.ecog + pspline(age), lung)
Call:
coxph(formula = Surv(time, status) ~ ph.ecog + pspline(age),
data = lung)

 coef   se(coef) se2 Chisq DF   p
ph.ecog  0.4505 0.11766  0.11723 14.66 1.00 0.00013
pspline(age), linear 0.0112 0.00927  0.00927  1.45 1.00 0.23000
pspline(age), nonlin  2.96 3.08 0.41000

Iterations: 4 outer, 10 Newton-Raphson
 Theta= 0.797 
Degrees of freedom for terms= 1.0 4.1
Likelihood ratio test=22.7  on 5.07 df, p=0.000412
  n=227 (1 observation deleted due to missingness)

__
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.


Re: [R] how to check linearity in Cox regression

2008-05-12 Thread Terry Therneau
  Use pspline within a Cox model.  It includes a fairly general test for 
nonlinearity, that is similar to GAM models.
Terry Therneau


 coxph(Surv(time, status) ~ ph.ecog + pspline(age), lung)
Call:
coxph(formula = Surv(time, status) ~ ph.ecog + pspline(age), 
data = lung)

 coef   se(coef) se2 Chisq DF   p  
ph.ecog  0.4505 0.11766  0.11723 14.66 1.00 0.00013
pspline(age), linear 0.0112 0.00927  0.00927  1.45 1.00 0.23000
pspline(age), nonlin  2.96 3.08 0.41000

Iterations: 4 outer, 10 Newton-Raphson
 Theta= 0.797 
Degrees of freedom for terms= 1.0 4.1 
Likelihood ratio test=22.7  on 5.07 df, p=0.000412
  n=227 (1 observation deleted due to missingness)

__
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.


[R] how to check linearity in Cox regression

2008-05-09 Thread array chip
Hi, I am just wondering if there is a test available for testing if a linear 
fit of an independent variable in a Cox regression is enough? Thanks for any 
suggestions.

John Zhang


  


[[elided Yahoo spam]]

__
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.


Re: [R] how to check linearity in Cox regression

2008-05-09 Thread Steven McKinney

 -Original Message-
 From: [EMAIL PROTECTED] on behalf of array chip
 Sent: Fri 5/9/2008 2:54 PM
 To: [EMAIL PROTECTED]
 Subject: [R] how to check linearity in Cox regression
 
 Hi, I am just wondering if there is a test available for testing if a linear 
 fit of an independent variable in a Cox regression is enough? Thanks for any 
 suggestions.

 John Zhang


 require(survival)
 require(splines)
 
 #
 # You could fit a model with the covariate of interest,
 # then an additional fit with other functional forms.
 # I'll use the stanford2 data and age as the covariate.
 # The first fit will have age as a linear term.
 # The second fit will have a quadratic term as well.
 #
 
 stanford2df - stanford2
 stanford2df$agectr - stanford2$age - mean(stanford2$age)
 stanford2df$age2 - stanford2df$agectr^2
 fit1 - coxph(Surv(time, status) ~ agectr, data = stanford2df)
 fit2 - coxph(Surv(time, status) ~ agectr + age2, data = stanford2df)
 summary(fit1)
Call:
coxph(formula = Surv(time, status) ~ agectr, data = stanford2df)

  n= 184 
 coef exp(coef) se(coef)z  p
agectr 0.0292  1.03   0.0106 2.74 0.0061

   exp(coef) exp(-coef) lower .95 upper .95
agectr  1.03  0.971  1.01  1.05

Rsquare= 0.044   (max possible= 0.996 )
Likelihood ratio test= 8.27  on 1 df,   p=0.00403
Wald test= 7.51  on 1 df,   p=0.00613
Score (logrank) test = 7.6  on 1 df,   p=0.00583

 summary(fit2)
Call:
coxph(formula = Surv(time, status) ~ agectr + age2, data = stanford2df)

  n= 184 
  coef exp(coef) se(coef)z   p
agectr 0.04100  1.04 0.010237 4.01 6.2e-05
age2   0.00197  1.00 0.000689 2.86 4.2e-03

   exp(coef) exp(-coef) lower .95 upper .95
agectr  1.04  0.960  1.02  1.06
age21.00  0.998  1.00  1.00

Rsquare= 0.08   (max possible= 0.996 )
Likelihood ratio test= 15.4  on 2 df,   p=0.00046
Wald test= 17.4  on 2 df,   p=0.00017
Score (logrank) test = 17.3  on 2 df,   p=0.000175

 anova(fit1, fit2, test = Chisq)
Analysis of Deviance Table

Model 1: Surv(time, status) ~ agectr
Model 2: Surv(time, status) ~ agectr + age2
  Resid. Df Resid. Dev  Df Deviance P(|Chi|)
1   1831017.94   
2   1821010.84   1 7.10  0.01
 # the likelihood ratio test suggests that a linear term for
 # age alone may not be adequate.  The quadratic form yields
 # a better fit.
 
 # You can also examine functional forms using spline fits.
 fit0 - coxph(Surv(time, status) ~ age, data = stanford2df)
 fit3 - coxph(Surv(time, status) ~ ns(age, df = 4), data = stanford2df)
 pred3 - predict(fit3, type=terms, se=TRUE)
 
 hfit - pred3$fit[,1]
 hse - pred3$se[,1]
 hmat=cbind(hfit, hfit+2*hse,hfit-2*hse)
 o - order(stanford2$age)
 matplot(stanford2$age[o], hmat[o, ], pch=*,col=c(2,4,4), xlab = age,
+ ylab=Log Relative Risk,main=Cox Model: Survival,type=l)
 # The plot of the spline fit for age shows a non-linear form
 
 anova(fit0, fit3, test = Chisq)
Analysis of Deviance Table

Model 1: Surv(time, status) ~ age
Model 2: Surv(time, status) ~ ns(age, df = 4)
  Resid. Df Resid. Dev  Df Deviance P(|Chi|)
1   1831017.94   
2   1801009.45   3 8.49  0.04
 # The likelihood ratio test suggests the linear model
 # can be improved upon.

Hope this helps

Steve McKinney


  


[[elided Yahoo spam]]

__
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.

__
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.