我正在开发一个 RcppArma 包,我正在集中和标准化提升算法的设计矩阵,这是剥离的代码:
// [[Rcpp::depends(RcppArmadillo)]]
#include <RcppArmadillo.h>
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
List centering(arma::mat & X) {
int p = X.n_cols;
rowvec meanx(p);
rowvec sigmax(p);
meanx=mean(X,0);
sigmax=stddev(X,0);
for(int j=0;j<p;j++)
{
X.col(j)=X.col(j)-meanx(j);
X.col(j)=X.col(j)/sigmax(j);
}
return List::create(Named("sigma") = sigmax, Named("X") = X);
}
居中工作正常,但是居中后向量"sigmax"的所有值都等于"1",因此向量以某种方式将自身更新为居中矩阵 X 的新标准偏差,而不会重新分配。我需要原始值来反向转换系数。为什么要这样做?我怎样才能避免这种情况?
我使用以下代码在 R 中测试了该函数:
set.seed(42)
X <- replicate(10, rnorm(100, 5, 3))
res <- centering(X)
res <- centering(X)
当我第二次调用它时,问题出现了。第一次有效。
简单地说:不要在函数定义中的参数X
旁边使用引用(&
)。这激活了 RcppArmadillo 对arma::mat
重用 R 对象内存 (c.f. include/RcppArmadilloWrap.h
)
因此,要解决此问题,请从:
List centering_reuse_memory(arma::mat & X) {
# ^ reference/reuse
# Routine given in OP
}
自:
List centering_new_memory(arma::mat X) {
# ^ Direct copy
# Routine given in OP
}
了解共享...
让我们来看看对象是如何变化的。
# Create the original object
set.seed(42)
X <- replicate(3, rnorm(5, 5, 3))
# Create a duplicate object not sharing memory with X
set.seed(42)
X_clone <- replicate(3, rnorm(5, 5, 3))
# View object
X
# [,1] [,2] [,3]
# [1,] 9.112875 4.681626 8.9146090
# [2,] 3.305905 9.534566 11.8599362
# [3,] 6.089385 4.716023 0.8334179
# [4,] 6.898588 11.055271 4.1636337
# [5,] 6.212805 4.811858 4.6000360
# Check equality
all.equal(X, X_clone)
# [1] TRUE
现在,使用参数传递运行函数 arma::mat & X
res <- centering_reuse_memory(X)
# Verify results are the same.
all.equal(X, X_clone)
# [1] "Mean relative difference: 8.387859"
# Check X manually to see what changed...
X
# [,1] [,2] [,3]
# [1,] 1.34167459 -0.7368308 0.6566715
# [2,] -1.45185917 0.8327104 1.3376293
# [3,] -0.11282266 -0.7257062 -1.2116948
# [4,] 0.27645691 1.3245379 -0.4417510
# [5,] -0.05344967 -0.6947113 -0.3408550
为什么会有区别?好吧,通过使用引用,C++函数中的修改传播回驻留在 R 中的 X
变量,该变量与存储在 res$X
处的对象匹配。
# Verify R's X matches the saved C++ routine X modification
all.equal(X, res$X)
# [1] TRUE