我正在尝试对大于 788000×788000 的巨大矩阵进行特征值分解以进行医学图像分析。矩阵不是稀疏的，矩阵中的每个元素都有一个实值。并且，例如，我想获得前 20 个最大的特征值和它们对应的 20 个特征向量。
虽然我的电脑配置非常出色，但是电脑无法对巨大的矩阵和内存溢出进行特征值分解。我用 Python 语言和其他相关包（如 NumPy、OpenCV、Matplotlib 等）编写计算机代码。是否有任何其他 Python 库或相关包可以进行特征值分解并解决计算问题？或者，有没有其他方法可以用 Python 解决这个问题？
我现在处境艰难，希望有人能帮助我。太感谢了。
非常抱歉，我写错了，巨大的矩阵也是对称的。
看看对面部识别做类似事情的文献——例如，搜索术语“eigenface”。
在这种情况下要说明的是，您正在寻找的信息实际上并不需要您考虑高分辨率图像。您可能有像素，在您可能想到的时间范围内，任何重要的计算都是不可行的。但是，当您查看特征向量分解（“特征脸”分解）时，您真的只关心事物的大规模特征，因此人们通常将图像缩小到或类似的低分辨率，并对这些进行计算。对于这些大小，特征向量分解、主成分分析等是可能的，因为您获得的矩阵相对较小。
几十年来，其他社区已经知道基本原则。例如，在 PDE 求解器社区中，这被称为“多级”分解，而在其他社区中，它被称为“分层分解”。这个想法是你可能有很多数据，但其中的信息可以用少得多的位数很好地近似。例如，如果您对全身的 X 射线图像感兴趣，如果您的目标是在个体之间进行比较（您可能想要计算“特征骨架”），则通常不需要极高的分辨率。如果你对单个人肺的特定分支模式感兴趣，那就不一样了，，而不是您在对 PCA 或特征向量分解感兴趣时考虑的数百个想象的统计属性。
如上所述，Thijs Steel 随机 svd 是一种解决方案，但数字 78800000 超出了我们计算机的计算能力。因此，您可以通过以下方式继续 rsvd 算法：
import numpy as np
n = 788
mu = 0
sigma = 1
A = np.random.normal(mu, sigma, (n,n))
Omega = np.random.normal(mu, sigma, (n,n))
def rsvd(A, Omega):
Y = A @ Omega
Q, _ = np.linalg.qr(Y)
B = Q.T @ A
u_tilde, s, v = np.linalg.svd(B, full_matrices = 0)
u = Q @ u_tilde
return u, s, v
u,s,v = rsvd(A,Omega)
让我们系统地进行：
如果您的矩阵确实是 ~10^6 x 10^6 元素，并且它来自医学成像，那么了解您拥有多少有效数字（如果您更喜欢您的准确度/噪声级别）非常重要。有效数字太少，那么您的问题在很大程度上取决于主要成分的强度（这就是全部，不是吗？）与其他特征值相比。如果它们没有很好地分开，那么您就无法在不丢失信息的情况下将最大的 20 个与其他 20 个分开。我还要在这里提醒您，特征向量中的误差取决于特征值的分离（即它们在绝对值中的差异）。如果您的可用数字很少，并且前 20 个特征值与测试之间的距离很小，那么错误将是巨大的，并且任何解决方案都可能会产生误导。还，假设您的数据精确到小数点后 4 位，那么，考虑到您的问题的规模，您将拥有非常高的特征值“密度”，这将使选择前 20 名充其量只是彩票。你需要比 20 更多的特征值。
这种矩阵的标准方法需要~10^13 字节（~10 TB）的双精度（你需要：单精度是没用的），这当然很多。这是标准库无法处理的。所以这条路很近。此外，非常重要的是，标准方法在计算特征值之前依赖于一些正交归约。这些方法需要 O(N^3) 次操作，并且因为矩阵是对称的（哎呀！），这些缩减不是非常可并行化的（每行/列一个全局同步）。因此，您将需要 ~10^18 次浮点运算（所有其他阶段都需要相同数量级的运算），大约为 1,000 Petaflops。最快的计算机系统总共大约 90 PFlops，在那个系统上，由于正交化的并行性差，如果可行，您需要等待几天/几周/几个月才能找到解决方案。在能够达到 1 TFlop 的系统上，您需要等待 1,000,000 小时，即查理曼大帝（约 800 广告）和现在之间的时间。因此，无法使用标准方法 tghrough 库。这是不可能的，走这条路只会导致挫败感。
那么，我们还剩下什么：可能的核外方法？我只能想到一个替代方案：一种基于矩阵向量乘积的核外方法，并且每次需要计算乘积时都需要重新加载部分矩阵。这需要大量的专业知识来最大化吞吐量，但请继续阅读。我已经使用了很好的效果（但也可以使用其他方法）结合 Lanczos 方法来生成合适的向量子空间（请参阅更多内容），然后进行一些子空间幂迭代。子空间：许多向量（> 20），它们近似于最大特征向量的向量空间（最大：最大特征值）。我不会在这里详细介绍。所以首先通过一系列 Krylov 方法步骤（Lanczos？CG？）来获得一些起点，比如 m > 20 个向量。然后使用子空间迭代从子空间给定的“池”中得到最大的特征向量。它确实有效，只要特征值之间有足够的差距（见上文），并且还允许退化（即多个特征值）。成本：每个矩阵向量乘积约 10^13 字节传输（当然是核外矩阵）。它需要〜m10 ^ 7 内存（如果 m 是 100，说一些千兆字节的顺序，非常微不足道）。如果需要 L (L>m) 次 Lanczos (CG) 迭代，则每次迭代都需要至少一个矩阵向量乘积，则成本将 > ~L 10^12 次浮点数，大约 100 次浮点数；如果子空间迭代需要 S 次迭代，则成本将是 ~ S m10^12，如果 S 再次为 100，则总共大约 10 Petaflops。因此，如果您有一个提供 1TFlops 的系统，您可以在 100 (Lacnzos) + 10000 (迭代) 秒内执行计算，或者大约几个小时。但是，它需要 ~(L+S) 矩阵加载。所以它需要，如果我们可以以 1Gbyte/sec (L+S)*10^13/10^9 ~ 10^6 秒的速度加载矩阵。但是，这里有一个有用的点：它可以并行化：如果你将矩阵分成“块”，每个处理器可以加载矩阵的一块，并且一些向量可以被复制（但它们包含更少的元素），尽管应该添加向量并跨处理器同步 - 但它们比矩阵元素少得多！所以假设，我们将批次分成 100 个处理器，忽略同步时间，整个过程可以在 (10^4+10^6)/100 ~ 10^4 秒内完成，只需几个小时。所以，这样做是可行的，但是......它需要相当多的专业知识来编写必要的定制代码。
总结一下：
第 1 点是绝对必要的。如果数字不支持计算，请不要这样做。任何结果都是无用的。
忘记标准库。无法使用标准方法。
它可以通过核外迭代方法来完成，但这将是一条非常陡峭的路径。
最重要的是上面的第 1 点。
顺便说一句，恒定图像当然只有一个非零特征值。具有随机像素的图像将具有 >> 20 个特征值。
希望这个对你有帮助。
一切顺利
可以办到。早在 1990 年代，某些小组（想想雷达横截面计算）就在使用 1,000,000 x 1,000,000 密集系统。与其他回复一致，这是在核心之外反复进行的。