在R中为ggplot绘制减少主轴线和置信区间

Question

有没有办法在ggplot中添加缩小的主轴线（理想情况下还有CI）？我知道我可以使用method =“lm”来获得OLS拟合，但似乎没有默认方法用于RMA。我可以从lmodel2软件包中获取RMA系数和CI区间，但是使用geom_abline()将它们添加到图表中似乎不起作用。这里是虚拟数据和代码。我只想用RMA线和CI替换OLS线和CI：

dat <- data.frame(a=log10(rnorm(50, 30, 10)), b=log10(rnorm(50, 20, 2)))

ggplot(dat, aes(x=a, y=b) ) + 
    geom_point(shape=1) +       
    geom_smooth(method="lm") 

编辑1：下面的代码获取RMA（此处称为SMA-标准化主轴）系数和CI。lmodel2包提供了更详细的输出，而smatr包仅返回系数和CI，如果有帮助的话：

library(lmodel2)
fit1 <- lmodel2(b ~ a, data=dat)

library(smatr)
fit2 <- line.cis(b, a, data=dat)

- Steve

---

请添加您正在使用的 lmodel2() 代码。 - Chase

Answer 1

在Joran的回答基础上，我认为将整个数据框传递给geom_abline会更容易些:

library(ggplot2)
library(lmodel2)

dat <- data.frame(a=log10(rnorm(50, 30, 10)), b=log10(rnorm(50, 20, 2)))
mod <- lmodel2(a ~ b, data=dat,"interval", "interval", 99)

reg <- mod$regression.results
names(reg) <- c("method", "intercept", "slope", "angle", "p-value")

ggplot(dat) + 
  geom_point(aes(b, a)) +
  geom_abline(data = reg, aes(intercept = intercept, slope = slope, colour = method))

- hadley

Answer 2

正如Chase所评论的那样，实际使用的lmodel2()代码和ggplot代码会很有帮助。但是这里有一个例子可能能指引你朝着正确的方向前进：

dat <- data.frame(a=log10(rnorm(50, 30, 10)), b=log10(rnorm(50, 20, 2)))
mod <- lmodel2(a ~ b, data=dat,"interval", "interval", 99)

#EDIT: mod is a list, with components (data.frames) regression.results and 
#        confidence.intervals containing the intercepts+slopes for different 
#        estimation methods; just put the right values into geom_abline
ggplot(dat,aes(x=b,y=a)) + geom_point() + 
   geom_abline(intercept=mod$regression.results[4,2],
            slope=mod$regression.results[4,3],colour="blue") + 
   geom_abline(intercept=mod$confidence.intervals[4,2],
            slope=mod$confidence.intervals[4,4],colour="red") + 
   geom_abline(intercept=mod$confidence.intervals[4,3],
            slope=mod$confidence.intervals[4,5],colour="red") + 
   xlim(c(-10,10)) + ylim(c(-10,10))

完整披露：我对RMA回归一无所知，因此只选取了相关斜率和截距并将它们放入geom_abline()，使用了来自lmodel2的示例代码作为指南。在这个玩具示例中生成的置信区间似乎没有太多意义，因为我不得不使用xlim()和ylim()迫使ggplot缩小以便看到CI线（红色）。

但是也许这会帮助您构建一个在ggplot()中工作的示例。

编辑2：通过OP添加的代码提取系数，ggplot()应该是这样的：

ggplot(dat,aes(x=b,y=a)) + geom_point() + 
geom_abline(intercept=fit2[1,1],slope=fit2[2,1],colour="blue") + 
geom_abline(intercept=fit2[1,2],slope=fit2[2,2],colour="red") + 
geom_abline(intercept=fit2[1,3],slope=fit2[2,3],colour="red")

- joran

---

正是我要写的内容，点赞。需要注意的是，lmodel2()返回一个列表，而regression.results是该列表中的数据框。从那里开始，只需要指定你想从该数据框中提取什么即可。 - Chase

@Chase - 已编辑以包含您的观点，谢谢。令人烦恼的是，regression.results数据框具有一些奇怪的列名，可能是为了视觉效果而设定的。 "Slope"列实际上应该是 " Slope"。这使我稍微有点困扰，并且也是我最终使用"["的原因。 - joran

谢谢，这对于实际行非常有效，但是有什么想法为什么置信区间如此奇怪？此外，在此包中，RMA=范围主轴，SMA=减少（标准化）主轴，因此行号将为3而不是4。我在问题的编辑中提供了代码的smatr包实际上更容易获得系数。 - Steve

没关系 - 一旦我把行号更正为3，置信区间就正常工作了，非常感谢。 - Steve

@Steve - 我一点也不清楚；对这种回归类型不是很了解。你可以去stats.stackexchange.com上问这个问题，那里可能会有更好的答案。 - joran

Answer 3

我发现自己陷入了同样的境地。

1. 使用 ggpmisc 包获得拟合值及其置信区间：
   library(ggpmisc)
ci <- predict.lm(fit1, method= 'RMA', interval= "confidence")

2. 将模型预测添加到您的数据中：
   datci <- cbind(dat, ci)

3. 使用 geom_smooth 参数，如透明度和线宽（当然，您可以自定义它们）：
   p <- ggplot(datci, aes(x= b, y= a)) + geom_point() + geom_line(aes(x= b, y= a), lwd= 1.1, alpha= 0.6)

4. 如果要添加置信区间，请使用 geom_ribbon ：
   p + geom_ribbon(aes(ymin= lwr, ymax= upr, fill= feather), alpha= 0.3, color= NA)

- Santi XGR