[Rcpp-devel] How to modifying the elements of a list in the global environment?

[Davor Cubranic]

Giovanni,

Anything you pass from R into your C++ routine should be considered 
read-only. It's typically a bad idea to do in-place updates of the R 
object you passed in as an argument.

In addition, I doubt that the code you wrote (addition/subtraction of 
matrix columns) is going to be that much faster in C++ over R simply 
because R already does those operations with optimized C-level code. 
Maybe not as you wrote it, but you can use available R's functionality a 
lot more effectively with something like this:

   sweep(array(x, dim=c(N, M, K)), 2:3, mu)

This should do exactly what you're trying to do with your C++ function, 
except storing it as a 3D array instead of a list of 2D matrices. (It's 
more efficient and works better with everything else in R.)

If it's the large number of iterations in your Gibbs sampler that you're 
worried about, then write the entire algorithm in C++ and return the 
final result.

Davor


On 12-10-18 07:28 PM, Giovanni Petris wrote:
> Thanks Dirk!
>
> By the way, I just figured out a working solution for what I wanted to achieve - 'code_f' below.
>
> Best,
> Giovanni
> =====
>
> ### Small reproducible example
> K <- 2; N <- 5; M <- 3
> x <- matrix(rnorm(N * M), N, M)
> mu <- matrix(rnorm(M * K, mean = 1 : K), M, K, TRUE)
>
> ### What I would do in R
> f_R <- function() {
>      ## center observations in 'x' using columns of 'mu'
>      w <- lapply(1L:K, function(i) x - rep(mu[, i], each = N))
>      ## eventually more stuff here...
>      w
> }
>
> ### For the Rcpp part, set aside working space in the global environment
> workspace <- lapply(1:K, function(k) matrix(0.0, N, M))
>
> code_f <- '
>      Rcpp::NumericMatrix cpp_x(x);
>      Rcpp::NumericMatrix cpp_mu(mu);
>      int K = cpp_mu.ncol();
>      int M = cpp_mu.nrow();
>      int N = cpp_x.nrow();
>      Environment glob = Environment::global_env();
>      List cpp_work(glob.get("workspace"));
>      NumericMatrix tmpMatrix;
>
>      for (int k=0; k < K; k++) {
>          tmpMatrix = as<SEXP>(cpp_work[k]);
>          for (int i=0; i < M; i++)
>              tmpMatrix(_, i) = cpp_x(_, i) - cpp_mu(i, k);
>      }
>
>      return R_NilValue;
> '
> ## library(Rcpp); library(inline)
> f_cpp <- cxxfunction(signature(x = "numeric", mu = "numeric"), code_f, plugin = "Rcpp")
>
> ## try it
> f_cpp(x, mu)
> all.equal(workspace, f_R())
>
>
> ________________________________________
> From: Dirk Eddelbuettel [edd at debian.org]
> Sent: Thursday, October 18, 2012 9:15 PM
> To: Giovanni Petris
> Cc: Dirk Eddelbuettel; rcpp-devel at lists.r-forge.r-project.org
> Subject: RE: [Rcpp-devel] How to modifying the elements of a list in the        global  environment?
>
> On 19 October 2012 at 01:32, Giovanni Petris wrote:
> | Hi Dirk,
> |
> | Thank you for the quick reply. I will look for more examples on the net.
> |
> | About your suggestion of allocating scrap space inside the C++ routine, am I wrong to think that when the matrices are large and the function is called repeatedly within a Gibbs sampler loop, this is not a very efficient approach?
>
> Yes, if you want to use it several times you can pass it around for reuse.
> As I mentioned, a reference is just a pointer and hence unrelated to the size
> of the object.
>
> Dirk
>
> | Thanks again!
> |
> | Best,
> | Giovanni
> |
> | ________________________________________
> | From: Dirk Eddelbuettel [edd at debian.org]
> | Sent: Thursday, October 18, 2012 8:06 PM
> | To: Dirk Eddelbuettel
> | Cc: Giovanni Petris; rcpp-devel at lists.r-forge.r-project.org
> | Subject: Re: [Rcpp-devel] How to modifying the elements of a list in the        global  environment?
> |
> | Giovanni,
> |
> | On 18 October 2012 at 17:46, Dirk Eddelbuettel wrote:
> | | | > code_f <- '
> | | | +     Rcpp::NumericMatrix cpp_x(x);
> | | | +     Rcpp::NumericMatrix cpp_mu(mu);
> | | | +     int K = cpp_mu.ncol();
> | | | +     int M = cpp_mu.nrow();
> | | | +     int N = cpp_x.nrow();
> | | | +     Environment glob = Environment::global_env();
> | | | +     List cpp_work(glob.get("workspace"));
> | | | +
> | | | +     for (int k=0; k < K; k++) {
> | | | +         cpp_work[k] = clone(cpp_x);
> | | | +         for (int i=0; i < M; i++)
> | | | +             cpp_work(_, i) = cpp_work(_, i) - cpp_mu(i, k);
> |
> | That can't work.  You have the types confused.
> |
> | You could just pass the list object down, and then pick from it.  Working
> | examples for that are in the list archives, in the RcppExamples package, and
> | in other places on the net.
> |
> | Or, as I suggested, just allocate the scrap space inside the C++ routine.
> |
> | Hth,  Dirk
> |
> | --
> | Dirk Eddelbuettel | edd at debian.org | http://dirk.eddelbuettel.com
>
> --
> Dirk Eddelbuettel | edd at debian.org | http://dirk.eddelbuettel.com
> _______________________________________________
> Rcpp-devel mailing list
> Rcpp-devel at lists.r-forge.r-project.org
> https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
>