《线性代数及其应用》笔记17.矩阵分解-QR分解.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# QR分解\n",
    "分解要求：矩阵的列线性无关。\n",
    "## 1.介绍\n",
    "$\\quad$如果$m\\times n$矩阵$A$的列线性无关，那么$A$可以分解为$A=QR$，其中$Q$是一个$m\\times n$矩阵，其列形成$ColA$的一个[标准正交基](《线性代数及其应用》笔记13.正交与正交基.ipynb)，$R$是一个$n\\times n$上三角可逆矩阵且在对角线上的元素为正数。$QR$分解常用于求解[特征值](《线性代数及其应用》笔记10.特征值和特征向量.ipynb)以及[最小二乘求解](《线性代数及其应用》笔记19.最小二乘法.ipynb)。  \n",
    "\n",
    "<img src=\"./_image/17_1.png\" width=\"300\" height=\"250\" />  \n",
    "\n",
    "$\\quad$对$m\\times n$矩阵$A$做$QR$分解步骤如下：\n",
    "* 对$A$的列向量应用[格拉姆-施密特方法](《线性代数及其应用》笔记15.正交矩阵与格拉姆-施密特方法.ipynb)得到一组正交基，再单位化后得到单位正交基，单位正交基构成矩阵Q的列；\n",
    "* 可知$Q^{T}Q=I$（正交矩阵的定理），所以有$Q^{T}A=Q^{T}QR=R$，因此计算$Q^{T}A$即可得到$R$。  \n",
    "\n",
    "举一个例子：\n",
    "\n",
    "<img src=\"./_image/17_2.png\" width=\"300\" height=\"250\" />  \n",
    "\n",
    "## 2.QR分解求特征值\n",
    "$\\quad$$QR$分解被广泛地用来估计一般矩阵$A$的特征值。算法流程如下：  \n",
    "1. $A$被分解为$Q_{1},R_{1}$，计算$R_{1}Q_{1}$得到$A_{1}$；\n",
    "2. $A_{1}$被分解为$Q_{2},R_{2}$，计算$R_{2}Q_{2}$得到$A_{2}$；  \n",
    "3. 以此类推，计算到第$k$次时，$A_{k}$几乎是上三角的，且对角线上单元素近似于$A$的特征值。  \n",
    "\n",
    "可以简单地证明，$A$是相似于$A_k$的：\n",
    "\n",
    "<img src=\"./_image/17_3.png\" width=\"300\" height=\"250\" />  \n",
    "\n",
    "$\\quad$我们用Python来实现该算法，取k=20或100来对比，$A=\n",
    "\\begin{bmatrix}\n",
    "   2 & 3\\\\\n",
    "   3 & -6\n",
    "  \\end{bmatrix}$，代码如下所示："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "without using qr, the eigenvalues is \n",
      "[ 3.+0.j -7.+0.j]\n",
      "without using qr, the eigenvectors is \n",
      "[[ 0.9486833  -0.31622777]\n",
      " [ 0.31622777  0.9486833 ]]\n",
      "\n",
      "after 20 qr, the mat is \n",
      "[[-7.00000000e+00  1.31094825e-06]\n",
      " [ 1.31094825e-06  3.00000000e+00]]\n",
      "after using qr, the eigenvalues is \n",
      "[-7.+0.j  3.+0.j]\n",
      "after using qr, the eigenvectors is \n",
      "[[-1.00000000e+00 -1.31094825e-07]\n",
      " [ 1.31094825e-07 -1.00000000e+00]]\n",
      "\n",
      "after 100 qr, the mat is \n",
      "[[-7.00000000e+00 -1.41631858e-15]\n",
      " [ 2.08884887e-43  3.00000000e+00]]\n",
      "after using qr, the eigenvalues is \n",
      "[-7.+0.j  3.+0.j]\n",
      "after using qr, the eigenvectors is \n",
      "[[ 1.00000000e+00 -1.41631858e-16]\n",
      " [ 0.00000000e+00  1.00000000e+00]]\n",
      "\n"
     ]
    }
   ],
   "source": [
    "import numpy as np  \n",
    "from scipy.linalg import eig  \n",
    "from scipy.linalg import qr  \n",
    "def qr_eig(A, times):  \n",
    "    for i in range(times):  \n",
    "        q, r = qr(A)  \n",
    "        A = np.dot(r, q)  \n",
    "    return A  \n",
    "  \n",
    "if __name__ == \"__main__\":  \n",
    "  \n",
    "    mat = np.asarray([[2, 3], [3, -6]])  \n",
    "    eigenvalues, eigenvectors = eig(mat)  \n",
    "    print('without using qr, the eigenvalues is \\n%s'%eigenvalues)  \n",
    "    print('without using qr, the eigenvectors is \\n%s\\n'%eigenvectors)  \n",
    "  \n",
    "    times = 20  \n",
    "    mat = qr_eig(mat, times)  \n",
    "    print('after %d qr, the mat is \\n%s'%(times, mat))  \n",
    "    eigenvalues, eigenvectors = eig(mat)  \n",
    "    print('after using qr, the eigenvalues is \\n%s'%eigenvalues)  \n",
    "    print('after using qr, the eigenvectors is \\n%s\\n'%eigenvectors)  \n",
    "  \n",
    "    times = 100  \n",
    "    mat = qr_eig(mat, times)  \n",
    "    print('after %d qr, the mat is \\n%s' % (times, mat))  \n",
    "    eigenvalues, eigenvectors = eig(mat)  \n",
    "    print('after using qr, the eigenvalues is \\n%s' % eigenvalues)  \n",
    "    print('after using qr, the eigenvectors is \\n%s\\n' % eigenvectors) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "k=100与k=20相比，通过$QR$分解计算出来的矩阵除对角线外元素更趋于0；通过$QR$分解计算得到的矩阵对角线上元素均趋于特征值。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}