从QR对象中获取R矩阵

Get the R matrix from QR object

本文关键字：矩阵 获取 QR 对象      更新时间：2023-10-16

我有以下R代码

P <- matrix(...)
qrP <- qr(t(P))
qR <- qr.R(qrP)

其中CCD_ 2被给出作为输入。

我正在尝试使用Eigen:在C++中编写相同的代码

auto qrP = P.transpose().fullPivHouseholderQr();
auto qr = qrP.matrixQR().template triangularView<Upper>();

但问题是矩阵不同（R与C++）。我是不是用错误的方法计算qr矩阵？

这是我打印qR对角线时得到的结果：

diag(qR)
# -1.0000000 -2.1718017 -0.4788378  0.0000000  0.0000000
cout << qr.diagonal();
// -370.247 1.37452 1 -1.5099e-14 -1.16018e-14

在Eigen版本中，您使用的是具有完全枢轴的QR因子分解，而R调用Lapack的DGEQP3例程，该例程对应于具有列枢轴的QR。在Eigen中，它可以通过colPivHouseholderQr方法或ColPivHouseholderQR类获得。

此外，您在一定程度上滥用了auto关键字。请看这张纸条。因此，一个更安全、更接近R的实现是：

ColPivHouseholderQR<MatrixXd> qrT(T.transpose());
MatrixXd q = qrT.matrixQR().triangularView<Upper>();
std::cout << q.diagonal().transpose() << std::endl;
相关文章：