我正在尝试在R中编写一个批量梯度下降函数,用于训练和测试数据集。目前我已经有了下面的代码。然而,当我运行它时,它只打印出最后的参数和运行的迭代次数。我希望能够存储每次迭代、测试误差,并能够可视化成本收敛过程。我不确定如何将代码添加到下面的函数中,或者在哪里添加。
GradD <- function(x, y, alpha = 0.006, epsilon = 10^-10){ iter <- 0 i <- 0 x <- cbind(rep(1, nrow(x)), x) theta <- matrix(c(1,1), ncol(x), 1) cost <- (1/(2*nrow(x)))* t(x%*% theta - y) %*% (x %*% theta - y) delta <- 1 while (delta > epsilon){ i <- i + 1 theta <- theta - (alpha / nrow(x)) * (t(x) %*% (x %*% theta - y)) cval <- (1/(2*nrow(x))) * t(x %*% theta - y) %*% (x %*% theta - y) cost <- append(cost, cval) delta <- abs(cost[i+1] - cost[i]) if((cost[i+1] - cost[i]) > 0){ print("The cost is increasing. Try reducing alpha.") return() } iter <- append(iter, i) } print(sprintf("Completed in %i iterations.", i)) return(theta)}TPredict <- function(theta, x){ x <- cbind(rep(1,nrow(x)), x) return(x %*% theta)}
EDIT我尝试创建一个列表来存储每次迭代… 但是现在运行代码时出现了错误
error.cost <- function(x, y, theta){ sum( (X %*% theta - y)^2 ) / (2*length(y))}num_iters <- 2000cost_history <- double(num_iters)theta_history <- list(num_iters)GradD <- function(x, y, alpha = 0.006, epsilon = 10^-10){ iter <- 2000 i <- 0 x <- cbind(rep(1,nrow(x)), x) theta <- matrix(c(1,1),ncol(x),1) cost <- (1/(2*nrow(x))) * t(x %*% theta - y) %*% (x %*% theta - y) delta <- 1 while(delta > epsilon){ i <- i + 1 theta <- theta - (alpha / nrow(x)) * (t(x) %*% (x %*% theta - y)) cval <- (1/(2*nrow(x))) * t(x %*% theta - y) %*% (x %*% theta - y) cost <- append(cost, cval) delta <- abs(cost[i+1] - cost[i]) cost_history[i] <- error.cost(x, y, theta) theta_history[[i]] <- theta if((cost[i+1] - cost[i]) > 0){ print("The cost is increasing. Try reducing alpha.") return() } iter <- append(iter, i) } print(sprintf("Completed in %i iterations.", i)) return(theta)}
我得到了error in nrow(x) %% theta: non-conformable arguments。如果我删除这个函数中的nrow():
error.cost <- function(x, y, theta){ sum( (x %*% theta - y)^2 ) / (2*length(y))}
然后它会打印出结果,但这些结果是错误的最终结果,而且我根本没有存储迭代
回答:
下面的代码主要使用了你的代码,仅作为捕获历史记录的概念验证。alpha值是通过尝试找到一个能通过if cost increasing
以创建超过一两个记录的值来确定的,这似乎是最初的问题:
GradD2 <- function(x, y, alpha = 0.0000056, epsilon = 10^-10) {cost <- vector(mode = 'numeric')iter <- vector(mode = 'integer')delta_hist <- vector(mode = 'numeric')i <- 0iter <- 0x <- cbind(rep(1, nrow(x)), x)theta <- matrix(c(1,1), ncol(x), 1)cost <- (1/(2*nrow(x))) * t(x%*% theta -y) %*% (x %*% theta -y)delta <- 1.0while(length(iter) < 1000 ) { #todo - change back to while(delta>epsilon)i <- i +1theta <- theta - (alpha / nrow(x)) * (t(x) %*% (x %*% theta -y))cval <- (1/(2*nrow(x))) * t(x %*% theta -y) %*% (x %*% theta -y)cost <- append(cost, cval, after = length(cost))delta <- abs(cost[i+1] - cost[i])delta_hist <- append(delta_hist, delta, after = length(delta_hist))if((cost[i+1] - cost[i]) > 0) {print('The cost is increasing. Try reducing alpha.')return(list(theta = theta,cost = cost, delta_hist = delta_hist, iter = iter))}iter <- append(iter, i, after = length(iter))}print(sprintf('Completed in %i iterations.', i))return(list(theta = theta,cost = cost, delta_hist = delta_hist, iter = iter))}
结果:
> str(GradD2_tst)List of 4 $ theta : num [1:2, 1] 1.693 0.707 $ cost : num [1:1000] 88564 87541 86587 85698 84870 ... $ delta_hist: num [1:999] 1024 954 889 828 771 ... $ iter : num [1:1000] 0 1 2 3 4 5 6 7 8 9 ...> GradD2_tst$delta_hist[994:999][1] 0.08580117 0.08580094 0.08580071 0.08580048 0.08580025 0.08580002> GradD2_tst$cost[994:999][1] 73493.72 73493.63 73493.54 73493.46 73493.37 73493.29>
除非另有说明,我的x和y是:
> x <- as.matrix(sample(1:1000, 400), ncol =1)> y <- sample(1:1000, 400)
当使用正确的while(delta > epsilon) {
与x的选择时,会出现cval的匹配错误,产生警告,现在没有了。该死。
再次尝试寻找alpha,epsilon为10^-1:
> GradD2_tst2 <- GradD2(x,y, alpha=0.006, epsilon = 10^-1)[1] "The cost is increasing. Try reducing alpha."> GradD2_tst2 <- GradD2(x,y, alpha=0.001, epsilon = 10^-1)[1] "The cost is increasing. Try reducing alpha."> GradD2_tst2 <- GradD2(x,y, alpha=0.0006, epsilon = 10^-1)[1] "The cost is increasing. Try reducing alpha."> GradD2_tst2 <- GradD2(x,y, alpha=0.00006, epsilon = 10^-1)[1] "The cost is increasing. Try reducing alpha."> GradD2_tst2 <- GradD2(x,y, alpha=0.000056, epsilon = 10^-1)[1] "The cost is increasing. Try reducing alpha."> GradD2_tst2 <- GradD2(x,y, alpha=0.000009, epsilon = 10^-1)[1] "The cost is increasing. Try reducing alpha."> GradD2_tst2 <- GradD2(x,y, alpha=0.0000056, epsilon = 10^-1)[1] "Completed in 160 iterations."> > str(GradD2_tst2)List of 4 $ theta : num [1:2, 1] 1.111 0.709 $ cost : num [1:161] 88564 87541 86587 85698 84870 ... $ delta_hist: num [1:160] 1024 954 889 828 771 ... $ iter : num [1:161] 0 1 2 3 4 5 6 7 8 9 ...> GradD2_tst2$cost[155:161][1] 73566.06 73565.96 73565.85 73565.75 73565.65 73565.55 73565.45> GradD2_tst2$delta_hist[155:161][1] 0.10493825 0.10364385 0.10243786 0.10131425 0.10026738 0.09929201 NA>
在epsilon = 10^-10 时运行一个核心需要很长时间。昨晚我的没能完成,Ctl-C。希望这对历史记录有所帮助。