2014-04-02 94 views
5

隨機梯度下降我想與R.隨機下降斜率比如我跟了一個名爲安德魯·吳的例子後勤迴歸程序:「ex2data1.txt」。編程Logistic迴歸與R中

的一點是,該算法正常工作,但θ驅動估計不正是我的預期。所以我試圖改變整個算法來解決這個問題。但是,我幾乎不可能。我無法檢測到導致此問題的錯誤。因此,如果有人可以檢查這個例子,並告訴我爲什麼計算不正確,那將是非常有用的。對此,我真的非常感激。

關於編程,我沒有使用任何功能在R 或矩陣計算。因爲我想使用hadoop中的代碼,所以我只是使用循環中的和和減法,我不能使用矩陣運算或甚至已經在R中編程的函數,例如「sum」,「sqrt」等。

隨機梯度下降是:

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

和Logistic迴歸:link to image

我的代碼是:

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 

我的θ驅動

Intercept: 0.4624024 
x[,2]: 1.3650234 

回答

5

我已經實現了其他Ng的例子集R中的解決方案: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)) 
+0

嗨!謝謝您的回答,但我只是用資金和減法循環中,由於我想使用Hadoop的代碼,我不能用矩陣演算,甚至被R中已經錄製諸如「求和」,「開方」的功能,」優化「等。 – user3488416

0

在少數情況下不應該是*%*%嗎?例如。 h<-1/(1+exp((-theta) %*% x[i,]))代替h<-1/(1+exp((-theta)*x[i,]))