5.5 QR分解一：施密特正交化

醒过来摸鱼

已于 2023-05-19 16:25:06 修改

阅读量1.7k

点赞数 3

分类专栏：线性代数【更新中】文章标签：线性代数矩阵算法

于 2023-01-31 20:27:25 首次发布

本文链接：https://blog.csdn.net/m0_66201040/article/details/128822415

版权

线性代数【更新中】专栏收录该内容

71 篇文章 38 订阅

订阅专栏

文章目录

计算方法
python实现
数据测试

计算方法

矩阵的正交三角分解QR decomposition，简称QR分解，是找到两个矩阵Q和R使得 $A = QR$ ,Q是正交矩阵（就是 $Q^TQ=E$ 的矩阵，也就是列向量之间互相正交的矩阵），R是对角线上没有0的上三角矩阵。
QR分解主要有三种方法：

施密特正交化
Householder变换
Givens旋转

本文只讲施密特正交化Gram-Schmidt orthogonalization，Gram-Schmit正交分解步骤很简单，是先根据原矩阵获取一组正交基，再求原矩阵在这组正交基下的坐标。以 $3\times3$ 矩阵 $w_1,w_2,w_3)$ 为例子。
$\\ u_1=w_1\\ u_2=w_2-\frac{(u_1, w_2)}{(u_1,u_1)} u_1\\ u_3=w_3-\frac{(u_1, w_3)}{(u_1,u_1)} u_1-\frac{(u_2,w_3)}{(u_2,u_2)} u_2\\ u_i=w_i-\sum_{j=1}^{j<i} \frac{(u_j,w_i)}{(u_j,u_j)} u_j\\ (w_1,w_2 ,w_3)= (u_1, u_2,u_3)\begin{pmatrix} 1& \frac{(u_1, w_2)}{(u_1,u_1)} & \frac{(u_1, w_3)}{(u_1,u_1)} \\ 0 & 1 & \frac{u_2\cdot w_3}{(u_2,u_2)}\\ 0 & 0 & 1 \end{pmatrix}\\ Q= (u_1,u_2,u_3)\\ R=\begin{pmatrix} 1& \frac{(u_1, w_2)}{(u_1,u_1)} & \frac{(u_1, w_3)}{(u_1,u_1)} \\ 0 & 1 & \frac{(u_2,w_3)}{(u_2,u_2)}\\ 0 & 0 & 1 \end{pmatrix}$
公式里的 $u_1,u_2,\cdots,u_i$ 是一组正交基。Gram-Schmit正交分解获得的两个矩阵，右边的是上三角矩阵没错，但一般在上面的步骤做完之后，还需要对得到的Q和R两个矩阵进行单位化。以 $3\times3$ 矩阵为例子，单位化是做如下处理：
$Q=\begin{pmatrix} \frac{u_1}{\sqrt{(u_1,u_1)}} &\frac{ u_2}{\sqrt{(u_2,u_2)}} & \frac{u_3}{\sqrt{u_3\cdot u_3}} \end{pmatrix}\\ R=\begin{pmatrix} \sqrt{(u_1,u_1)} & \frac{(u_1, w_2)}{\sqrt{(u_1,u_1)}}& \frac{(u_1, w_3)}{\sqrt{(u_1,u_1)}} \\ 0 & \sqrt{(u_2,u_2)} & \frac{(u_2,w_3)}{\sqrt{(u_2,u_2)}}\\ 0 & 0 & \sqrt{u_3\cdot u_3} \end{pmatrix}$
以下列矩阵为例子：
$\begin{pmatrix} 2& 1& -1\\ 1& 2& -1\\ -1& -1& 2 \end{pmatrix}\\ Q= \begin{pmatrix}2 & -0.667 & 0.364\\ 1 & 1.167 & 0.364\\ -1 & -0.167 & 1.091\\ \end{pmatrix}\\ R= \begin{pmatrix}1 & 0.833 & -0.833\\ 0 & 1 & -0.455\\ 0 & 0 & 1\\ \end{pmatrix}$
单位化之后，就是下面的结果：
$\begin{pmatrix}0.816 & -0.492 & 0.302\\ 0.408 & 0.862 & 0.302\\ -0.408 & -0.123 & 0.905\\ \end{pmatrix}\\ R= \begin{pmatrix}2.449 & 2.041 & -2.041\\ 0.0 & 1.354 & -0.615\\ 0.0 & 0.0 & 1.206\\ \end{pmatrix}$

python实现

实现代码就是严格按照公式：


    def gram_schmidt(self, matrix = None):
        # q矩阵先复制原矩阵
        q = [[i for i in vector] for vector in self.__vectors]
        # 从第二个向量开始做减法
        # q的长度数组
        q_square = [0 for _ in q]
        q_square[0] = inner_product(q[0], q[0], matrix)
        q_length = [0 for _ in q]
        q_length[0] = math.sqrt(q_square[0])
        columns = len(self.__vectors)
        # r 初始化为全部为0
        r = [[0 for _ in vector] for vector in self.__vectors]
        # r对角线为1
        for i, vector in enumerate(r):
            vector[i] = 1
        # i 是列
        for i in range(1, columns):
            q[i] = self.__vectors[i]
            # sum其实是一个向量
            for j in range(0, i):
                frac = inner_product(q[j], self.__vectors[i], matrix) / q_square[j]
                # 把系数存储在r中
                r[i][j] = frac
                q[i] = sub(q[i], mul_num(q[j], frac))
            q_square[i] = inner_product(q[i], q[i], matrix)
            q_length[i] = math.sqrt(q_square[i])
        # Q单位化
        # R乘回来
        for i, vector in enumerate(q):
            for j, e in enumerate(vector):
                vector[j] = e / q_length[i]
                # vector[j] = round(vector[j] * 10000) / 10000.0
        # i 代表向量，也就是列
        for i, vector in enumerate(r):
            # j才代表行
            for j, e in enumerate(vector):
                # 每一行乘以长度
                vector[j] = e * q_length[j]
                # vector[j] = round(vector[j] * 10000) / 10000.0
        return Matrix(q), Matrix(r)

数据测试

以下矩阵，在标准内积下分解为：
$\begin{pmatrix}2 & 1 & -1\\ 1 & 2 & -1\\ -1 & -1 & 2\\ \end{pmatrix}\\= \begin{pmatrix}0.816 & -0.493 & 0.303\\ 0.408 & 0.862 & 0.301\\ -0.408 & -0.123 & 0.906\\ \end{pmatrix} \begin{pmatrix}2.449 & 2.041 & -2.041\\ 0 & 1.353 & -0.614\\ 0 & 0 & 1.204\\ \end{pmatrix}$
定义内积为：
$(x,y)=x^T\begin{pmatrix}1 & 0 & 0\\ 0 & 2 & 0\\ 0 & 0 & 3\\ \end{pmatrix}y$
在这个内积下，QR分解为：
$\begin{pmatrix}2 & 1 & -1\\ 1 & 2 & -1\\ -1 & -1 & 2\\ \end{pmatrix} \\=\begin{pmatrix}0.667 & -0.577 & 0.471\\ 0.333 & 0.577 & 0.236\\ -0.333 & 0 & 0.471\\ \end{pmatrix} \begin{pmatrix}3 & 3 & -3.333\\ 0 & 1.732 & -0.577\\ 0 & 0 & 1.887\\ \end{pmatrix}$