r - 逻辑回归梯度下降算法从 R 内置的 GLM 函数返回不同的系数

Question

我一直在尝试在 R 中编写用于逻辑回归的梯度下降算法，以便更好地理解它。在 Andrew NG 的 ML 类(class)中，他们似乎跳过了这部分，而是展示了高级优化方法。但是，我想自己重新创建梯度下降法。这是我的尝试:

###my data

X <- c(34.62366, 30.28671, 35.84741, 60.18260, 79.03274)
X <- cbind(rep(1,5),X)
y <- c(0, 0, 0, 1, 1)

###sigmoid function to calculate predicted probabilities

sigmoid <- function(z) {
  #SIGMOID Compute sigmoid function
  z <- as.matrix(z)
  g <- matrix(0,dim(z)[1],dim(z)[2])
  g <- 1 / (1 + exp(-1 * z))
  g
}


###Gradient Descent 

theta <-  c(0,0)
iterations <- 15000
alpha <- 0.02
m <- length(y)

for (i in 1:iterations) {

  theta_prev = theta

  p = dim(X)[2]

  for (j in 1:p) {
    h <- sigmoid(X %*% theta_prev)

    #sigmoid derivative
    deriv <- (t(h - y) %*% X[,j]) / m

    theta[j] = theta_prev[j] - (alpha * deriv)
  }
}

这给出了 -11.95 和 0.24 的最终系数，而在 R 中使用 GLM 函数我得到 -90.87 和 1.89。有谁知道我的代码哪里出错了？

这是 GLM 模型的代码:

X <- X[,2]
mod <- glm(y ~ X, family = 'binomial')
coef(mod)

提前致谢!

编辑:对于这个没有完美分离的较大数据集，系数之间的差异仍然存在。此外，对于包含 100 个观测值的更大数据集，差异仍然存在。

X <- c(34.62366, 30.28671, 35.84741, 60.18260, 79.03274, 45.08328, 61.10666,
   75.02475, 76.09879, 84.43282, 95.86156, 75.01366, 82.30705, 69.36459, 39.53834)
X <- cbind(rep(1,5),X)
y <- c(0, 0, 0, 1, 1, 0, 1, 1, 1, 1, 0, 0, 1, 1, 0)

对于这个稍微大一点的数据集，我的尝试返回 -18,46 和 0.15 的系数，而 R 的 GLM 返回 -4.12 和 0.07。

Answer 1

您看到的问题是由您的数据引起的。您拥有可以被多个平面分隔的数据。查看此讨论 http://r.789695.n4.nabble.com/glm-fit-quot-fitted-probabilities-numerically-0-or-1-occurred-quot-td849242.html

请注意，当我尝试 glm() 时，我会收到警告

glm.fit: glm.fit: "fitted probabilities numerically 0 or 1 occurred"

这应该会提示您有些地方不正确。基本上你会发现有无限的平面可以分开你的点(你可以说轴的左边全是 0，轴的右边全是 1)。在我引用的讨论中，链接得到了很好的解释。您自行开发的 GD 根据您的起始值(尝试!)返回不同的值，因为有几个是可以的...从

theta <-  c(20,20)

给你

> theta
[1] -18.6533438   0.3883605

在图中，您可以看到我从不同起始条件的方法中得到的三行，正如您所看到的，它们都很好地分隔了您的点...

希望对你有所帮助。最好的，翁贝托

编辑:查看您的数据后，我会说您的数据不是线性可分的(与您的初始数据所建议的相反)。 glm 给出的模型并没有真正起作用。检查 summary(mod)

Coefficients:
            Estimate Std. Error z value Pr(>|z|)  
(Intercept) -4.11494    2.32945  -1.766   0.0773 .
X[, 2]       0.06759    0.03527   1.916   0.0553 .

检查错误和 z 值...所以就我个人而言，我不会对您从 glm 获得的结果给予太多重视...并且您的代码给出的结果取决于(如预期的那样)初始值...顺便说一下，要使用您的代码和超参数获得稳定的结果，您需要更多的迭代......仍在检查一下。一旦发现更多，将尽快更新答案。

编辑 2:得到了一些东西。如果使用以下参数

theta <-  c(-4,0.05)
iterations <- 1000000
alpha <- 0.001

从你的方法中得到

> theta
[1] -4.11500250  0.06758884

从 glm 你得到

> coef(mod)
(Intercept)      X[, 2] 
-4.11493568  0.06758787 

所以相同的值(好吧，彼此非常非常接近)。现在请注意，如果您使用初始参数 c(0,0)，您仍然会得到相同的结果……学习率也是一个问题(如果您选择的太大，您的参数不会收敛)。我检查了 theta 值的行为，发现参数在两个值之间振荡，这是学习率太大的明显迹象。此外，您需要更多迭代才能收敛...

在图中截距的行为与迭代次数的关系给你一个想法......

希望对您有所帮助，翁贝托