14. Estimation with external information

AAA data set for adaptation method case study

R code for adaptation method case study

The R code for the core analyses is included in the file below.
r code for adaptation.doc

# 7 predictor model in individual patient data
full <- lrm(STATUS~SEX+AGE10+MI+CHF+ISCHEMIA+LUNG+RENAL,data=AAA,x=T,y=T)

# Estimate shrinkage factor
full.validate <- validate(full, B=200)
# Apply shrinkage factor on model coefficients
full.shrunk.coefficients	<- full.validate["Slope","index.corrected"] * full$coefficients
# Adjust intercept; alternatively, use as offset variable 
full.shrunk.coefficients[1] <- full.shrunk.coefficients[1] + full.validate["Intercept","index.corrected"]

# Full model, penalized estimation
penalty	<- pentrace(full,penalty=c(0.5,1,2,3,4,6,8,12,16,24),maxit=25)
full.penalized <- update(full, penalty=penalty$penalty)

# Adaptation method; start with univariate literature coefficients
# analyzed with random effects meta-analysis
lit.coefficients	<- c(0.3606,0.788,1.034,1.590,1.514,0.8502,1.302)
lit.se		<- c(0.176,0.112,0.317,0.4109,0.378,0.2367,0.2495)
# Univariate coefs in individual patient data, followed by simple adaptation method
fit.ind.x	<- rep(0,7)
for (i in 1:7) {fit.ind.x	<- lrm.fit(y=full$y,x=full$x[,i]) 
# Adaptation method 1, identical to proposal Greenland
adapt.coef[i] 	<- full$coefficients[i+1] + 1 * (lit.coefficients[i] - fit.ind.x$coefficients[2])
adapt.var[i]	<- full$var[i+1,i+1] - fit.ind.x$var[2,2] + lit.se[i]^2 }
# We note that the variance has decreased considerably
sqrt(adapt.var) / sqrt(diag(full$var)[2:8])

# Round the coefficients after shrinkage 0.9
cat("\n",round(9*adapt.coef,0))
# we estimate the intercept, first scale age around 70 years (7 decades)
X	<- as.matrix(full$x)
X[,2] <- X[,2] – 7	# age rescaled
# offset with scores / 10
offset.AAA	<- X %*% c(.3,.3,.2,.8,.7,.6,1.1)
fit.offset  <- lrm.fit(y=full$y, offset=offset.AAA)
# intercept given these scores for predictors
fit.offset$coef[1]	# result: -3.48
fit.offset$coef[1] + log((0.05/0.95)/(18/220))		# result: -3.92


# there is a function ‘bootcor.uni.mult’ which performs bootstrapping to obtain the correlation between univariate and multivariable logistic regression coefficients in the individual patient data; this correlation is required for Adaptation method 2 
# Obtain correlations uni – mult, and shrinkage factor
full.uni.mult.cor	<- bootcor.uni.mult(full,B=200,maxit=10,trim=0.1,save.indices = T)
# Adaptation method 2, with c optimal
for (i in 1:7) {
fit.ind.x		<- lrm.fit(y=full$y,x=full$x[,i],maxit=25)
adapt.factors[i]	<- (full.uni.mult.cor$r[i] * sqrt(full$var[i+1,i+1]) *
 sqrt(fit.ind.x$var[2,2]) ) / (lit.se[i]^2 + fit.ind.x$var[2,2])
adapt.coefficients[i]	<- full$coefficients[i+1] + adapt.factors[i] * 
(lit.coefficients[i] - fit.ind.x$coefficients[2])
adapt.var[i]		<- full$var[i+1,i+1] * (1-(full.uni.mult.cor$r[i]^2 * fit.ind.x$var[2,2] / 
				(lit.se[i]^2 + fit.ind.x$var[2,2]))) } 
# end loop over 7 predictors

## Recalibrate literature coefficients with one calibration factor ##
# make offset based on univariate literature coefficients
offset.AAA2	<- X %*% lit.coefficients
fit.offset2  <- lrm.fit(y=full$y, x=offset.AAA2)
fit.offset2$stats
round(fit.offset2$coef[2],2)	# result: 0.69 
# scores with this overall recalibration
cat("\n",round(10*fit.offset2$coef[2]*lit.coefficients,0),"\n")
# use these scores / 10 as offset for intercept
offset.AAA3	<- X %*% c(.3,.5,.7,1.1,1.0,.6,.9)
fit.offset3  <- lrm.fit(y=full$y, offset=offset.AAA3)
fit.offset3$coef[1] # -4.05
fit.offset3$coef[1] + log((0.05/0.95)/(18/220)) # -4.49 

rcode_and_data/chapter14.txt · Last modified: 2008/07/17 15:33 by afcgerritsen
chimeric.de = chi`s home Creative Commons License Valid CSS Driven by DokuWiki do yourself a favour and use a real browser - get firefox!! Recent changes RSS feed Valid XHTML 1.0