R help - Aug 2010 - How to compare the effect of a variable across regression models?

Hello,

I would like, if it is possible, to compare the effect of a variable across 
regression models. I have looked around but I haven't found anything. Maybe 
someone could help? Here is the problem:

I am studying the effect of a variable (age) on an outcome (local recurrence: 
lr). I have built 3 models:
- model 1: lr ~ age      y = \beta_(a1).age
- model 2: lr ~ age +  presentation variables (X_p)        y = \beta_(a2).age + 
\BETA_(p2).X_p
- model 3: lr ~ age + presentation variables + treatment variables( X_t) 
       y = \beta_(a3).age  + \BETA_(p3).X_(p) + \BETA_(t3).X_t
 
Presentation variables include variables such as tumor grade, tumor size, etc...
the physician cannot interfer with these variables.
Treatment variables include variables such as chemotherapy, radiation, surgical 
margins (a surrogate for adequate surgery).

I have used cph for the models and restricted cubic splines (Design library) for
age. I have noted that the effect of age decreases from model 1 to 3.

I would like to compare the effect of age on the outcome across the different 
models. A test of \beta_(a1) = \beta_(a2) = \beta_(a3) and then two by two 
comparisons or a global trend test maybe? Is that possible?

Thank you for your help,


David Biau.



      
	[[alternative HTML version deleted]]

David,

In the Cox and many other regression models, the effect of a variable 
is context-dependent.  There is an identifiability problem in what you 
are doing, as discussed by

@ARTICLE{for95mod,
   author = {Ford, Ian and Norrie, John and Ahmadi, Susan},
   year = 1995,
   title = {Model inconsistency, illustrated by the {Cox} proportional 
hazards
           model},
   journal = Stat in Med,
   volume = 14,
   pages = {735-746},
   annote = {covariable adjustment; adjusted estimates; baseline 
imbalances;
            RCT; model misspecification; model identification}
}

One possible remedy, which may not work for your goals, is to embed 
all models in a grand model that is used for inference.

When coefficients ARE comparable in some sense, you can use the 
bootstrap to get confidence bands for differences in regressor effects 
between models.

Frank

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

On Fri, 13 Aug 2010, Biau David wrote:
> Hello,
>
> I would like, if it is possible, to compare the effect of a variable across
> regression models. I have looked around but I haven't found anything.
Maybe
> someone could help? Here is the problem:
>
> I am studying the effect of a variable (age) on an outcome (local
recurrence:
> lr). I have built 3 models:
> - model 1: lr ~ age      y = \beta_(a1).age
> - model 2: lr ~ age +  presentation variables (X_p)        y =
\beta_(a2).age +
> \BETA_(p2).X_p
> - model 3: lr ~ age + presentation variables + treatment variables( X_t)
>       y = \beta_(a3).age  + \BETA_(p3).X_(p) + \BETA_(t3).X_t
>
> Presentation variables include variables such as tumor grade, tumor size,
etc...
> the physician cannot interfer with these variables.
> Treatment variables include variables such as chemotherapy, radiation,
surgical
> margins (a surrogate for adequate surgery).
>
> I have used cph for the models and restricted cubic splines (Design
library) for
> age. I have noted that the effect of age decreases from model 1 to 3.
>
> I would like to compare the effect of age on the outcome across the
different
> models. A test of \beta_(a1) = \beta_(a2) = \beta_(a3) and then two by two
> comparisons or a global trend test maybe? Is that possible?
>
> Thank you for your help,
>
>
> David Biau.
>
>
>
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at 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.
>

If you just want to visualize the effect on one variable on the response from
some different models then you might try Predict.Plot from the TeachingDemos
package.  It takes a little tweaking to get it to work with cph objects, but
here is a basic example (partly stolen from the help page for cph):

library(rms)
library(TeachingDemos)


     set.seed(731)
     age <- 50 + 12*rnorm(n)
     label(age) <- "Age"
     sex <- factor(sample(c('Male','Female'), n, 
                   rep=TRUE, prob=c(.6, .4)))
     cens <- 15*runif(n)
     h <- .02*exp(.04*(age-50)+.8*(sex=='Female'))
     dt <- -log(runif(n))/h
     label(dt) <- 'Follow-up Time'
     e <- ifelse(dt <= cens,1,0)
     dt <- pmin(dt, cens)
     units(dt) <- "Year"
  

tmp.df <- data.frame(dt=dt, e=e, age=age, sex=sex)
f <- cph(Surv(dt,e) ~ rcs(age,4) + sex, data=tmp.df )
f$data <- tmp.df
Predict.Plot(f, 'age', age=50, sex='Male', type='lp', 
	plot.args=list(col='blue',ylim=c(-1,2)))
Predict.Plot(f, 'age', age=50, sex='Female', type='lp',
add=TRUE,
	plot.args=list(col='red'))

What that all means depends on the things mentioned in the other replies.  Hope
this helps,

-- 
Gregory (Greg) L. Snow Ph.D.
Statistical Data Center
Intermountain Healthcare
greg.snow at imail.org
801.408.8111

> -----Original Message-----
> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-
> project.org] On Behalf Of Biau David
> Sent: Friday, August 13, 2010 7:42 AM
> To: r help list
> Subject: [R] How to compare the effect of a variable across regression
> models?
> 
> Hello,
> 
> I would like, if it is possible, to compare the effect of a variable
> across
> regression models. I have looked around but I haven't found anything.
> Maybe
> someone could help? Here is the problem:
> 
> I am studying the effect of a variable (age) on an outcome (local
> recurrence:
> lr). I have built 3 models:
> - model 1: lr ~ age      y = \beta_(a1).age
> - model 2: lr ~ age +  presentation variables (X_p)        y >
\beta_(a2).age +
> \BETA_(p2).X_p
> - model 3: lr ~ age + presentation variables + treatment variables(
> X_t)
>        y = \beta_(a3).age  + \BETA_(p3).X_(p) + \BETA_(t3).X_t
> 
> Presentation variables include variables such as tumor grade, tumor
> size, etc...
> the physician cannot interfer with these variables.
> Treatment variables include variables such as chemotherapy, radiation,
> surgical
> margins (a surrogate for adequate surgery).
> 
> I have used cph for the models and restricted cubic splines (Design
> library) for
> age. I have noted that the effect of age decreases from model 1 to 3.
> 
> I would like to compare the effect of age on the outcome across the
> different
> models. A test of \beta_(a1) = \beta_(a2) = \beta_(a3) and then two by
> two
> comparisons or a global trend test maybe? Is that possible?
> 
> Thank you for your help,
> 
> 
> David Biau.
> 
> 
> 
> 
> 	[[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help at 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.