使用 Numpy (np.linalg.svd) 进行奇异值分解

Posted

技术标签:

【中文标题】使用 Numpy (np.linalg.svd) 进行奇异值分解【英文标题】:Using Numpy (np.linalg.svd) for Singular Value Decomposition 【发布时间】:2014-09-14 19:11:31 【问题描述】:

我正在阅读 Abdi & Williams (2010) “主成分分析”,我正在尝试重做 SVD 以获得进一步 PCA 的值。

文章指出以下SVD:

X = P D Q^t

我将数据加载到 np.array X 中。

X = np.array(data)
P, D, Q = np.linalg.svd(X, full_matrices=False)
D = np.diag(D)

但我在检查时没有得到上述相等性

X_a = np.dot(np.dot(P, D), Q.T)

X_a 和 X 是相同的维度,但值不同。我是否遗漏了什么,或者 np.linalg.svd 函数的功能与论文中的公式不兼容?

【问题讨论】:

【参考方案1】:

我认为对于那些在 Python/linalg 库中使用 SVD 的人来说,仍然有一些要点。首先,https://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.svd.html 是 SVD 计算函数的一个很好的参考。

以 SVD 计算为 A= U D (V^T), 对于 U, D, V = np.linalg.svd(A),此函数已经以 V^T 形式返回 V。 D也只包含特征值,因此它必须被塑造成矩阵形式。因此重建可以用

import numpy as np
U, D, V = np.linalg.svd(A)
A_reconstructed = U @ np.diag(D) @ V

关键是,如果 A 矩阵不是方阵而是矩形矩阵,这行不通,你可以改用它

import numpy as np
U, D, V = np.linalg.svd(A)
m, n = A.shape
A_reconstructed = U[:,:n] @ np.diag(D) @ V[:m,:]

或者您可以在 SVD 函数中使用 'full_matrices=False' 选项;

import numpy as np
U, D, V = np.linalg.svd(A,full_matrices=False)
A_reconstructed = U @ np.diag(D) @ V

【讨论】:

【参考方案2】:

TL;DR:numpy 的 SVD 计算 X = PDQ,所以 Q 已经转置。

SVD 有效地将矩阵X 分解为旋转PQ 以及对角矩阵Dlinalg.svd() 的版本我已经返回了 PQ 的正向旋转。在计算 X_a 时,您不想转换 Q

import numpy as np
X = np.random.normal(size=[20,18])
P, D, Q = np.linalg.svd(X, full_matrices=False)
X_a = np.matmul(np.matmul(P, np.diag(D)), Q)
print(np.std(X), np.std(X_a), np.std(X - X_a))

我得到:1.02、1.02、1.8e-15,表明X_a 非常准确地重构了X

如果您使用的是 Python 3,@ 运算符实现矩阵乘法并使代码更易于理解:

import numpy as np
X = np.random.normal(size=[20,18])
P, D, Q = np.linalg.svd(X, full_matrices=False)
X_a = P @ diag(D) @ Q
print(np.std(X), np.std(X_a), np.std(X - X_a))
print('Is X close to X_a?', np.isclose(X, X_a).all())

【讨论】:

根据np.dot's documentation,矩阵乘法首选np.matmul 根据罗德里戈的评论更新了答案。还添加了较新的“@”符号。【参考方案3】:

来自 scipy.linalg.svd 文档字符串,其中 (M,N) 是输入矩阵的形状,K 是两者中较小的一个:

Returns
-------
U : ndarray
    Unitary matrix having left singular vectors as columns.
    Of shape ``(M,M)`` or ``(M,K)``, depending on `full_matrices`.
s : ndarray
    The singular values, sorted in non-increasing order.
    Of shape (K,), with ``K = min(M, N)``.
Vh : ndarray
    Unitary matrix having right singular vectors as rows.
    Of shape ``(N,N)`` or ``(K,N)`` depending on `full_matrices`.

如上所述,Vh 是 Abdi 和 Williams 论文中使用的 Q 的转置。所以只是

X_a = P.dot(D).dot(Q)

应该给你答案。

【讨论】:

以上是关于使用 Numpy (np.linalg.svd) 进行奇异值分解的主要内容,如果未能解决你的问题,请参考以下文章

基于最小二乘法的线性回归

python numpy访问行列元素的方法

python vs R中的矩阵乘法没有为SVD白化返回相同的结果

多处理队列已满

numpy 使用

numpy使用np.printoptions函数抑制numpy对于长numpy数组的自动省略显示numpy数组中的所有数据