假设我想画一个类似于这里的图,这里使用R,即y轴是风险比,x轴是基于样条项的Cox模型的一些预测变量。唯一的例外是我想设置自己的x轴点;termplot似乎从数据集中挑选了所有唯一的值,但我想使用某种网格。这是因为我在做多次imputation,每一轮都会产生不同的唯一值。否则,我可以很容易地进行组合推断,但在每个估算轮中对相同的预测值进行预测会容易得多。
所以,我需要找到一种方法来使用termplot函数,这样我就可以修复预测值或找到一个解决方案。我试图使用预测函数,但它的newdata参数也需要所有其他(调整)变量的值,这会增加标准误差。这是一个问题,因为我也在绘制置信区间。我想我可以手动完成这个,不需要任何函数,除了样条项在这个意义上超出了我的范围。
这是一个说明性的例子。
library(survival)
data(diabetic)
diabetic<-diabetic[diabetic$eye=="right",]
# Model with spline term and two adjusting variables
mod<-coxph(Surv(time,status)~+pspline(age,df=3)+risk+laser,data=diabetic)
summary(mod)
# Let's pretend this is the grid
# These are in the data but it's easier for comparison in this example
ages<-20:25
# Right SEs but what if I want to use different age values?
termplot(mod,term=1,se=TRUE,plot=F)$age[20:25,]
# This does something else
termplot(mod,data=data.frame(age=ages),term=1,se=TRUE,plot=F)$age
# This produces an error
predict(mod,newdata=data.frame(age=ages),se.fit=T)
# This inflates variance
# May actually work with models without categorical variables: what to do with them?
# Actual predictions are different but all that matters is the difference between ages
# and they are in line
predict(mod,newdata=data.frame(age=ages,risk=mean(diabetic$risk),laser="xenon"),se.fit=T)
请告诉我是否没有充分解释我的问题。我尽量让它保持简单。
最后,我是这样做的。首先,我使用termplot函数进行预测和se,然后我使用线性插值为我的自定义网格获得近似正确的预测和se。
ptemp<-termplot(mod,term=1,se=TRUE,plot=F)
ptemp<-data.frame(ptemp[1]) # Picking up age and corresponding estimates and SEs
x<-ptemp[,1]; y<-ptemp[,2]; se<-ptemp[,3]
f<-approxfun(x,y) # Linear interpolation function
x2<-seq(from=10,to=50,by=0.5) # You can make a finer grid if you want
y2<-f(x2) # Interpolation itself
f_se<-approxfun(x,se) # Same for SEs
se2<-f_se(x2)
dat<-data.frame(x2,y2,se2) # The wanted result