绘图和 beta 回归输出不匹配
Plot and beta- regression output do not match
在拟合似乎没有意义的 beta 回归模型后，我使用 tab_model() 得到了一个奇怪的输出。截距和系数似乎与情节不符。例如 soil_moisture 的 CI 似乎太高了。另外，因为我必须回归线，所以不应该有两个截距吗？我错过了什么？
这是我的数据：
df <- structure(list(weed_coverage = c(0.002, 0.01, 0.015, 0.01, 0.001,
0.03, 0.006, 0.012, 0.03, 0.01, 0.002, 0.05, 0.004, 0.02, 0.02,
0.006, 0.03, 0.01, 0.015, 0.01), soil_moisture = c(0.1125, 0.1343,
0.1662, 0.3402, 0.2195, 0.1923, 0.2277, 0.2577, 0.148, 0.2715,
0.104, 0.1495, 0.2788, 0.3477, 0.1835, 0.3175, 0.134, 0.3488,
0.3477, 0.1097), distance = structure(c(2L, 2L, 1L, 2L, 1L, 1L,
1L, 2L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L), .Label = c("2",
"5"), class = "factor"), wc_pr = c(`1` = 0.0144421159096612,
`2` = 0.0148173851570077, `3` = 0.0146600960637327, `4` = 0.0188698067335207,
`5` = 0.0132256975788894, `6` = 0.0139395025623511, `7` = 0.0130176686719618,
`8` = 0.0171297102952414, `9` = 0.0150581171360966, `10` = 0.0119600057973879,
`11` = 0.0142983464947494, `12` = 0.0150847074475541, `13` = 0.0117921770637613,
`14` = 0.019036340203918, `15` = 0.0141784723499839, `16` = 0.0109405819057696,
`17` = 0.0148121562892363, `18` = 0.0190608859962305, `19` = 0.0103185336737258,
`20` = 0.0163480105406738)), class = "data.frame", row.names = c(NA,
-20L))
这是我的代码
library(sjPlot)
library(betareg)
library(ggplot2)
betareg (weed_coverage ~ soil_moisture * distance, data = df) -> model_b # fit beta regression model
tab_model(model_b)
df %>% mutate(wc_pr= predict(model_b , type = "response")) -> df # create column with prediction values for weed_coverage
ggplot(df, aes(x = soil_moisture, y = weed_coverage, color = distance)) + # Plot the model
geom_point(size = 2, shape = 21) +
geom_line(aes(y = wc_pr, color = distance), data = df)+
theme_bw()
tab_model(model_b)
的输出
这是我得到的情节
好的，所以这里发生了很多事情，我将要经历这些事情。第一个是最简单的，就是明确这个模型只有一个截距，就是你的自变量等于0时weed_coverage的估计值。因此，这些值并不总是可以直接解释的。
图中的两条回归线是 soil_moisture 在两个不同距离值处的预测斜率。这是因为存在交互作用，所以 soil_moisture 的直线斜率取决于 distance.
的水平
现在，tab_model() 正确报告了 CI，但 sjPlot 包中大多数函数的默认行为是使用 exp() 函数在对数刻度上转换值，以获得更接近概率值的东西。如果您想要模型适合的日志 space 中的原始模型值，您必须像这样指定 transform=NULL：
tab_model(model_b,
transform = NULL)
您现在会注意到模型报告的数字与您调用以下命令时的数字相同：
summary(model_b)
Call:
betareg(formula = weed_coverage ~ soil_moisture * distance, data = df)
Standardized weighted residuals 2:
Min 1Q Median 3Q Max
-2.2157 -0.5334 -0.0546 0.8254 1.7902
Coefficients (mean model with logit link):
Estimate Std. Error z value Pr(>|z|)
(Intercept) -3.8822 0.7840 -4.952 7.35e-07 ***
soil_moisture -1.9591 3.3122 -0.591 0.554
distance5 -0.4752 0.9303 -0.511 0.610
soil_moisture:distance5 3.1533 3.9372 0.801 0.423
Phi coefficients (precision model with identity link):
Estimate Std. Error z value Pr(>|z|)
(phi) 102.30 34.72 2.946 0.00322 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Type of estimator: ML (maximum likelihood)
Log-likelihood: 65.73 on 5 Df
Pseudo R-squared: 0.05632
Number of iterations: 399 (BFGS) + 6 (Fisher scoring)
CI 非常大，但这是因为 weed_coverage 的所有值都非常低。在log-oddsspace中，小的百分比或比例变化在接近0或1时变得非常大。在这个示例模型中，模型非常不确定，这就是为什么back-transformed概率尺度上的 CI 非常大。还值得注意的是，这些估计是针对 weed_coverage 中的变化 soil_moisture 中的 1 的变化。如果你缩放这个变量，模型不会改变，但数字会更容易在你的数据的上下文中解释，而不是估计回归线中远远超出你的数据范围的不确定性，一切都会报告soil_moisture 的 1 个标准偏差的变化：
df$scaled_soil_moisture <- scale(df$soil_moisture)
model_b2 <- betareg(weed_coverage ~ scaled_soil_moisture * distance, data = df)
# Notice that I have let it back-transform the variable again.
tab_model(model_b2)
如果有什么不明白的地方请告诉我。
编辑：为了澄清您在评论中提出的关于截距的问题，我绘制了带有缩放自变量的回归的两条回归线。此模型中的截距对应于 weed_coverage 的值，当 scaled_soil_moisture 等于零且距离等于 2（因为距离是一个 two-level 因子，默认行为通常是创建一个虚拟变量，第一级由截距表示）。如您所见，在 x-axis 等于 0 时的红线上，y-axis 大约为 0.01（已向下舍入）。这与 non-scaled 模型形成对比，其中截距未在图表上显示，因为数据中没有零值。如果您将红线继续穿过数据到达 x=0 的位置，您会看到它穿过 y-axis 大约 0.02.
This page 说明这是一种相当直观的方式。
在拟合似乎没有意义的 beta 回归模型后，我使用 tab_model() 得到了一个奇怪的输出。截距和系数似乎与情节不符。例如 soil_moisture 的 CI 似乎太高了。另外，因为我必须回归线，所以不应该有两个截距吗？我错过了什么？
这是我的数据：
df <- structure(list(weed_coverage = c(0.002, 0.01, 0.015, 0.01, 0.001,
0.03, 0.006, 0.012, 0.03, 0.01, 0.002, 0.05, 0.004, 0.02, 0.02,
0.006, 0.03, 0.01, 0.015, 0.01), soil_moisture = c(0.1125, 0.1343,
0.1662, 0.3402, 0.2195, 0.1923, 0.2277, 0.2577, 0.148, 0.2715,
0.104, 0.1495, 0.2788, 0.3477, 0.1835, 0.3175, 0.134, 0.3488,
0.3477, 0.1097), distance = structure(c(2L, 2L, 1L, 2L, 1L, 1L,
1L, 2L, 2L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L), .Label = c("2",
"5"), class = "factor"), wc_pr = c(`1` = 0.0144421159096612,
`2` = 0.0148173851570077, `3` = 0.0146600960637327, `4` = 0.0188698067335207,
`5` = 0.0132256975788894, `6` = 0.0139395025623511, `7` = 0.0130176686719618,
`8` = 0.0171297102952414, `9` = 0.0150581171360966, `10` = 0.0119600057973879,
`11` = 0.0142983464947494, `12` = 0.0150847074475541, `13` = 0.0117921770637613,
`14` = 0.019036340203918, `15` = 0.0141784723499839, `16` = 0.0109405819057696,
`17` = 0.0148121562892363, `18` = 0.0190608859962305, `19` = 0.0103185336737258,
`20` = 0.0163480105406738)), class = "data.frame", row.names = c(NA,
-20L))
这是我的代码
library(sjPlot)
library(betareg)
library(ggplot2)
betareg (weed_coverage ~ soil_moisture * distance, data = df) -> model_b # fit beta regression model
tab_model(model_b)
df %>% mutate(wc_pr= predict(model_b , type = "response")) -> df # create column with prediction values for weed_coverage
ggplot(df, aes(x = soil_moisture, y = weed_coverage, color = distance)) + # Plot the model
geom_point(size = 2, shape = 21) +
geom_line(aes(y = wc_pr, color = distance), data = df)+
theme_bw()
tab_model(model_b)
的输出这是我得到的情节
好的，所以这里发生了很多事情，我将要经历这些事情。第一个是最简单的，就是明确这个模型只有一个截距，就是你的自变量等于0时weed_coverage的估计值。因此，这些值并不总是可以直接解释的。
图中的两条回归线是 soil_moisture 在两个不同距离值处的预测斜率。这是因为存在交互作用，所以 soil_moisture 的直线斜率取决于 distance.
现在，tab_model() 正确报告了 CI，但 sjPlot 包中大多数函数的默认行为是使用 exp() 函数在对数刻度上转换值，以获得更接近概率值的东西。如果您想要模型适合的日志 space 中的原始模型值，您必须像这样指定 transform=NULL：
tab_model(model_b,
transform = NULL)
您现在会注意到模型报告的数字与您调用以下命令时的数字相同：
summary(model_b)
Call:
betareg(formula = weed_coverage ~ soil_moisture * distance, data = df)
Standardized weighted residuals 2:
Min 1Q Median 3Q Max
-2.2157 -0.5334 -0.0546 0.8254 1.7902
Coefficients (mean model with logit link):
Estimate Std. Error z value Pr(>|z|)
(Intercept) -3.8822 0.7840 -4.952 7.35e-07 ***
soil_moisture -1.9591 3.3122 -0.591 0.554
distance5 -0.4752 0.9303 -0.511 0.610
soil_moisture:distance5 3.1533 3.9372 0.801 0.423
Phi coefficients (precision model with identity link):
Estimate Std. Error z value Pr(>|z|)
(phi) 102.30 34.72 2.946 0.00322 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Type of estimator: ML (maximum likelihood)
Log-likelihood: 65.73 on 5 Df
Pseudo R-squared: 0.05632
Number of iterations: 399 (BFGS) + 6 (Fisher scoring)
CI 非常大，但这是因为 weed_coverage 的所有值都非常低。在log-oddsspace中，小的百分比或比例变化在接近0或1时变得非常大。在这个示例模型中，模型非常不确定，这就是为什么back-transformed概率尺度上的 CI 非常大。还值得注意的是，这些估计是针对 weed_coverage 中的变化 soil_moisture 中的 1 的变化。如果你缩放这个变量，模型不会改变，但数字会更容易在你的数据的上下文中解释，而不是估计回归线中远远超出你的数据范围的不确定性，一切都会报告soil_moisture 的 1 个标准偏差的变化：
df$scaled_soil_moisture <- scale(df$soil_moisture)
model_b2 <- betareg(weed_coverage ~ scaled_soil_moisture * distance, data = df)
# Notice that I have let it back-transform the variable again.
tab_model(model_b2)
如果有什么不明白的地方请告诉我。
编辑：为了澄清您在评论中提出的关于截距的问题，我绘制了带有缩放自变量的回归的两条回归线。此模型中的截距对应于 weed_coverage 的值，当 scaled_soil_moisture 等于零且距离等于 2（因为距离是一个 two-level 因子，默认行为通常是创建一个虚拟变量，第一级由截距表示）。如您所见，在 x-axis 等于 0 时的红线上，y-axis 大约为 0.01（已向下舍入）。这与 non-scaled 模型形成对比，其中截距未在图表上显示，因为数据中没有零值。如果您将红线继续穿过数据到达 x=0 的位置，您会看到它穿过 y-axis 大约 0.02.
This page 说明这是一种相当直观的方式。