Rcpp矩阵:循环遍历行,一次一列

Vin*_*ent 20 r rcpp

这是我第一次尝试Rcpp,这个非常简单的问题给了我麻烦.我想使用嵌套for循环来操作矩阵的各个值,一次一列.我瞄准的脚本看起来像这样:

src <- '
    Rcpp::NumericMatrix Am(A);
    int nrows = Am.nrow();
    int ncolumns = Am.ncol();
    for (int i = 0; i < ncolumns; i++){
        for (int j = 1; j < nrows; j++){
            Am[j,i] = Am[j,i] + Am[j-1,i];
        }
    }
    return Am;
'
fun <- cxxfunction(signature(A = "numeric"), body = src, plugin="Rcpp")
fun(matrix(1,4,4))
Run Code Online (Sandbox Code Playgroud)

期望的输出是这样的:

     [,1] [,2] [,3] [,4]
[1,]    1    1    1    1
[2,]    2    2    2    2
[3,]    3    3    3    3
[4,]    4    4    4    4
Run Code Online (Sandbox Code Playgroud)

问题显然在这一行,我不知道如何引用矩阵的各个元素.

Am[j,i] = Am[j,i] + Am[j-1,i];
Run Code Online (Sandbox Code Playgroud)

如果这是一个愚蠢的新手问题,请道歉.任何提示将不胜感激!

Dir*_*tel 30

不能在单个[ ]表达式中使用多个索引.这是一个C语言的限制,即没有 C++矩阵类系统或图书馆中,我知道克服.所以请( )改用.

修复那个以及你实际没有传递src给的错误cxxfunction(),我们得到了这个:

R> src <- '
+     Rcpp::NumericMatrix Am(A);
+     int nrows = Am.nrow();
+     int ncolumns = Am.ncol();
+     for (int i = 0; i < ncolumns; i++) {
+         for (int j = 1; j < nrows; j++) {
+             Am(j,i) = Am(j,i) + Am(j-1,i);
+         }
+     }
+     return Am;
+ '
R> fun <- cxxfunction(signature(A = "numeric"), body = src, plugin="Rcpp")
R> fun(matrix(1,4,4))
     [,1] [,2] [,3] [,4]
[1,]    1    1    1    1
[2,]    2    2    2    2
[3,]    3    3    3    3
[4,]    4    4    4    4
R> 
Run Code Online (Sandbox Code Playgroud)

最后,请注意Rcpp sugar有一次处理整行或列的示例,请参阅邮件列表存档和插图.

编辑:只是为了明确,这里只使用一个循环和Rcpp糖的逐列索引:

R> src <- '
+     Rcpp::NumericMatrix Am(A);
+     int nrows = Am.nrow();
+     for (int j = 1; j < nrows; j++) {
+         Am(j,_) = Am(j,_) + Am(j-1,_);
+     }
+     return Am;
+ '
R> fun <- cxxfunction(signature(A = "numeric"), body = src, plugin="Rcpp")
R> fun(matrix(1,4,4))
     [,1] [,2] [,3] [,4]
[1,]    1    1    1    1
[2,]    2    2    2    2
[3,]    3    3    3    3
[4,]    4    4    4    4
R> 
Run Code Online (Sandbox Code Playgroud)

  • 相反,'Am(j,_)`你可以使用'Am.row(j)`. (4认同)
  • 美丽!感谢Dirk的回应.我将把未来的问题(如果有的话)引导到Rcpp-devel列表.再次感谢! (3认同)