迭代非零元素本徵

Question

我傳遞小稀疏矩陣（用於測試），以一個C++從R.函數矩陣屬於類dgCMatrix，如下所示：

5 x 5 sparse Matrix of class "dgCMatrix" 

[1,] . . . . . 
[2,] 1 1 . . . 
[3,] . . . . . 
[4,] . . 1 . . 
[5,] . 1 . . . 

我迭代這個矩陣在文檔here中提到。 我的函數輸出迭代器的值和行索引，列索引。

C++函數被定義如下：

#include <RcppEigen.h> 
// [[Rcpp::depends(RcppEigen)]] 
using Eigen::MappedSparseMatrix; 
using Eigen::SparseMatrix; 
using Eigen::VectorXi; 
using Eigen::Map; 
using namespace Rcpp; 
using namespace std; 

// [[Rcpp::export]] 
void createRec(RObject sparse_mat, IntegerVector sparse_vec) { 

const MappedSparseMatrix<int> spmat(as<MappedSparseMatrix<int> >(sparse_mat)); 
long int nrow = spmat.rows(); 
long int ncol = spmat.cols(); 
NumericVector sim(nrow); 

for(int k=0;k<spmat.outerSize();k++){ 
    for(SparseMatrix<int,Eigen::ColMajor>::InnerIterator it(spmat,k);it;++it){ 
     cout<<"k="<<k<<endl; 
     cout<<"value="<<it.value()<<endl; 
     cout<<"it.row="<<it.row()<<endl; 
     cout<<"it.col="<<it.col()<<endl; 
     cout<<"index="<<it.index()<<endl; 
    } 
} 
} 

對於上面的結果如下給出的矩陣被印刷：

k=0 
value=156148016 
it.row=66211520 
it.col=0 
index=66211520 
k=1 
value=0 
it.row=0 
it.col=1 
index=0 
k=1 
value=1 
it.row=4 
it.col=1 
index=4 
k=2 
value=1 
it.row=3 
it.col=2 
index=3 

1）不限於對應於k = 0的值的解釋？這些可能是由於以錯誤的方式傳遞矩陣？

2.）k正在迭代outerSize，它等於5，爲什麼它不迭代k = 3,4？考慮到它是一個sparseMatrix，這個行爲是來自迭代器的。

Answer 1

無論何時您看到非常大的數字，如156148016或66211520，您可能有未定義的行爲（UB）或未正確初始化值。在這種情況下，它是後者。具體而言，dgCMatrix class的基礎類型是double而非int。

dgCMatrix類是壓縮的稀疏列式格式中的一類稀疏數字矩陣。在這個實現中，列中的非零元素被排序爲增加的行順序。 dgCMatrix是Matrix包中稀疏數字矩陣的「標準」類。

因此，當你正在試圖創建一個映射到底層RObject的存儲位置存在要重新重新創建對象在被請求不同類型所需的附加步驟。在添加const條款之後，我願意打賭這些條目是如預期的，因爲編譯器可能會記住中間對象。

所以，變化如下：

MappedSparseMatrix<int> spmat(as<MappedSparseMatrix<int> >(sparse_mat)); 

到：

MappedSparseMatrix<double> spmat(as<MappedSparseMatrix<double> >(sparse_mat)); 

應該足夠了。

---

的linked example使用SparseMatrix矩陣，在這裏你使用的是MappedSparseMatrix但不要設置爲第二循環的適當MappedSparseMatrix::InnerIterator。

因此，我們有：

for(SparseMatrix<int,Eigen::ColMajor>::InnerIterator it(spmat,k);it;++it){ 

要：

for(MappedSparseMatrix<double>::InnerIterator it(spmat,k);it;++it){ 

另外請注意，不需要SparseMatrix<int, Eigen::ColMajor>::InnerIterator內使用Eigen::ColMajor爲that is the default initialization。所以，我已經刪除了這個聲明。

---

關於你的第二個問題，在迭代k。

k確實迭代在k=3,4上，但在這些列中沒有元素。因此，輸出k的內循環確實是而不是被調用。

如果我們在外部和內部循環中放置兩個k聲明性輸出語句，很容易看出。

例如

for(int k = 0; k < spmat.outerSize(); ++k) { 
    Rcpp::Rcout << "Overall k = " << k << std::endl << std::endl; 
    for(MappedSparseMatrix<double>::InnerIterator it(spmat,k); it; ++it) { 
    Rcpp::Rcout << "Inner k = " << k << std::endl; 
    } 
} 

---

避免using namespace std;

添加在命名空間中有時也有意想不到的後果，尤其是一個大std。

---

以從上面稍點簡化你的榜樣，我們下面光禿禿的骨頭工作示例：

#include <RcppEigen.h> 
// [[Rcpp::depends(RcppEigen)]] 
using Eigen::MappedSparseMatrix; 
using Eigen::SparseMatrix; 
using Eigen::VectorXi; 
using Eigen::Map; 

// [[Rcpp::export]] 
void createRec(Rcpp::RObject sparse_mat) { 

    MappedSparseMatrix<double> spmat(Rcpp::as<MappedSparseMatrix<double> >(sparse_mat)); 

    long int nrow = spmat.rows(); 
    Rcpp::NumericVector sim(nrow); 

    for(int k = 0; k < spmat.outerSize(); ++k) { 
    Rcpp::Rcout << "Overall k = " << k << std::endl << std::endl; 
    for(MappedSparseMatrix<double>::InnerIterator it(spmat,k); it; ++it) { 
     Rcpp::Rcout << "Inner k = " << k << std::endl 
        << "value = " << it.value() << std::endl 
        << "it.row = " << it.row() << std::endl 
        << "it.col = " << it.col() << std::endl 
        << "index = " << it.index() << std::endl; 
    } 


    } 
} 

/***R 

# Setup values 
id_row = c(2, 2, 4, 5) 
id_col = c(1, 2, 3, 2) 
vals = rep(1,4) 

# Make the matrix 
x = sparseMatrix(id_row, id_col, x = vals, dims = c(5, 5)) 

# Test the function 
createRec(x) 
*/ 

輸出：

Overall k = 0 

Inner k = 0 
value = 1 
it.row = 1 
it.col = 0 
index = 1 
Overall k = 1 

Inner k = 1 
value = 1 
it.row = 1 
it.col = 1 
index = 1 
Inner k = 1 
value = 1 
it.row = 4 
it.col = 1 
index = 4 
Overall k = 2 

Inner k = 2 
value = 1 
it.row = 3 
it.col = 2 
index = 3 
Overall k = 3 

Overall k = 4 

---

有關更多詳細信息Eigen和Rcpp中的稀疏矩陣，您可能希望閱讀Soren Hojsgaard和Dou的Rcpp Gallery: Using iterators for sparse vectors and matrices克貝茨。