Python之QR分解的简单学习

@lik-F

于 2024-10-12 22:08:43 发布

阅读量1.3k

点赞数 24

文章标签： python 学习开发语言

本文链接：https://blog.csdn.net/weixin_75064522/article/details/142886256

版权

1. QR分解的基本原理

QR分解是线性代数中一种重要的矩阵分解技术，它将一个矩阵分解为一个正交矩阵和一个上三角矩阵的乘积。其基本原理如下：

给定一个m×n矩阵A，QR分解将A分解为：

$A=QR$
其中：

Q是一个m×m正交矩阵，即Q的转置等于其逆矩阵。

R是一个m×n上三角矩阵，即R的下方元素均为0。

2. QR分解的数学性质

2.1 QR分解的几何意义

QR分解将一个矩阵分解为一个正交矩阵和一个上三角矩阵。正交矩阵的列向量是单位正交向量，即它们相互垂直且具有单位长度。上三角矩阵是一个对角线以上元素为零的矩阵。

从几何上讲，QR分解将一个矩阵分解为一个旋转和平移的组合。正交矩阵代表旋转，它将矩阵的列向量旋转到标准正交基中。上三角矩阵代表平移，它将旋转后的列向量平移到目标位置。

2.2 QR分解的正交性

QR分解的一个重要性质是正交性。正交矩阵的转置等于其逆矩阵，即正交矩阵的列向量是正交的。

import numpy as np
A = np.array([[1, 2], [3, 4]])
Q, R = np.linalg.qr(A)
print("正交矩阵 Q：")
print(Q)
print("上三角矩阵 R：")
print(R)
print("Q 的转置：")
print(Q.T)
print("Q 的逆矩阵：")
print(np.linalg.inv(Q))

2.3 QR分解的稳定性

QR分解的另一个重要性质是稳定性。对于一个给定的矩阵，QR分解的结果对于矩阵的微小扰动是稳定的。即QR分解可以鲁棒地处理输入数据中的噪声和误差。

import numpy as np
# 原始矩阵 A
A = np.array([[1, 2], [3, 4]])
# 加入噪声
A_noise = A + 0.1 * np.random.randn(*A.shape)
# 对原始矩阵和加入噪声的矩阵进行 QR 分解
Q1, R1 = np.linalg.qr(A)
Q2, R2 = np.linalg.qr(A_noise)
# 计算两个 QR 分解结果之间的差异
diff = np.linalg.norm(Q1 - Q2) + np.linalg.norm(R1 - R2)
print("QR 分解结果差异：", diff)

3.QR分解的实际应用

3.1 线性方程组求解

QR分解在求解线性方程组方面有着广泛的应用。给定一个线性方程组：

$Ax=b$

其中 A 是一个 m×n 矩阵，x 是一个 n×1 列向量，b 是一个 m×1 列向量。

QR分解可以将矩阵 A 分解为：

$A = QR$

其中 Q 是一个 m×n 正交矩阵，R 是一个 n×n 上三角矩阵。

利用QR分解，我们可以将求解线性方程组转化为求解上三角方程组：

$Rx = Q^Tb$

上三角方程组可以通过回代法轻松求解。

import numpy as np
# 随机生成一个矩阵 A
A = np.random.rand(5, 3)
# QR分解
Q, R = np.linalg.qr(A)
# 求解线性方程组
b = np.random.rand(5, 1)
x = np.linalg.solve(R, Q.T @ b)
print(x)

分析：

np.linalg.qr(A) 函数执行 QR 分解，返回正交矩阵 Q 和上三角矩阵 R。

Q.T @ b 计算 Q 的转置与 b 的乘积，得到上三角方程组的右端项。

np.linalg.solve(R, Q.T @ b) 函数求解上三角方程组，得到解向量 x。

3.2 最小二乘问题求解

最小二乘问题是指求解一个函数 f(x) 的最小值，其中 f(x) 是一个非线性函数。QR分解可以将最小二乘问题转化为一个线性方程组求解问题。

给定一个最小二乘问题：

$minf(x)= \frac{1}{2}\left \| Ax-b \right \|^{2}$

其中 A 是一个 m×n 矩阵，x 是一个 n×1 列向量，b 是一个 m×1 列向量。

QR分解可以将矩阵 A 分解为：

$A=QR$

其中 Q 是一个 m×n 正交矩阵，R 是一个 n×n 上三角矩阵。

利用QR分解，我们可以将最小二乘问题转化为求解上三角方程组：

$R^Tx = Q^Tb$

import numpy as np
# 随机生成一个矩阵 A
A = np.random.rand(5, 3)
# QR分解
Q, R = np.linalg.qr(A)
# 求解最小二乘问题
b = np.random.rand(5, 1)
x = np.linalg.solve(R.T, Q.T @ b)
print(x)

分析：

np.linalg.qr(A) 函数执行 QR 分解，返回正交矩阵 Q 和上三角矩阵 R。

Q.T @ b 计算 Q 的转置与 b 的乘积，得到上三角方程组的右端项。

np.linalg.solve(R.T, Q.T @ b) 函数求解上三角方程组，得到最小二乘问题的解向量 x。

3.3 QR分解推导列向量线性无关的矩阵在值域上的投影

QR分解是一种将矩阵分解为一个正交矩阵和一个上三角矩阵的算法。对于列向量线性无关的矩阵，我们可以通过QR分解来找到这个矩阵在值域上的投影。

首先，我们需要了解QR分解的基本步骤。对于一个m×n的矩阵A，其列向量线性无关，我们可以通过Gram-Schmidt正交化过程来构造QR分解。这个过程包括以下几个步骤：

1.将第一个向量 $a_{1}$ 单位化，得到 $q_{1}= \frac{a_{1}}{\left \| a_{1} \right \|}$ ，其中 $\left \| a_{1} \right \|$ 是 $a_{1}$ 的范数。

2.对于每个后续的向量 $a_{i}$ ，计算它在之前所有 $q_{j}$ 上的投影，然后从 $a_{i}$ 中减去这些投影，得到正交向量。

3.将每个正交向量单位化，得到 $q_{i}$

通过这个过程，我们可以得到矩阵A的QR分解，其中Q是一个正交矩阵，其列向量是单位正交的，而R是一个上三角矩阵。

对于列向量线性无关的矩阵A，其值域就是由这些列向量张成的空间。矩阵A在值域上的投影可以表示为 $A(A^{T}A)^{-1}A^{T}$ ，这是一个投影矩阵，它将任何向量投影到A的值域上。

在QR分解中，由于Q是正交矩阵，它保持了向量的长度和方向，而R的上三角性质意味着它包含了A的列向量在值域上的“缩放”信息。因此，我们可以利用QR分解来简化计算，特别是当我们需要多次进行这种投影操作时。

import numpy as np

# 定义或生成一个矩阵A，其列向量线性无关
A = np.array([[1, 2, 3],
              [4, 5, 6],
              [7, 8, 10]])

# 定义一个向量x
x = np.array([1, 2, 3])

# 使用NumPy的qr函数进行QR分解
Q, R = np.linalg.qr(A)

# 计算向量x在矩阵A值域上的投影
# 首先计算Q的转置与x的乘积
Qt_x = np.dot(Q.T, x)

# 然后计算(R^T * R)^(-1) * R^T * Q^T * x
# 由于R是上三角矩阵，我们可以直接使用numpy的inv函数来计算逆
R_inv = np.linalg.inv(R)
projection = np.dot(A, np.dot(R_inv, Qt_x))

print("QR分解的结果:")
print("Q:")
print(Q)
print("R:")
print(R)
print("向量x在矩阵A值域上的投影:")
print(projection)

分析：代码首先定义了一个矩阵A和一个向量x，并numpy.linalg.qr函数进行QR分解，最后计算向量x在矩阵A值域上的投影。

3.4 QR分解在图像处理和计算机视觉中的应用

1）图像去噪：QR分解可以用于图像去噪。通过将图像矩阵分解为 QR 矩阵，可以将噪声成分分离出来，从而实现图像去噪。

2）图像压缩：QR分解可以用于图像压缩。通过将图像矩阵分解为 QR 矩阵，可以将图像信息压缩到 R 矩阵中，从而实现图像压缩。

3）计算机视觉：QR分解在计算机视觉中用于特征提取和匹配。通过对图像进行 QR 分解，可以提取图像中的特征，并通过匹配这些特征来进行图像识别和匹配。

import numpy as np
from PIL import Image
import matplotlib.pyplot as plt

# 读取图像
image = Image.open('image.jpg')
# 转换为灰度图像
image = image.convert('L')
# 将图像转换为矩阵
image_matrix = np.array(image)

# QR分解
Q, R = np.linalg.qr(image_matrix)

# 读取第二张图像
image1 = Image.open('image1.jpg')
# 确保两张图像大小相同
image1 = image1.resize(image.size)
image1 = image1.convert('L')
image1_matrix = np.array(image1)

# QR分解
Q1, R1 = np.linalg.qr(image1_matrix)

# 匹配
# 由于R和R1的形状可能不同，我们不能直接比较它们
# 这里我们使用一种简单的方法来比较图像的特征
# 例如，我们可以比较它们的行数和列数是否相同
if R.shape == R1.shape:
    print('图像匹配')
else:
    print('图像不匹配')

# 显示原始图像和第二张图像
plt.figure(figsize=(10, 5))
plt.subplot(1, 2, 1)
plt.imshow(image, cmap='gray')
plt.title('Original Image')
plt.subplot(1, 2, 2)
plt.imshow(image1, cmap='gray')
plt.title('Second Image')
plt.show()