文章来源于网络收集而来，版权归原创者所有，如有侵权请及时联系！

QR 算法

发布于 2025-01-01 12:38:42 字数 10618 浏览 0 评论 0 收藏 0

我们使用幂方法，找到对应维基百科链接矩阵的最大特征值的特征向量。这个特征向量给了我们每个维基百科页面的相对重要性（就像简化的 PageRank）。

接下来，让我们看一个方法，查找对称正定矩阵的所有特征值。该方法包括数值线性代数中的 2 种基本算法，并且是许多更复杂方法的基础。

Google Matrix 的第二个特征值：具有“implications for the convergence rate of the standard PageRank algorithm as the web scales, for the stability of PageRank to perturbations to the link structure of the web, for the detection of Google spammers, and for the design of algorithms to speed up PageRank”。

避免混淆：QR 算法和 QR 分解

QR 算法使用称为 QR 分解的东西。两者都很重要，所以不要混淆他们。 QR 分解将矩阵 A = QR 分解为一组正交列 Q 和三角矩阵 R 。我们将在未来的课程中查看几种计算 QR 分解的方法。现在，只要知道它给了我们一个正交矩阵和一个三角矩阵。

线性代数

如果存在非奇异矩阵 X ：

$B = X^{-1}AX$

则两个矩阵 A 和 B 是相似的。

观看这个：基变换

定理：如果 X 是非奇异的，那么 A 和 $X^{-1}AX$ 具有相同的特征值。

其他资源

回顾：线性组合，跨度和基向量。

上述定理的证明（以及更多！）请参阅第 24 讲。

算法

QR 算法的最基本版本：

for k=1,2,...
    Q, R = A
    A = R @ Q

在适当的假设下，该算法收敛于 A 的 Schur 形式！

工作原理

再写一次，只有下标：

$A_0 = A$

对于 $k=1,2,\ldots$

$Q_k, R_k = A_{k-1} \\ A_k = R_kQ_k$

我们可以将其视为构建 $A_k$ ， $Q_k$ 和 $R_k$ 的序列。

$A_k = Q_k \, R_k \\ Q_k^{-1} \, A_k = R_k$

所以：

$R_k Q_k = Q_k^{-1} \, A_k \, Q_k \\ A_k = Q_k^{-1} \ldots Q_2^{-1} Q_1^{-1} A Q_1 Q_2 \dots Q_k$

Trefethen 在第 216-217 页证明了以下内容：

$A^k = Q_1 Q_2 \dots Q_k R_k R_{k-1}\dots R_1$

要点：QR 算法为连续的幂构造 $A_k$ 的正交基。记住 A 的幂与特征值分解之间的密切关系。

要了解更多信息，请阅读瑞利熵。

纯 QR

n = 6
A = np.random.rand(n,n)
AT = A @ A.T

def pure_qr(A, max_iter=1000):
    Ak = np.copy(A)
    n = A.shape[0]
    QQ = np.eye(n)
    for k in range(max_iter):
        Q, R = np.linalg.qr(Ak)
        Ak = R @ Q
        QQ = QQ @ Q
        if k % 100 == 0:
            print(Ak)
            print("\n")
    return Ak, QQ

A

'''
array([[ 0.62878,  0.23258,  0.63909,  0.90223,  0.94772,  0.80247],
       [ 0.64361,  0.52469,  0.92231,  0.32869,  0.58532,  0.75104],
       [ 0.44363,  0.00427,  0.62418,  0.47093,  0.6762 ,  0.28078],
       [ 0.14629,  0.76324,  0.23316,  0.55208,  0.21712,  0.20255],
       [ 0.56122,  0.08282,  0.12788,  0.10419,  0.40358,  0.69272],
       [ 0.41172,  0.06411,  0.92162,  0.53139,  0.27901,  0.61592]])
'''

Ak, Q = pure_qr(A)

'''
[[ 2.65646  0.21386  0.16765  0.75256  0.61251  0.93518]
 [ 0.52744  0.47579  0.17052 -0.41086 -0.21182 -0.01876]
 [ 0.29923  0.06964  0.11173  0.1879  -0.29101  0.60032]
 [ 0.2274   0.46162 -0.26654  0.08899  0.24631  0.26447]
 [-0.06093  0.02892  0.34162  0.07533  0.02393 -0.05456]
 [-0.06025  0.02694 -0.11675 -0.00927 -0.11939 -0.00767]]


[[ 2.78023  0.52642  0.0395  -0.11135  0.1569   1.15184]
 [ 0.       0.18624 -0.297   -0.07256 -0.04537  0.27907]
 [ 0.       0.69328  0.34105 -0.12198  0.11029  0.0992 ]
 [-0.      -0.0494  -0.02057  0.09461  0.59466  0.09115]
 [-0.       0.00008 -0.02659 -0.40372  0.06542  0.38612]
 [-0.       0.       0.       0.      -0.      -0.11832]]


[[ 2.78023 -0.12185 -0.51401  0.17625 -0.07467  1.15184]
 [ 0.       0.2117  -0.70351  0.09974 -0.02986  0.00172]
 [ 0.       0.28284  0.32635 -0.0847  -0.08488 -0.29104]
 [-0.      -0.00068 -0.00088 -0.01282  0.54836  0.13447]
 [-0.       0.      -0.00102 -0.45718  0.16208 -0.37726]
 [-0.       0.       0.       0.      -0.      -0.11832]]


[[ 2.78023 -0.33997  0.4049   0.17949  0.06291  1.15184]
 [ 0.       0.48719 -0.48788 -0.05831 -0.12286 -0.23486]
 [ 0.       0.49874  0.05104 -0.07191  0.03638  0.17261]
 [ 0.       0.00002  0.       0.02128  0.41958  0.3531 ]
 [ 0.      -0.       0.00002 -0.58571  0.1278  -0.18838]
 [ 0.       0.       0.      -0.       0.      -0.11832]]


[[ 2.78023  0.35761  0.38941  0.07462  0.17493  1.15184]
 [ 0.       0.0597  -0.55441 -0.0681  -0.04456  0.14084]
 [ 0.       0.43221  0.47853 -0.06068  0.12117  0.25519]
 [-0.      -0.      -0.       0.16206  0.45708  0.37724]
 [-0.       0.      -0.      -0.54821 -0.01298  0.1336 ]
 [ 0.       0.       0.       0.      -0.      -0.11832]]


[[ 2.78023  0.06853 -0.52424  0.05224 -0.18287  1.15184]
 [ 0.       0.36572 -0.6889   0.07864 -0.09263  0.105  ]
 [ 0.       0.29772  0.17251 -0.09836 -0.02347 -0.27191]
 [ 0.      -0.      -0.       0.13719  0.57888 -0.20884]
 [ 0.       0.      -0.      -0.42642  0.01189 -0.34139]
 [ 0.       0.       0.       0.      -0.      -0.11832]]


[[ 2.78023 -0.52782  0.03045 -0.14389 -0.12436  1.15184]
 [ 0.       0.25091 -0.27593  0.08994 -0.06581 -0.28672]
 [ 0.       0.7107   0.28732  0.10154  0.04751 -0.05245]
 [ 0.      -0.      -0.       0.0297  -0.59054 -0.01234]
 [ 0.      -0.       0.       0.41475  0.11938 -0.40001]
 [ 0.       0.       0.       0.       0.      -0.11832]]


[[ 2.78023  0.1533   0.50599  0.18983  0.01158  1.15184]
 [ 0.       0.18627 -0.69511 -0.0991  -0.00189  0.01621]
 [ 0.       0.29151  0.35196  0.05638 -0.10949  0.29102]
 [ 0.      -0.      -0.      -0.02207 -0.48261  0.25246]
 [ 0.      -0.       0.       0.52268  0.17116  0.31053]
 [ 0.       0.       0.       0.       0.      -0.11832]]


[[ 2.78023  0.29683 -0.43751 -0.13852  0.13032  1.15184]
 [ 0.       0.48375 -0.53231 -0.01164  0.13482  0.216  ]
 [ 0.       0.45431  0.05448 -0.07972 -0.01795 -0.19571]
 [ 0.      -0.       0.       0.10042 -0.40743 -0.39915]
 [ 0.      -0.       0.       0.59786  0.04867 -0.02893]
 [ 0.       0.       0.       0.       0.      -0.11832]]


[[ 2.78023 -0.39373 -0.35284  0.00838 -0.19     1.15184]
 [ 0.       0.05184 -0.51278  0.05752 -0.0564  -0.16496]
 [ 0.       0.47384  0.48639  0.09426  0.09806 -0.24031]
 [ 0.      -0.      -0.       0.17043 -0.52593  0.30622]
 [ 0.      -0.       0.       0.47936 -0.02134 -0.25766]
 [ 0.       0.       0.       0.       0.      -0.11832]]
'''

让我们和特征值比较。

np.linalg.eigvals(A)

'''
array([ 2.78023+0.j     , -0.11832+0.j     ,  0.26911+0.44246j,
        0.26911-0.44246j,  0.07454+0.49287j,  0.07454-0.49287j])
'''

检查 Q 是正交的。

np.allclose(np.eye(n), Q @ Q.T), np.allclose(np.eye(n), Q.T @ Q)

# (True, True)

非常非常慢。

实际的 QR（带移位的 QR）

不是将 $A_k$ 分解为 $Q_kR_k$ ，

1）获取 QR 分解：

$A_k - s_k I = Q_k R_k$

2）设置：

$A_{k+1} = R_k Q_k + s_k I$

选择 $s_k$ 来近似 A 的特征值。我们将使用 $s_k = A_k(m,m)$ 。

在数值线性代数（包括幂方法，逆迭代和瑞利商迭代）中的许多算法中都出现了添加移位以加速收敛的想法。

作业：向 QR 算法添加移位。

# 练习：向 QR 算法添加移位
# 练习：def practical_qr(A, iters=10):
# 练习：    return Ak, Q

Ak, Q = practical_qr(A, 10)

'''
[ 5.16264  2.53603  0.31923  0.35315  0.97569  0.43615]
[ 7.99381  0.05922  0.34478  0.29482  0.79026  0.29999]
[ 8.00491  0.04358  0.89735  0.26386  0.26182  0.31135]
[ 8.00493  0.13648  0.91881  0.14839  0.24313  0.33115]
[ 8.00493  0.43377  0.62809  0.13429  0.24592  0.33589]
[ 8.00493  0.81058  0.25128  0.13297  0.24722  0.3359 ]
[ 8.00493  0.98945  0.07221  0.13292  0.24747  0.3359 ]
[ 8.00493  1.0366   0.02497  0.13296  0.24751  0.3359 ]
[ 8.00493  1.04688  0.01465  0.13299  0.24752  0.3359 ]
[ 8.00493  1.04902  0.0125   0.13301  0.24753  0.3359 ]
'''

检查 Q 是正交的。

np.allclose(np.eye(n), Q @ Q.T), np.allclose(np.eye(n), Q.T @ Q)

# (True, True)

让我们和特征值比较。

np.linalg.eigvals(A)

'''
array([ 2.68500+0.j     ,  0.19274+0.41647j,  0.19274-0.41647j,
       -0.35988+0.43753j, -0.35988-0.43753j, -0.18346+0.j     ])
'''

问题：这比未移位的版本更好（它甚至不能保证收敛），但仍然很慢！事实上，它是 O(n^4) ，这是非常糟糕的。

在对称矩阵的情况下，它是 O(n^3)

但是，如果从海森堡矩阵（第一个子对角线下面是零）开始，它会更快： O(n^3) ，如果是对称的则是 O(n^2) 。

分享到QQ

分享到微博