Rcpp:从矩阵中消除一列和一行

Rcpp: Eliminating a column and a row from a matrix

我正在尝试创建一个接受矩阵 n x p 和索引 e 和 [=39= 的函数] 从X中消去e-th列和e-th行得到的子矩阵。我认为最简单的事情是创建一个 n -1 x p - 1 矩阵并在由 e-th[=31] 形成的 cross 周围插入角=] 其中的行和列。使用 Corner 语法,同样的方法也适用于 EigenRcpp。 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;
}

感谢您的帮助,Marco

关于这个错误,

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

我认为您(或客户端的任何人)对此无能为力,因为问题来自 SubMatrixclass definition。根据编译器,

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

m 应该是 static 成员,但它不是,因此出现错误。


虽然我无法为您提供以这种方式使用 Range 作业的一般解决方法,但对于您关于 "crossing out" [=] 的一部分的具体问题,这是一种可能的解决方法18=]:

#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 指定的列或行中,我将把 "source" 迭代器的当前值复制到 "destination" 迭代器.得到行索引为(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;
}