R help - Apr 2004 - Regression models w/ splines

Hi - I am fitting various Cox PH models with spline predictors.  After
fitting the model, I would like to use termplot() to examine the
functional form of the fitted model (e.g., to obtain a plot of the
relative risk (or log r.r.) versus the predictors). 

When there is only 1 predictor in the model, termplot returns a "?". 
In this case, I have not been able to figure out how to create the
desired plot.  Otherwise,termplot is easy to use. 

Example:
library(survival)
library(splines)

data(pbc)

fit1<-coxph(Surv(time,status)~ns(bili,4)+ns(age,4), data=pbc)
fit2<-coxph(Surv(time,status)~ns(bili,4), data=pbc)

par(mfrow=c(2,2))
termplot(fit1) # returns two plots (one for bili and one for age)
termplot(fit2) # returns a "?" prompt

Does anyone have any suggestions for obtaining a plot of the log
relative risk as a function of bili in the second example (or for this
type of model in general)?

On a similar note, I would like to be able to obtain predicted values
for "new" subjects.  I tried using:

temp<-data.frame(bili=0.8)
predict(fit2,newdata=temp, safe=T)

This gives the following error message:
Error in qr(t(const)) : NA/NaN/Inf in foreign function call (arg 1)

Is there a simple method available for obtaining predicted values for
new subjects when using these types of models?

I have encountered similar problems when trying to use smoothing
splines w/ the above examples (i.e., 'pspline' in the survival package).
 I am currently using R version 1.8.1 for windows. Any help would be
greatly appreciated!

John

John Fieberg, Ph.D.
Wildlife Biometrician, MN DNR
5463-C W. Broadway
Forest Lake, MN 55434
Phone: (651) 296-2704

"John Fieberg" <John.Fieberg at dnr.state.mn.us> writes:
> Hi - I am fitting various Cox PH models with spline predictors.  After
> fitting the model, I would like to use termplot() to examine the
> functional form of the fitted model (e.g., to obtain a plot of the
> relative risk (or log r.r.) versus the predictors). 
> 
> When there is only 1 predictor in the model, termplot returns a
"?".
> In this case, I have not been able to figure out how to create the
> desired plot.  Otherwise,termplot is easy to use. 
Kasper Hansen and Per Jensen from my department reported this
recently, AFAIR with a suggested fix. You might want to go over the
archives for the last month or so.

-- 
   O__  ---- Peter Dalgaard             Blegdamsvej 3  
  c/ /'_ --- Dept. of Biostatistics     2200 Cph. N   
 (*) \(*) -- University of Copenhagen   Denmark      Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk)             FAX: (+45) 35327907

On Fri, 09 Apr 2004 15:16:34 -0500
"John Fieberg" <John.Fieberg at dnr.state.mn.us> wrote:
> Hi - I am fitting various Cox PH models with spline predictors.  After
> fitting the model, I would like to use termplot() to examine the
> functional form of the fitted model (e.g., to obtain a plot of the
> relative risk (or log r.r.) versus the predictors). 
This kind of application is why the Design package was created.  It is
very easy to obtain such plots:

library(Design)
dd <- datadist(mydata); options(datadist='dd')
f <- cph(Surv(y,d) ~ sex*rcs(age,4)+rcs(bp,5))
plot(f, age=NA, sex='female')

---
Frank E Harrell Jr   Professor and Chair           School of Medicine
                     Department of Biostatistics   Vanderbilt University