【问题标题】:add confidence interval to splines from quantile regression从分位数回归向样条曲线添加置信区间
【发布时间】:2018-05-07 23:05:32
【问题描述】:

我正在使用 quantreg 来拟合分位数回归,并且还添加了一些基于 x 值(分位数)的节点。我现在想绘制这个并且还有置信区间。不知道该怎么做。

这是一个可重现的例子:

  #create data    
  x <- seq(0,100,length.out = 100)        
  sig <- 0.1 + 0.05*x 
  b_0 <- 6                                
  b_1 <- 0.1                              
 set.seed(1)                             
 e <- rnorm(100,mean = 0, sd = sig)      
 y <- b_0 + b_1*x + e 

mydata <- data.frame(x,y, age=sample(30:70,100,replace=TRUE), sex=sample(c("Male","Female"),100, replace=TRUE))

#run regression
library(quantreg)
library(splines)
model <- rq(y ~ ns(x, knots=c(25,50,75))+age+sex, tau=0.5, data=mydata ) 

#plot
sp <- c(25,50,75)
ggplot(mydata, aes(x=x,y=y))+ geom_point()+ geom_quantile(formula=y~ns(x,knots=sp), quantiles=0.5, se=T)

这不显示置信区间?? 另外,这个情节没有考虑协变量吗?有没有办法做到这一点?

【问题讨论】:

  • 您可以添加到 quantiles 参数以获得快速结果,geom_quantile(formula = y ~ ns(x, knots = sp), quantiles = c(0.05, 0.5, 0.95)),但是,如果要使用功能区,则需要不同的方法

标签: r plot spline quantile


【解决方案1】:

您实际上并没有在情节中使用来自model 的结果。所以我不完全确定你想要做什么。

您可以使用predict 获得基于分位数回归模型的预测,包括置信区间,并使用geom_ribbon 绘制这些预测:

# Get model-based predictions
pred <- as.data.frame(predict(model, data.frame(x = mydata$x, age = mydata$age, sex = mydata$sex), interval = "confidence"));
pred$x <- mydata$x;

#plot
sp <- c(25,50,75);
ggplot(mydata, aes(x=x,y=y)) +
    geom_point() +
    geom_line(data = pred, aes(x = x, y = fit)) +
    geom_ribbon(data = pred, aes(ymin = lower, ymax = higher, x = x), alpha = 0.4) +
    geom_quantile(formula = y ~ ns(x, knots = sp), quantiles = 0.5, se = T);

【讨论】:

  • 这太棒了!谢谢!是否可以让置信区间更平滑?
  • @user6121484 您需要评估模型拟合的质量;在model 中,您将xagesex 视为预测变量。稍后,当 (gg) 绘图时,您从 不同的 分位数回归模型添加拟合,您只在 x 上回归。我不确定您要做什么,但是如果您从 model 中删除 agesex,您将重现更平滑的蓝色拟合曲线,包括置信区间。
  • 在模型中,我正在校正年龄和性别。所以我想绘制 x 和 y 之间的关系,并根据年龄和性别进行校正。
  • @user6121484 在这种情况下,您可能想要删除geom_quantile(...) 位,因为您只是在x 上回归。就个人而言,我会使用 GLM 而不是分位数回归来探索不同预测变量对响应 y 的作用。或者,如果您只对平滑趋势的可视化感兴趣,也许是非参数 LOESS 拟合。也存在带有平滑的分位数回归模型,参见例如here.
猜你喜欢
  • 2018-10-30
  • 1970-01-01
  • 2018-03-05
  • 1970-01-01
  • 2018-05-05
  • 1970-01-01
  • 1970-01-01
  • 2010-11-13
  • 2019-08-29
相关资源
最近更新 更多