RCPP：从矩阵消除列和行

Question

我想创建一个函数，矩阵，NXP和索引，ē，并返回消除Ë获得的子矩阵第柱和从X的第e个行我认为最简单的事情是创造一个n -1个XP - 1个矩阵，并插入由第e个形成的横周围的角落行和列中。 EigenRcpp使用Corner语法可以很好地工作。看起来Rcpp不喜欢赋值给范围。我得到了以下错误消息：

error: non-static reference member 'Rcpp::SubMatrix<14>::MATRIX& Rcpp::SubMatrix<14>::m', can't use default assignment operator

我下面的代码复制。我的问题： 如何为矩阵块分配值？ 有没有更好的方法来做到这一点？

#include <Rcpp.h> 
using namespace Rcpp; 

NumericMatrix elimMat(NumericMatrix X, int e){ 

int p = X.cols() - 1; 
int n = X.rows() - 1; 
int f = e - 1; 
NumericMatrix M = NumericMatrix(n, p); 

M(Range(0, f - 1), Range(0, f - 1)) = X(Range(0, f - 1), Range(0, f - 1)); //TopLeft same 
M(Range(0, f - 1), Range(f, p - 1)) = X(Range(0, f - 1), Range(f + 1, p)); //TopRight 
M(Range(f, n - 1), Range(0, f - 1)) = X(Range(f + 1, n), Range(0, f - 1)); //BottomLeft 
M(Range(f, n - 1), Range(f, p - 1)) = X(Range(f + 1, n), Range(f + 1, p)); //BottomRight 

return M; 
} 

感谢您的帮助，马可

Answer 1

关于此错误，

error: non-static reference member ‘Rcpp::SubMatrix<14>::MATRIX& Rcpp::SubMatrix<14>::m’, can’t use default assignment operator

我不认为这有什么，你（或客户端上的任何人），因为可以做这个该问题来自SubMatrix的class definition。据编译器，

private: 
    MATRIX& m ; 
    vec_iterator iter ; 
    int m_nr, nc, nr ; 

m应该是一个static成员，它不是，因此错误。

---

虽然我不能给你一个一般的工作，周围以这种方式使用Range分配，这里是一个可行的方法，以您的具体问题有关“划掉”一NumericMatrix的部分：

#include <Rcpp.h> 

// [[Rcpp::export]] 
Rcpp::NumericMatrix crossout(Rcpp::NumericMatrix X, int e) { 
    Rcpp::NumericMatrix result = Rcpp::no_init_matrix(X.nrow() - 1, X.ncol() - 1); 
    Rcpp::NumericMatrix::iterator src = X.begin(), dst = result.begin(); 

    while (src != X.end()) { 
    if (((src - X.begin()) % X.nrow()) != e && ((src - X.begin())/X.nrow()) != e) { 
     *dst = *src; 
     ++dst; 
    } 
    ++src; 
    } 

    return result; 
} 

---

/*** R 

M <- matrix(1:9 + .5, nrow = 3) 
R> all.equal(M[c(1, 3), c(1, 3)], crossout(M, 1)) 
#[1] TRUE 

R> crossout(M, 1) 
#  [,1] [,2] 
#[1,] 1.5 7.5 
#[2,] 3.5 9.5 

M2 <- matrix(1:20 + .5, nrow = 4) 
R> all.equal(M2[c(1:2, 4), c(1:2, 4:5)], crossout(M2, 2)) 
#[1] TRUE 

R> crossout(M2, 2) 
#  [,1] [,2] [,3] [,4] 
#[1,] 1.5 5.5 13.5 17.5 
#[2,] 2.5 6.5 14.5 18.5 
#[3,] 4.5 8.5 16.5 20.5 

R> any(unique(c(M2[3,], M2[,3])) %in% crossout(M2, 2)) 
#[1] FALSE 

*/ 

---

上面，我复制 “源” 伊特拉的当前值如果源迭代器的位置不落在e指定的列或行内，则返回“目标”迭代器。行索引获得为(src - X.begin()) % X.nrow()，列索引为(src - X.begin())/X.nrow()。

---

编辑： 对于下方的评论，这将是非常简单的，以适应非方阵/不同的行和列维度指数：

#include <Rcpp.h> 

// [[Rcpp::export]] 
Rcpp::NumericMatrix crossout2(Rcpp::NumericMatrix X, int r, int c) { 
    if (r >= X.nrow() || c >= X.ncol()) { 
    Rf_error("Invalid row or column index.\n"); 
    } 

    Rcpp::NumericMatrix result = Rcpp::no_init_matrix(X.nrow() - 1, X.ncol() - 1); 
    Rcpp::NumericMatrix::iterator src = X.begin(), dst = result.begin(); 

    while (src != X.end()) { 
    if (((src - X.begin()) % X.nrow()) != r && ((src - X.begin())/X.nrow()) != c) { 
     *dst = *src; 
     ++dst; 
    } 
    ++src; 
    } 
    return result; 
} 

---