逻辑回归如何使用二项分布?

机器算法验证 物流 二项分布 直觉
2022-02-04 02:26:40

我试图了解逻辑回归如何使用二项分布。

假设我正在研究鸟类的筑巢成功率。筑巢成功的概率是 0.6。使用二项分布,我可以计算给定 n 次试验(研究的巢数)的 r 次成功的概率。

但是在建模环境中如何使用二项分布呢?假设我想知道平均每日温度如何影响筑巢成功,我使用逻辑回归来探讨这个问题。

在我所描述的上下文中,逻辑回归如何使用二项分布?

我正在寻找一个直观的答案,因此是一个没有方程式的答案!我认为只有在直观的层面上实现了理解后,方程才有用。

3个回答

没有方程?哎呀。让我们来看看:

逻辑回归模型实际上是一个模型p二项分布的参数;使用连续预测器,每个点都可以有自己的分布。(在观测值为 0-1 的情况下,我们处理伯努利特例;这是一种常见情况。)

n是给定的,不是建模的。所以结果是,有一个与pi和一个已知的ni,我们可以根据预测变量对二项式数据进行建模,该预测变量通过其模型描述均值(和方差)p. 该模型可能通过最大似然估计来拟合,但由于其特殊形式(指数族),ML 相对“不错”。

因为逻辑链接对于二项式家族来说是规范的,所以它甚至更好,因为足够的统计数据具有非常简单的形式——这使得处理大样本甚至开发“在线”算法变得方便。

当然,p,作为概率,位于 0 和 1 之间。这自然意味着当我们根据其他变量为其编写模型时,该模型不应崩溃这些限制,因此自变量变得足够大或足够小,关系必须弯曲以保持在界限内。

对于逻辑回归,该曲线(链接函数)是逻辑函数。其他功能是可能的,并且许多包实现了几个(如果我没记错的话,R 在其功能中内置了三个合适的glm功能)。


在这篇文章的制作过程中没有损害任何平等符号。

假设您在不同的日平均温度下观察到几个巢穴t. 概率如何π(t)筑巢成功取决于温度t? (如果巢是独立的,在温度下成功的巢数t然后是二项式分布n等于观察到的巢数和成功概率π(t).)

逻辑回归是一种方法(使用逻辑函数),通过拉伸和移动逻辑曲线将成功概率指定为温度的函数,需要从数据中估计拉伸和移动的量。

您的模型假设一个巢穴的成功可以被视为一场赌博:上帝掷硬币,硬币的两面标有“成功”和“失败”。一个巢的翻转结果与任何其他巢的翻转结果无关。

不过,这些鸟儿确实对它们有好处:与其他温度相比,硬币在某些温度下可能非常有利于成功。因此,当您有机会在给定温度下观察巢穴时,成功的次数等于成功翻转同一硬币的次数——即该温度下的硬币。相应的二项分布描述了成功的机会。也就是说,它通过嵌套的数量建立了零成功、一、二、...等等的概率。

温度与上帝如何装载硬币之间的关系的一个合理估计是由在该温度下观察到的成功比例给出的。这是最大似然估计 (MLE)。

例如,假设您观察7在温度下筑巢10学位和3那些巢是成功的。MLE 是3/7. 也就是说,我们估计上帝的硬币有3/7显示成功的机会。相应的二项分布绘制在图的第一行(见下文)的标题“10 度”下。它表示垂直线段高度的机会。红色部分对应于观察值3成功。

您的数据中的温度必须有所不同。作为一个运行示例,让我们假设在温度5,10,15,20你观察到的度数0,3,2,3之间的成功2,7,5,3巢。该数据集由图中“拟合”面板中的灰色圆圈绘制。圆圈的高度代表它的成功率。圆圈面积与嵌套数量成正比(从而强调具有更多嵌套的数据)。

该图的第一行显示了四个观察到的温度中的每个温度下的 MLE。 “适合”面板中的红色曲线描绘了硬币的加载方式,具体取决于温度。通过构造,此跟踪通过每个数据点。(它在中间温度下的作用是未知的;我粗略地连接了这些值以强调这一点。)

这种“饱和”模型并不是很有用,正是因为它没有给我们提供依据来估计上帝将如何在中间温度下装载硬币。为此,我们需要假设存在某种将硬币装载量与温度相关联的“趋势”曲线。

数字

该图的底行符合这样的趋势。 趋势的作用是有限的:当绘制在适当的(“对数赔率”)坐标中时,如左侧的“Logit 响应”面板所示,它只能沿着一条直线。任何这样的直线都决定了硬币在所有温度下的装载量,如“适合”面板中的相应曲线所示。反过来,该载荷决定了所有温度下的二项式分布。底行绘制了观察巢穴的温度分布。(黑色虚线标记了分布的预期值,有助于相当精确地识别它们。您在图中的第一行看不到这些线,因为它们与红色线段重合。)

现在必须做出权衡:这条线可能会靠近某些数据点,但会偏离其他点。这会导致相应的二项式分布为大多数观测值分配比以前更低的概率。您可以在 10 度和 15 度处清楚地看到这一点:观测值的概率不是可能的最高概率,也不接近上排分配的值。

逻辑回归滑动和摆动可能的线(在“Logit 响应”面板使用的坐标系中),将它们的高度转换为二项式概率(“拟合”面板),评估分配给观察的机会(四个右侧面板),并选择能够提供这些机会的最佳组合的线路。

什么是“最好”? 简单地说,所有数据的组合概率尽可能大。通过这种方式,不允许单个概率(红色部分)非常小,但通常大多数概率不会像饱和模型中那样高。

这是逻辑回归搜索的一次迭代,其中线向下旋转:

图 2

首先,请注意保持不变: “拟合”散点图中的灰点是固定的,因为它们代表数据。同样,四个二项式图中红色线段的值范围水平位置也是固定的,因为它们也代表数据。然而,这条新线路以完全不同的方式加载硬币。这样做,它改变了四个二项分布(灰色部分)。例如,它在温度为10度,对应于 4 到 6 次成功概率最高的分布。这条线实际上在拟合数据方面做得很好15度,但拟合其他数据的工作很糟糕。(在 5 度和 20 度时,分配给数据的二项式概率非常小,您甚至看不到红色部分。)总体而言,这比第一个图中显示的拟合度差得多。


我希望这个讨论能帮助你形成二项式概率随着线的变化而变化的心理形象,同时保持数据不变。逻辑回归拟合线试图使这些红条总体上尽可能高。因此,逻辑回归与二项分布族之间的关系是深刻而密切的。


附录:R生成数字的代码

#
# Create example data.
#
X <- data.frame(temperature=c(5,10,15,20),
                nests=c(2,7,5,3),
                successes=c(0,3,2,3))
#
# A function to plot a Binomial(n,p) distribution and highlight the value `k0`.
#
plot.binom <- function(n, p, k0, highlight="#f02020", ...) {
  plot(0:n, dbinom(0:n, n, p), type="h", yaxt="n",
       xlab="Trials", ylab="Probability", ...)
  abline(v = p*n, lty=3, lwd=2)
  if(!missing(k0)) lines(rep(k0,2), c(0, dbinom(k0,n,p)), lwd=2, col=highlight)
}
#
# A function to convert from probability to log odds.
#
logit <- function(p) log(p) - log(1-p)
#
# Fit a saturated model, then the intended model.
#
# Ordinarily the formula for the saturated model would be in the form
# `... ~ factor(temperature)`, but the following method makes it possible to  
# plot the predicted values in a visually effective way.
#
fit.0 <- glm(cbind(successes, nests-successes) ~ factor(round(temperature/5)), 
             data=X, family=binomial)
summary(fit.0)

fit <- glm(cbind(successes, nests-successes) ~ temperature, 
           data=X, family=binomial)
summary(fit)
#
# Plot both fits, one per row.
#
lfits <- list(fit.0, fit)
par.old <- par(mfrow=c(length(lfits), nrow(X)+2))
for (fit in lfits) {
  #
  # Construct arrays of plotting points.
  #
  X$p.hat <- predict(fit, type="response")
  Y <- data.frame(temperature = seq(min(X$temperature), max(X$temperature), 
                                    length.out=101))
  Y$p.hat <- predict(fit, type="response", newdata=Y)  # Probability
  Y$lambda.hat <- predict(fit, type="link", newdata=Y) # Log odds
  #
  # Plot the fit in terms of log odds.
  #
  with(Y, plot(temperature, lambda.hat, type="n", 
               yaxt="n", bty="n", main="Logit Response",
               ylab=expression(hat(lambda))))
  if (isTRUE(diff(range(Y$lambda.hat)) < 6)) {
    # Draw gridlines and y-axis labels
    p <- c( .10, .25, .5, .75, .9)
    q <- logit(p)
    suppressWarnings(rug(q, side=2))
    abline(h=q, col="#d0d0d0")
    mtext(signif(p, 2), at=q, side=2, cex=0.6)
  }
  with(Y, lines(temperature, lambda.hat, lwd=2, col="#f02020"))
  #
  # Plot the data and the fit in terms of probability.
  #
  with(X, plot(temperature, successes/nests, ylim=0:1,
               cex=sqrt(nests), pch=21, bg="Gray",
               main="Fit"))
  with(Y, lines(temperature, p.hat, col="#f02020", lwd=2))
  #
  # Plot the Binomial distributions associated with each row of the data.
  #
  apply(X, 1, function(x) plot.binom(x[2], x[4], x[3], bty="n", lwd=2, col="Gray",
                                     main=paste(x[1], "Degrees")))
}
par(mfrow=par.old)