<div dir="ltr">Hi there,<div><br></div><div>I am new to RcppArmadillo and very excited about the sparse matrix functionality it offers. I have a large optim that I wish to break down into stochastic descent in Rcpp to speed things up.</div>
<div><br></div><div>To start small, I used a simple logisitic regression function that I have working in R -</div><div><br></div><div><div>Norm <- function(w1, w2){</div><div> sqrt(sum((w2 - w1)^2))</div><div>}</div><div>
<br></div><div>xentropy <- function(x, y, w) {</div><div> - ((y * x) / (1 + exp(y * t(w) %*% x)) )</div><div>}</div><div><br></div><div>stochDesc <- function(param, eta = 0.01, maxit = 10, reltol = .01, x, y) {</div>
<div> # Initialize</div><div> w = param</div><div> iter = 0</div><div> rel_err = reltol + 1</div><div> x_dash = cbind(rep(1, nrow(x)), x) # Add bias</div><div> </div><div> while(iter < maxit & rel_err > reltol ) { </div>
<div> iter = iter + 1</div><div> </div><div> prev_w <- w</div><div> for(curidx in sample(nrow(x))) { # Random shuffle</div><div> grad <- xentropy(x_dash[curidx, ], y[curidx], w)</div><div> w <- w - eta * grad # Update weights</div>
<div> }</div><div> rel_err = Norm(w, prev_w)<br></div><div> # could also use: </div><div> # rel_err = base::norm(matrix(prev_w - w, nrow = 1), type = "F")</div><div> cat(paste("Epoch ", iter, " Relative cost:", rel_err, "\n"))</div>
<div> }</div><div> list(w = w, iter = iter)</div><div>}</div></div><div><br></div><div>This does every operation that my more complex gradient descent will do, but it is pretty straightforward so I started converting it to Rcpp -</div>
<div><br></div><div><div>#include <RcppArmadillo.h></div><div>#include <iostream></div><div>using namespace std;</div><div>using namespace Rcpp ;</div><div>// [[Rcpp::depends(RcppArmadillo)]]</div><div><br></div>
<div>double Norm(NumericVector w1, NumericVector w2) { </div><div> // Compute the Forbius Norm</div><div> double op = pow(sum(pow(w2 - w1, 2)), .5); </div><div> return(op);</div><div>}</div><div><br></div><div><br></div>
<div>NumericVector xentropy(NumericVector x, NumericVector y, NumericVector w) { </div><div> // Return the cross entropy</div><div> // NumericVector op = - ((y * x) / (1 + exp(y * trans(w) %*% x)) ); // How to transpose?</div>
<div> NumericVector op = w; // Not correct</div><div> return (op);</div><div>}</div><div><br></div><div><br></div><div>// [[Rcpp::export]]</div><div>NumericVector stochDescCpp(NumericVector param, NumericMatrix x, NumericVector y,</div>
<div> double eta = 0.01, int maxit = 10, double reltol = .01) {</div><div> // Run stochastic Descent to compute logistic regression</div><div> // Initialize</div><div> NumericVector w (param);</div>
<div> NumericVector prev_w (w);</div><div> int epoch = 0;</div><div> NumericVector grad (param.length());</div><div> </div><div> double rel_err = reltol + 1;</div><div> // Don't know how to add bias, will do that in R before passing it here.</div>
<div> //NumericMatrix x_dash(x.nrow(), x.ncol() + 1, rep(1, x.nrow()), x.begin());</div><div> </div><div> // Loop through for each epoch</div><div> while(epoch < maxit && rel_err > reltol) { </div><div> epoch++;</div>
<div> </div><div> prev_w = w;</div><div> for(int curidx = 0; curidx < x.nrow(); ++curidx) { // How to random shuffle?</div><div> grad = xentropy(x(curidx, _), y[curidx], w);</div><div> w = w - eta * grad; // Update weights</div>
<div> }</div><div> rel_err = Norm(w, prev_w);</div><div> cout << "Epoch " << epoch << " Relative cost:" << rel_err << endl;</div><div> }</div><div> return(w);</div>
<div>}</div></div><div><br></div><div>I am running into many issues. Here they are ranked in order of importance) -</div><div>1) When I run this (incorrect) program -</div><div>stochDescCpp(param = rep(0,3), x = X, y = Y, maxit = 100)<br>
</div><div><br></div><div>I get -</div><div><div>Error in .Primitive(".Call")(<pointer: 0x0000000069d42560>, param, x, : </div><div> negative length vectors are not allowed</div></div><div><br></div><div>
All the inputs are valid and work fine with the R version. So I don't know where this is coming from? I put a bunch of couts and I know that it is happening inside the for loop but don't know why. In general, what is the best way to debug code written this way? I am using RStudio.</div>
<div><br></div><div>2) I don't know how to transpose the matrix in rCpp. I read that I could that using rCppArmadillo, but then I don't know how to convert from NumericMatrix to arma matrix. I don't want to copy as this operation will happen 1000x of times, so this has to be fast.</div>
<div><br></div><div>3) I don't know how to randomly shuffle the rows. I tried <span class="" style="background-color:transparent;color:inherit;font-family:Menlo,Monaco,'Courier New',monospace;font-size:12.025px;line-height:18px;white-space:pre-wrap">RcppArmadillo</span><span class="" style="background-color:transparent;color:inherit;font-family:Menlo,Monaco,'Courier New',monospace;font-size:12.025px;line-height:18px;white-space:pre-wrap;font-weight:bold">::</span><span class="" style="background-color:transparent;color:inherit;font-family:Menlo,Monaco,'Courier New',monospace;font-size:12.025px;line-height:18px;white-space:pre-wrap">sample but haven't gotten the output to work. E.g. do I access elements as x(rowOrd[curIdx],_)?</span></div>
<div><span class="" style="background-color:transparent;color:inherit;font-family:Menlo,Monaco,'Courier New',monospace;font-size:12.025px;line-height:18px;white-space:pre-wrap"><br></span></div><div><span class="" style="background-color:transparent;color:inherit;font-family:Menlo,Monaco,'Courier New',monospace;font-size:12.025px;line-height:18px;white-space:pre-wrap">4) (minor issue) How do I add the bias unit? In other words, add a new column to the input matrix filled with all 1s. Do I have to create a new matrix, copy it column by column? this is not a dealbreaker since I can do this in R before I send in the inputs.</span></div>
<div><br></div><div>Sorry if my questions are too basic and don't seem well researched. I am trying to go for a practical example. </div><div><br></div><div>Any help will be greatly appreciated.</div><div><br></div><div>
Sincerely,</div><div>Saurabh</div></div>