在 R 中使用随机梯度下降编程逻辑回归

Question

我正在尝试在 R 中使用随机下降梯度来编写逻辑回归。例如，我遵循了 Andrew Ng 的示例，命名为：“ex2data1.txt”。

关键是算法工作正常，但θ估计并不完全符合我的预期。所以我试图改变整个算法来解决这个问题。然而，这对我来说几乎是不可能的。我无法检测到导致此问题的错误。因此，如果有人可以检查示例并告诉我为什么未正确计算 thetas，那将非常有用。对此，我真的非常感激。

关于编程，我没有使用 R 中实现的任何函数或矩阵计算。我只是在循环中使用加法和减法，因为我想在 hadoop 中使用代码，我不能使用矩阵微积分，甚至不能使用 R 中已经编程的函数，例如“sum”、“sqrt”等

随机梯度下降是：

Loop {
   for i = 1 to m, {
     θj := θj + α(y(i) - hθ(x(i)))(xj)(i)
  }
}`

逻辑回归：

我的代码是：

data1 <- read.table("~/ex2data1.txt", sep = ",")
names(data1) <- c("Exam1", "Exam2", "Admit")

# Sample the data for stochastic gradient decent 

ss<-data1[sample(nrow(data1),size=nrow(data1),replace=FALSE),]

x <- with(ss, matrix(cbind(1, Exam1), nrow = nrow(ss)))
y <- c(ss$Admit)
m <- nrow(x)

# startup parameters

iterations<-1
j<-vector()
alpha<-0.05
theta<-c(0,0)



#My loop

while(iterations<=10){

    coste<-c(0,0)
    suma<-0

    for(i in 1:m){

        # h<-1/(1+exp(-Q*x)

        h<-1/(1+exp((-theta)*x[i,]))

        #Cost(hQ(x),y)=y(i)*log(hQ(x))+(1-y(i))*log(1-hQ(x))

            cost<-((y[i]*log(h))+((1-y[i])*log(1-h)))

        #sum(cost) i=1 to m

            suma<-suma+cost

        #Diferences=(hQ(x(i))-y(i))*x(i)

            difference<-(h-y[i])*x[i,]  

        #sum the differences 

            coste<-coste+difference

        #calculation thetas and upgrade = Qj:= Qj - alpha* sum((h-y[i])*x[i,]*x(i))

            theta[1]<-(theta[1]-alpha*1/m*(coste[1]))
            theta[2]<-(theta[2]-alpha*1/m*(coste[2]))

    }
        #J(Q)=(-1/m)* sum ( y(i)*log(hQ(x))+(1-y(i))*log(1-hQ(x)))

            j[iterations]<-(-1/m)*suma

            iterations=iterations+1

}



#If I compare my thetas with R glm 


Call:  glm(formula = y ~ x[, 2], family = binomial("logit"), data = data1)

Coefficients:

Intercept:-4.71816 

x[, 2]  :0.08091  

我的 thetas

Intercept: 0.4624024 
 x[,2]: 1.3650234

Answer 1

我已经在 R 中为另一个 Ng 的示例集实现了一个解决方案：ex2data2.txt。这是我的代码：

sigmoid <- function(z) {
return(1/(1 + exp(-z)))
}


mapFeature <- function(X1, X2) {
degree <- 6
out <- rep(1, length(X1))
for (i in 1:degree) {
for (j in 0:i) {
out <- cbind(out, (X1^(i - j)) * (X2^j))
}
}
return(out)
}


## Cost Function
fr <- function(theta, X, y, lambda) {
m <- length(y)
return(1/m * sum(-y * log(sigmoid(X %*% theta)) - (1 - y) *
log(1 - sigmoid(X %*% theta))) + lambda/2/m * sum(theta[-1]^2))
}


## Gradient
grr <- function(theta, X, y, lambda) {
return(1/m * t(X) %*% (sigmoid(X %*% theta) - y) + lambda/m *
c(0, theta[-1]))
}

data <- read.csv("ex2data2.txt", header = F)
X = as.matrix(data[,c(1,2)])
y = data[,3]
X = mapFeature(X[,1],X[,2])
m <- nrow(X)
n <- ncol(X)
initial_theta = rep(0, n)
lambda <- 1
res <- optim(initial_theta, fr, grr, X, y, lambda,
method = "BFGS", control = list(maxit = 100000))

Answer 2

在少数情况下 * 不应该是 %*% 吗？例如。 h<-1/(1+exp((-theta) %*% x[i,])) 而不是 h<-1/(1+exp((-theta)*x[i,]))