Eigen学习笔记22：求解线性最小二乘系统

一个超定方程组，例如Ax = b，没有解。

在这种情况下，在令Ax - b尽可能小的意义上，搜索最接近解的向量x是有意义的。该x称为最小二乘解（如果使用欧几里得范数）。

这里讨论三种方法： the SVD decomposition, the QR decomposition and normal equations.

当然, the SVD decomposition 通常最准确但最慢, normal equations 最快但最不准确, and the QR decomposition 介于两者之间。

使用SVD分解

The solve() method in the BDCSVD class can be directly used to solve linear squares systems. 

仅计算奇异值（此类的默认值）是不够的。您还需要 the singular vectors ，但是 the singular vectors 足以计算最小二乘解：

#include <iostream>
#include <Eigen/Dense>

using namespace std;
using namespace Eigen;

int main()
{
   MatrixXf A = MatrixXf::Random(3, 2);
   cout << "Here is the matrix A:\n" << A << endl;
   VectorXf b = VectorXf::Random(3);
   cout << "Here is the right hand side b:\n" << b << endl;
   cout << "The least-squares solution is:\n"
        << A.bdcSvd(ComputeThinU | ComputeThinV).solve(b) << endl;
}

输出：

Here is the matrix A:
 0.680375   0.59688
-0.211234  0.823295
 0.566198 -0.604897
Here is the right hand side b:
-0.329554
 0.536459
-0.444451
The least-squares solution is:
-0.669626
 0.314253

使用QR分解

QR分解类中的solve()方法还计算最小二乘解。有3种QR分解类： There are three QR decomposition classes: HouseholderQR (no pivoting, so fast but unstable), ColPivHouseholderQR (column pivoting, thus a bit slower but more accurate) and FullPivHouseholderQR (full pivoting, so slowest and most stable). 

#include <iostream>
#include <Eigen/Dense>

using namespace std;
using namespace Eigen;

int main()
{
    MatrixXf A = MatrixXf::Random(3, 2);
    VectorXf b = VectorXf::Random(3);
    cout << "The solution using the QR decomposition is:\n"
         << A.colPivHouseholderQr().solve(b) << endl;
}

输出：

The solution using the QR decomposition is:
-0.669626
 0.314253

使用 normal equations

找到Ax = b的最小二乘解等效于求解法线方程

#include <iostream>
#include <Eigen/Dense>

using namespace std;
using namespace Eigen;

int main()
{
    MatrixXf A = MatrixXf::Random(3, 2);
    VectorXf b = VectorXf::Random(3);
    cout << "The solution using normal equations is:\n"
         << (A.transpose() * A).ldlt().solve(A.transpose() * b) << endl;
}

输出：

The solution using normal equations is:
-0.669626
 0.314253

如果矩阵A是病态的，那么这不是一个好方法，因为的条件数是A的条件数的平方。This means that you lose twice as many digits using normal equation than if you use the other methods.