-1

我正在对高度偏斜的比例(中位数~2%,平均~13%)进行 beta 回归。然而,来自模型的平均预测值通常大约是样本均值的两倍。甚至最小预测值也大于样本均值!我已经尝试了 link 和 link.phi 函数的所有变体,结果保持不变。

关于如何在预测分布和样本分布之间获得更好匹配的任何建议?谢谢!

library(betareg)

data <- structure(list(y = c(0.09053397, 0.02396534, 0.0618883, 1.9e-07, 
                             0.07121946, 1.9e-07, 1.9e-07, 0.02395879, 0.01905196, 0.06781774, 
                             1.9e-07, 1.9e-07, 0.99999981, 1.9e-07, 1.9e-07, 0.19839486, 1.9e-07, 
                             0.05087847, 1.9e-07, 0.00257154, 0.08083247, 0.02506733, 1.9e-07, 
                             1.9e-07, 0.19138035, 0.11433343, 0.02341815, 0.99999981, 0.01890199, 
                             1.9e-07, 1.9e-07, 0.09462199, 1.9e-07, 0.62029209, 0.01332234, 
                             1.9e-07, 1.9e-07, 0.00053404, 1.9e-07, 0.03991713, 0.05449004, 
                             0.39588225, 1.9e-07, 1.9e-07, 0.21705814, 0.01063151, 0.08353221, 
                             0.00856108, 0.11422871, 0.99999981, 0.04512786, 0.12043445, 1.9e-07, 
                             0.01948717, 0.14161796, 0.00139597, 1.9e-07, 0.12162885, 0.79714543, 
                             0.24379762, 1.9e-07, 0.64810177, 0.21043641, 0.12024897, 1.9e-07, 
                             0.51773257, 0.99999981, 1.9e-07, 0.33720735, 1.9e-07, 0.03963053, 
                             0.00405092, 0.04025961, 0.16409692, 1.9e-07, 0.64428465, 1.9e-07, 
                             1.9e-07, 0.07384747, 1.9e-07, 0.01558129, 0.0366185, 0.03218959, 
                             0.01835648, 0.13076873, 1.9e-07, 0.03025186, 0.0914969, 0.00021478, 
                             1.9e-07, 1.9e-07, 1.9e-07, 0.02444584, 0.02189736, 0.097028, 
                             0.01347678, 1.9e-07, 0.63884538, 0.09762437, 0.5166029), 
                       x = c(0.07825506, 0.08064626, 0.08514504, 1.9e-07, 1.9e-07, 1.9e-07, 1.9e-07, 0.00990069, 
                              1.9e-07, 0.08254355, 0.87667832, 1.9e-07, 0.64975944, 1.9e-07, 
                              1.9e-07, 0.22367289, 1.9e-07, 0.04823323, 1.9e-07, 0.43635115, 
                              0.08165183, 0.02220256, 1.9e-07, 1.9e-07, 1.9e-07, 0.01418114, 
                              0.01097939, 1.9e-07, 1.9e-07, 0.03326785, 1.9e-07, 0.06957586, 
                              1.9e-07, 0.25541391, 0.32377859, 0.23855679, 0.02837065, 0.06388796, 
                              1.9e-07, 1.9e-07, 0.05840598, 0.01999456, 1.9e-07, 1.9e-07, 0.40921674, 
                              1.9e-07, 0.002533, 0.04552465, 0.1046386, 0.75826368, 0.40296593, 
                              1.9e-07, 0.02852656, 0.05545389, 1.9e-07, 0.00030394, 0.01199183, 
                              0.21549378, 0.05438194, 0.03442961, 0.15914659, 0.20376901, 0.09622411, 
                              0.38810661, 1.9e-07, 0.90504166, 1.9e-07, 0.00803977, 0.03232884, 
                              1.9e-07, 0.02517311, 0.00013941, 0.11415697, 1.9e-07, 1.9e-07, 
                              0.42301297, 1.9e-07, 1.9e-07, 0.05489811, 1.9e-07, 0.3773369, 
                              0.23680447, 0.01163705, 0.01191086, 0.28056447, 1.9e-07, 0.0030056, 
                              0.03705583, 0.0637212, 1.9e-07, 0.00752029, 1.9e-07, 0.10052723, 
                              0.04053147, 0.02817251, 0.0165664, 1.9e-07, 1.9e-07, 0.01264018, 
                              0.018953)), 
                  row.names = c(NA, -100L), 
                  class = c("tbl_df", "tbl", "data.frame"))


model <- betareg(formula = "y ~ x", data = data)
data$pred <- predict(model, newdata = data)

summary(data)
       y                   x                  pred       
 Min.   :0.0000002   Min.   :0.0000002   Min.   :0.1854  
 1st Qu.:0.0000002   1st Qu.:0.0000002   1st Qu.:0.1854  
 Median :0.0236885   Median :0.0153738   Median :0.1918  
 Mean   :0.1277322   Mean   :0.0957217   Mean   :0.2366  
 3rd Qu.:0.1142549   3rd Qu.:0.0818748   3rd Qu.:0.2216  
 Max.   :0.9999998   Max.   :0.9050417   Max.   :0.7300
4

2 回答 2

0

您可以使用诸如马氏距离之类的多元方法去除异常值。尝试这个:

library(dplyr)

Sx <- cov(data)
D2 <- mahalanobis(data, colMeans(data), Sx)
 
# Number of outliers to remove
sum(D2 > qchisq(0.001,ncol(data),lower.tail=FALSE)
 
data2 <- filter(data,D2 <= qchisq(0.001,ncol(data),lower.tail=FALSE))
model <- betareg(formula = "y ~ x", data = data2)

data$pred <- predict(model, newdata = data)

summary(data)

输出:

#       y                   x                  pred        
# Min.   :0.0000002   Min.   :0.0000002   Min.   :0.06280  
# 1st Qu.:0.0000002   1st Qu.:0.0000002   1st Qu.:0.06280  
# Median :0.0236885   Median :0.0153738   Median :0.06638  
# Mean   :0.1277322   Mean   :0.0957217   Mean   :0.10657  
# 3rd Qu.:0.1142549   3rd Qu.:0.0818748   3rd Qu.:0.08413  
# Max.   :0.9999998   Max.   :0.9050417   Max.   :0.68639 
于 2020-12-10T19:05:08.820 回答
0

为了帮助您回答问题,我们可以训练 5 个模型并绘制预测,如下所示:

library(quantreg)
library(betareg)
mod1 <- betareg(formula = y ~ x, data = data, link = "logit")
pred1 <- predict(object = mod1, newdata = data, type = "response")
mod2 = loess(y ~ x, data = data)
pred2 = predict(mod2, data$x)
mod3 = lm(y ~ x, data = data)
mod4 <- rq(y ~ x, data = data)
x = data$x; y = data$y
y2 = y[y < 0.95]; x2 = x[y < 0.95]
mod5 = betareg(formula = y ~ x, data = data.frame(x = x2, y = y2), link = "logit")
pred5 = predict(object = mod5, newdata = data.frame(x = data$x), type = "response")
plot(data$x, data$y, xlab = "x", ylab = "y", pch = 19)
lines(sort(data$x), pred1[order(data$x)], type = "l", lty = 2)
lines(sort(data$x), pred2[order(data$x)], type = "l", col = "red2", lty = 2)
abline(mod3, col = "blue", lty = 2)
abline(mod4, col = "green", lty = 2)
lines(sort(data$x), pred5[order(data$x)], type = "l", col = "brown", lty = 2)

输出图

蓝色是线性模型,红色是黄土模型。如您所见,这些往往遵循条件均值 E[Y|X]。在绿色中,您有分位数回归,它倾向于遵循条件中位数,中位数(Y|X)。最后,黑色是 betareg 模型,棕色是去掉 y > 0.95 后的 beta 回归。为什么只删除 4 个点会产生如此戏剧性的效果?在 beta 回归中,y 的方差是高度异方差的(在 y = 0.5 时最大方差,在 y = 0 和 y = 1 时方差等于零)。y = 0.99999981 的这 4 个观察值几乎没有方差。这些是如此确定,以至于将整个预测线向上推。通过删除它们,你得到你想要的。

于 2020-12-07T21:57:17.163 回答