RcppArmadillo中的QR分解 [英] QR decomposition in RcppArmadillo

问题描述

真的很困惑，为什么使用RcppArmadillo的QR输出不同于R的QR输出；为什么?犰狳的文件也没有给出明确的答案.本质上，当我给R矩阵n为n * q(例如1000 X 20)时，我得到的Q是1000 X 20和R 20 X1000.这就是我所需要的.但是当我在Armadillo中使用QR解算器时，它使我回想起Q 1000 X 1000和R 1000 X 20.我可以代替调用R的qr函数吗?我需要Q的维数为n x q，而不是q x q.下面的代码是我正在使用的代码(它是更大功能的一部分).

Really confused why the QR output using RcppArmadillo is different than QR output from R; Armadillo documentation doesnt give a clear answer either. Essentially when I give R a matrix Y that is n * q (say 1000 X 20 ) , I get back Q which is 1000 X 20 and R 20 X 1000. This is what I need. But when I use the QR solver in Armadillo, it throws me back Q 1000 X 1000 and R 1000 X 20. Can I call R's qr function instead? I need Q to have dimension n x q, not q x q. Code below is what I am using(its a part of a bigger function).

如果有人可以建议如何在RcppEigen中进行操作，那也将有所帮助.

If someone can suggest how to do it in RcppEigen, that'd be helpful too.

library(inline)
library(RcppArmadillo)

src <- '
    Rcpp::NumericMatrix Xr(Xs);
    int q = Rcpp::as<int>(ys);

    int n = Xr.nrow(), k = Xr.ncol();
    arma::mat X(Xr.begin(), n, k, false);

    arma::mat G, Y, B;

    G = arma::randn(n,q);

    Y = X*G;

    arma::mat Q, R;
    arma::qr(Q,R,Y);

    return Rcpp::List::create(Rcpp::Named("Q")=Q,Rcpp::Named("R")=R,Rcpp::Named("Y")=Y);'


rsvd <- cxxfunction(signature(Xs="numeric", ys="integer"), body=src, plugin="RcppArmadillo")

推荐答案

(注意:该答案解释了为什么R和RcppArmadillo返回具有不同维数的矩阵，而不是如何使RcppArmadillo返回R所做的1000 * 20矩阵. ，也许看看在qr.Q()函数定义底部附近使用的策略.)

(NOTE: This answer explains why R and RcppArmadillo return matrices with different dimensions, but not how to make RcppArmadillo return the 1000*20 matrix that R does. For that, perhaps take a look at the strategy used near the bottom of the qr.Q() function definition.)

R的qr()函数不会直接返回Q.为此，您需要使用qr.Q()，如下所示:

R's qr() function does not return Q directly. For that, you need to use qr.Q(), like this:

m <- matrix(rnorm(10), ncol=2)
qr.Q(qr(m))
#             [,1]        [,2]
# [1,] -0.40909444  0.05243591
# [2,]  0.08334031 -0.07158896
# [3,]  0.38411959 -0.83459079
# [4,] -0.69953918 -0.53945738
# [5,] -0.43450340  0.06759767

请注意，qr.Q()返回的矩阵尺寸与m相同，而不是完整的5 * 5 Q矩阵.您可以使用complete=参数来控制此行为，并获得完整尺寸的Q矩阵:

Notice that qr.Q() returns a matrix of the the same dimension as m, rather than a full 5*5 Q matrix. You can use the complete= argument to control this behavior, and get a Q matrix of full dimension:

qr.Q(qr(m), complete=TRUE)
#             [,1]        [,2]       [,3]       [,4]        [,5]
# [1,] -0.40909444  0.05243591  0.3603937 -0.7158951 -0.43301590
# [2,]  0.08334031 -0.07158896 -0.8416121 -0.5231477  0.07703927
# [3,]  0.38411959 -0.83459079  0.2720003 -0.2389826  0.15752300
# [4,] -0.69953918 -0.53945738 -0.2552198  0.3453161 -0.18775072
# [5,] -0.43450340  0.06759767  0.1506125 -0.1935326  0.86400136

在您的情况下，听起来RcppArmadillo正在返回完整的1000x1000 Q矩阵(如qr.Q(q(m, complete=FALSE)))，而不是仅返回其前20列(如qr.Q(q(m))).

In your case, it sounds like RcppArmadillo is returning the full 1000x1000 Q matrix (like qr.Q(q(m, complete=FALSE)) would), rather than just its 1st 20 columns (like qr.Q(q(m)) would).

这篇关于RcppArmadillo中的QR分解的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！