哪个链接函数可用于响应为百分比(0 - 100%)的 glm?

机器算法验证 r 广义线性模型 生态 链接功能 贝塔回归
2022-03-28 23:00:27

我正在考虑建立一个模型(glm),其中响应变量(y)是特定区域内植物物种的覆盖率(百分比),取决于环境变量。但是,我不认为它是二项式链接,因为输出不是介于 0 和 1 之间的概率。您建议将什么用作链接函数?

1个回答

好问题!

您可以通过转换覆盖变量来开始建模工作,使其表示为比例而不是百分比(即,只需将原始覆盖变量的值除以 100)。

一旦你得到转换后的变量,看看你有多少 0 和 1 值,如果有的话,除了落在 (0,1) 区间的值。

如果您没有任何 0 或 1 值,则可以使用beta 回归建模beta 回归模型不是 glm 模型,但它与 glm 模型有相似之处。beta回归模型具有logit链接函数;这意味着您正在将 logit 转换的平均覆盖建模为各种预测变量的函数。(回想一下,现在覆盖率表示为一个比例。)

如果您有大量 0 值但没有 1 值,则需要使用零膨胀 beta 回归模型进行建模。

如果您有大量 1 值但没有 0 值,则需要使用一次膨胀的 beta 回归模型

如果您有大量的 0 值和大量的 1 值,则需要使用零一一膨胀的 beta 回归模型

你用R吗?如果是,那么gamlss包是实现这些模型的最佳选择 - 例如,请参阅本文档的第 107 页:

https://www.gamlss.com/wp-content/uploads/2013/01/book-2010-Athens1.pdf

Beta 回归模型的膨胀版本具有多个链接,因为它们同时对多个参数进行建模。其中一些链接可以是 logit 链接,一些可以是 log 链接等。在gamlss中,这些模型的系列分布可以是以下之一:

  1. BE 用于 beta 回归;

  2. BEINF0 用于零膨胀 beta 回归;

  3. BEINF1 用于单膨胀 beta 回归;

  4. BEINF 用于零和一膨胀 beta 回归。

附录

betareg() 和 gamlss() 都能够将 beta 回归模型拟合到您的数据,然后使用该模型进行预测。下面的 R 代码提供了如何执行此操作的示例。请注意,来自同一数据集的两个函数生成的模型摘要之间的输出存在一些差异。

#===========================================================
# Load Gasoline data into R 
#===========================================================

library(betareg)

data("GasolineYield", package = "betareg")

#===========================================================
# Fit beta regression model to Gasoline data using betareg()
#===========================================================

model_betareg <- betareg(yield ~ temp, 
                         link = "logit", 
                         link.phi = "log", 
                         data = GasolineYield)

 summary(model_betareg)

#===========================================================
# Fit beta regression model to Gasoline data using gamlss()
#===========================================================

library(gamlss) 

model_gamlss <- gamlss(yield ~ temp, 
                       family = BE(mu.link = "logit", 
                                   sigma.link = "log"), 
                       data = GasolineYield)

summary(model_gamlss)

#===========================================================
# Predict from model_betareg for a new temp value
#===========================================================

predict_betareg <- predict(model_betareg, 
                           newdata = data.frame(temp = 349), 
                           type="response")
predict_betareg

#===========================================================
# Predict from model_gamlss for a new temp value
#===========================================================

predict_gamlss <- predict(model_gamlss, 
                          newdata = data.frame(temp = 349), 
                          type="response")
predict_gamlss

model_betareg预测的屈服值为0.2046289,model_gamlss产生的屈服值为0.204634。正如预期的那样,这些非常接近。

请注意,预测是在响应变量产量的尺度上执行的,因此预测值表示为区间 (0,1) 中的比例。

newdata 必须是一个数据框,它列出了包含在 beta 回归模型右侧的所有预测变量的值。这些值的类型必须与 R 通过 str() 命令看到的预测变量的类型相匹配:

str(Gasoline)

这是此 str() 命令的输出:

str(Gasoline)
Classes ‘nfnGroupedData’, ‘nfGroupedData’, ‘groupedData’ and 'data.frame':      
32 obs. of  6 variables:
$ yield   : num  6.9 14.4 7.4 8.5 8 2.8 5 12.2 10 15.2 ...
$ endpoint: num  235 307 212 365 218 235 285 205 267 300 ...
$ Sample  : Ord.factor w/ 10 levels "1"<"2"<"7"<"9"<..: 8 9 5 1 3 4 7 10 6 8 
...
$ API     : num  38.4 40.3 40 31.8 40.8 41.3 38.1 50.8 32.2 38.4 ...
$ vapor   : num  6.1 4.8 6.1 0.2 3.5 1.8 1.2 8.6 5.2 6.1 ...
$ ASTM    : num  220 231 217 316 210 267 274 190 236 220 ...

如果 temp 被列为数据集中的一个因素(例如,“High”、“Medium”、“Low”),您可以在 newdata 中指定它:

predict_gamlss <- predict(model_gamlss, 
                          newdata = data.frame(temp = c("High","Low")), 
                          type="response")
predict_gamlss