Hub*_*ast 10 python numpy linear-algebra scipy eigenvector
这里有两个关于方阵的特征向量和特征值的假设。我相信两者都是正确的:
如果一个矩阵是对称的并且只包含实数,那么它就是一个厄米矩阵,那么所有的特征值都应该是实数,所有特征向量的所有分量也应该是实数。当您从 Hermitian 矩阵计算特征向量和特征值时,结果中不应出现复数。
从给定矩阵计算的给定特征值的特征向量应始终指向仅由矩阵和特征值确定的方向。用于计算它的算法对结果没有影响,只要正确执行算法即可。
但是,当您在 Python 中使用标准库来计算特征向量和特征值时,这两个假设都不成立。这些方法是否包含错误?
有四种不同的方法可以从 Hermitian 矩阵计算特征值和特征向量:
#1 和 #2 可用于任何方阵(包括 Hermitian 矩阵)。
#3 和 #4 仅用于 Hermitian 矩阵。据我所知,他们的目的只是为了让他们跑得更快,但结果应该是一样的(只要输入真的是 Hermitian)。
但是对于相同的输入,这四种方法会产生三种不同的结果。这是我用来测试所有四种方法的程序:
#!/usr/bin/env python3
import numpy as np
import scipy.linalg as la
A = [
[19, -1, -1, -1, -1, -1, -1, -1],
[-1, 19, -1, -1, -1, -1, -1, -1],
[-1, -1, 19, -1, -1, -1, -1, -1],
[-1, -1, -1, 19, -1, -1, -1, -1],
[-1, -1, -1, -1, 19, -1, -1, -1],
[-1, -1, -1, -1, -1, 19, -1, -1],
[-1, -1, -1, -1, -1, -1, 19, -1],
[-1, -1, -1, -1, -1, -1, -1, 19]
]
A = np.array(A, dtype=np.float64)
delta = 1e-12
A[5,7] += delta
A[7,5] += delta
if np.array_equal(A, A.T):
print('input is symmetric')
else:
print('input is NOT symmetric')
methods = {
'np.linalg.eig' : np.linalg.eig,
'la.eig' : la.eig,
'np.linalg.eigh' : np.linalg.eigh,
'la.eigh' : la.eigh
}
for name, method in methods.items():
print('============================================================')
print(name)
print()
eigenValues, eigenVectors = method(A)
for i in range(len(eigenValues)):
print('{0:6.3f}{1:+6.3f}i '.format(eigenValues[i].real, eigenValues[i].imag), end=' | ')
line = eigenVectors[i]
for item in line:
print('{0:6.3f}{1:+6.3f}i '.format(item.real, item.imag), end='')
print()
print('---------------------')
for i in range(len(eigenValues)):
if eigenValues[i].imag == 0:
print('real ', end=' | ')
else:
print('COMPLEX ', end=' | ')
line = eigenVectors[i]
for item in line:
if item.imag == 0:
print('real ', end='')
else:
print('COMPLEX ', end='')
print()
print()
Run Code Online (Sandbox Code Playgroud)
这是它产生的输出:
#!/usr/bin/env python3
import numpy as np
import scipy.linalg as la
A = [
[19, -1, -1, -1, -1, -1, -1, -1],
[-1, 19, -1, -1, -1, -1, -1, -1],
[-1, -1, 19, -1, -1, -1, -1, -1],
[-1, -1, -1, 19, -1, -1, -1, -1],
[-1, -1, -1, -1, 19, -1, -1, -1],
[-1, -1, -1, -1, -1, 19, -1, -1],
[-1, -1, -1, -1, -1, -1, 19, -1],
[-1, -1, -1, -1, -1, -1, -1, 19]
]
A = np.array(A, dtype=np.float64)
delta = 1e-12
A[5,7] += delta
A[7,5] += delta
if np.array_equal(A, A.T):
print('input is symmetric')
else:
print('input is NOT symmetric')
methods = {
'np.linalg.eig' : np.linalg.eig,
'la.eig' : la.eig,
'np.linalg.eigh' : np.linalg.eigh,
'la.eigh' : la.eigh
}
for name, method in methods.items():
print('============================================================')
print(name)
print()
eigenValues, eigenVectors = method(A)
for i in range(len(eigenValues)):
print('{0:6.3f}{1:+6.3f}i '.format(eigenValues[i].real, eigenValues[i].imag), end=' | ')
line = eigenVectors[i]
for item in line:
print('{0:6.3f}{1:+6.3f}i '.format(item.real, item.imag), end='')
print()
print('---------------------')
for i in range(len(eigenValues)):
if eigenValues[i].imag == 0:
print('real ', end=' | ')
else:
print('COMPLEX ', end=' | ')
line = eigenVectors[i]
for item in line:
if item.imag == 0:
print('real ', end='')
else:
print('COMPLEX ', end='')
print()
print()
Run Code Online (Sandbox Code Playgroud)
如您所见,numpy.linalg.eig并scipy.linalg.eig在其输出中生成复数,但它们不应该。如果虚部的大小与实部的大小相比微不足道,则可以将其视为某种舍入误差。但这种情况并非如此。产生的数字之一是-0.013+0.016i。这里虚部的量级甚至高于实部。
更糟糕的是:这四种方法会产生三种不同的结果。
所有四种方法仅计算一次的20的特征值12倍和7倍和所有的特征向量的特征值总是具有长度1。这意味着,所有的四种方法应该产生非常相同的特征向量的特征值用于12.但是,只有numpy.linalg.eig与scipy.linalg.eig产生相同的输出.
以下是特征值 12 的特征向量的分量。仔细查看标有箭头 ( <==)的线。在这里您可以找到三个不同的值,但这些值应该完全相等。如果你有第二次看,你会看到,只有在1日线的所有三个值相等。在所有其他行中,您会发现 2 或 3 个不同的值。
input is symmetric
============================================================
np.linalg.eig
12.000+0.000i | -0.354+0.000i 0.913+0.000i 0.204+0.000i -0.013+0.016i -0.013-0.016i 0.160+0.000i -0.000+0.000i 0.130+0.000i
20.000+0.000i | -0.354+0.000i -0.183+0.000i 0.208+0.000i 0.379-0.171i 0.379+0.171i -0.607+0.000i 0.000+0.000i -0.138+0.000i
20.000+0.000i | -0.354+0.000i -0.182+0.000i 0.203+0.000i -0.468-0.048i -0.468+0.048i 0.153+0.000i 0.001+0.000i -0.271+0.000i
20.000+0.000i | -0.354+0.000i -0.182+0.000i 0.203+0.000i 0.657+0.000i 0.657-0.000i 0.672+0.000i -0.001+0.000i 0.617+0.000i
20.000-0.000i | -0.354+0.000i -0.182+0.000i 0.203+0.000i -0.276+0.101i -0.276-0.101i -0.361+0.000i 0.001+0.000i -0.644+0.000i
20.000+0.000i | -0.354+0.000i -0.001+0.000i -0.612+0.000i -0.001+0.000i -0.001-0.000i 0.001+0.000i 0.706+0.000i -0.000+0.000i
20.000+0.000i | -0.354+0.000i -0.182+0.000i 0.203+0.000i -0.276+0.101i -0.276-0.101i -0.018+0.000i -0.000+0.000i 0.306+0.000i
20.000+0.000i | -0.354+0.000i -0.001+0.000i -0.612+0.000i -0.001+0.000i -0.001-0.000i 0.001+0.000i -0.708+0.000i 0.000+0.000i
---------------------
real | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
COMPLEX | real real real real real real real real
COMPLEX | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
============================================================
la.eig
12.000+0.000i | -0.354+0.000i 0.913+0.000i 0.204+0.000i -0.013+0.016i -0.013-0.016i 0.160+0.000i -0.000+0.000i 0.130+0.000i
20.000+0.000i | -0.354+0.000i -0.183+0.000i 0.208+0.000i 0.379-0.171i 0.379+0.171i -0.607+0.000i 0.000+0.000i -0.138+0.000i
20.000+0.000i | -0.354+0.000i -0.182+0.000i 0.203+0.000i -0.468-0.048i -0.468+0.048i 0.153+0.000i 0.001+0.000i -0.271+0.000i
20.000+0.000i | -0.354+0.000i -0.182+0.000i 0.203+0.000i 0.657+0.000i 0.657-0.000i 0.672+0.000i -0.001+0.000i 0.617+0.000i
20.000-0.000i | -0.354+0.000i -0.182+0.000i 0.203+0.000i -0.276+0.101i -0.276-0.101i -0.361+0.000i 0.001+0.000i -0.644+0.000i
20.000+0.000i | -0.354+0.000i -0.001+0.000i -0.612+0.000i -0.001+0.000i -0.001-0.000i 0.001+0.000i 0.706+0.000i -0.000+0.000i
20.000+0.000i | -0.354+0.000i -0.182+0.000i 0.203+0.000i -0.276+0.101i -0.276-0.101i -0.018+0.000i -0.000+0.000i 0.306+0.000i
20.000+0.000i | -0.354+0.000i -0.001+0.000i -0.612+0.000i -0.001+0.000i -0.001-0.000i 0.001+0.000i -0.708+0.000i 0.000+0.000i
---------------------
real | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
COMPLEX | real real real real real real real real
COMPLEX | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
real | real real real COMPLEX COMPLEX real real real
============================================================
np.linalg.eigh
12.000+0.000i | -0.354+0.000i 0.000+0.000i 0.000+0.000i -0.086+0.000i 0.905+0.000i -0.025+0.000i 0.073+0.000i 0.205+0.000i
20.000+0.000i | -0.354+0.000i 0.000+0.000i -0.374+0.000i 0.149+0.000i -0.236+0.000i -0.388+0.000i 0.682+0.000i 0.206+0.000i
20.000+0.000i | -0.354+0.000i 0.001+0.000i 0.551+0.000i 0.136+0.000i -0.180+0.000i 0.616+0.000i 0.317+0.000i 0.201+0.000i
20.000+0.000i | -0.354+0.000i 0.001+0.000i -0.149+0.000i 0.719+0.000i -0.074+0.000i -0.042+0.000i -0.534+0.000i 0.207+0.000i
20.000+0.000i | -0.354+0.000i -0.005+0.000i 0.505+0.000i -0.386+0.000i -0.214+0.000i -0.556+0.000i -0.274+0.000i 0.203+0.000i
20.000+0.000i | -0.354+0.000i -0.707+0.000i -0.004+0.000i 0.002+0.000i 0.001+0.000i 0.002+0.000i -0.000+0.000i -0.612+0.000i
20.000+0.000i | -0.354+0.000i 0.003+0.000i -0.529+0.000i -0.535+0.000i -0.203+0.000i 0.398+0.000i -0.262+0.000i 0.203+0.000i
20.000+0.000i | -0.354+0.000i 0.707+0.000i 0.001+0.000i 0.001+0.000i 0.000+0.000i -0.005+0.000i -0.001+0.000i -0.612+0.000i
---------------------
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
============================================================
la.eigh
12.000+0.000i | -0.354+0.000i 0.000+0.000i 0.000+0.000i -0.225+0.000i 0.882+0.000i 0.000+0.000i 0.065+0.000i -0.205+0.000i
20.000+0.000i | -0.354+0.000i 0.000+0.000i -0.395+0.000i 0.332+0.000i -0.156+0.000i 0.227+0.000i 0.701+0.000i -0.205+0.000i
20.000+0.000i | -0.354+0.000i 0.001+0.000i 0.612+0.000i 0.011+0.000i -0.204+0.000i -0.597+0.000i 0.250+0.000i -0.200+0.000i
20.000+0.000i | -0.354+0.000i 0.001+0.000i -0.086+0.000i 0.689+0.000i 0.030+0.000i -0.054+0.000i -0.589+0.000i -0.205+0.000i
20.000+0.000i | -0.354+0.000i -0.005+0.000i 0.413+0.000i -0.264+0.000i -0.245+0.000i 0.711+0.000i -0.165+0.000i -0.205+0.000i
20.000+0.000i | -0.354+0.000i -0.707+0.000i -0.004+0.000i -0.000+0.000i 0.001+0.000i -0.002+0.000i -0.001+0.000i 0.612+0.000i
20.000+0.000i | -0.354+0.000i 0.003+0.000i -0.540+0.000i -0.542+0.000i -0.309+0.000i -0.290+0.000i -0.261+0.000i -0.205+0.000i
20.000+0.000i | -0.354+0.000i 0.707+0.000i 0.001+0.000i -0.000+0.000i 0.001+0.000i 0.005+0.000i -0.001+0.000i 0.612+0.000i
---------------------
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
real | real real real real real real real real
Run Code Online (Sandbox Code Playgroud)
以下是我的问题:
ps:以下是相关版本信息:
这是numpy.show_config()
(根据评论中的要求)的输出:
numpy.linalg.eig | |
scipy.linalg.eig | numpy.linalg.eigh | scipy.linalg.eigh
------------------+---------------------+-------------------
-0.354+0.000i | -0.354+0.000i | -0.354+0.000i
0.913+0.000i | 0.000+0.000i | 0.000+0.000i
0.204+0.000i | 0.000+0.000i | 0.000+0.000i
-0.013+0.016i | -0.086+0.000i | -0.225+0.000i <===
-0.013-0.016i | 0.905+0.000i | 0.882+0.000i <===
0.160+0.000i | -0.025+0.000i | 0.000+0.000i <===
-0.000+0.000i | 0.073+0.000i | 0.065+0.000i <===
0.130+0.000i | 0.205+0.000i | -0.205+0.000i
Run Code Online (Sandbox Code Playgroud)
对评论的反应:
实对称矩阵的复特征向量
@Rethipher:谢谢!我确实阅读并理解了您链接到的问题(真正的对称矩阵可以有复杂的特征向量吗?),我也确实阅读了答案,但我不明白。他们说“是”还是“不是”?(修辞问题,无需回答,见下一行)
@Mark Dickinson & @bnaecker:谢谢你澄清,我的假设是错误的。
实对称矩阵与 Hermitian 矩阵
@bnaecker:实数集是复数集的子集。那些等于它们自己的复共轭的复数称为实数。因此,实对称矩阵集是 Hermitian 矩阵的子集。这很重要,因为 numpy.linalg.eigh 和 scipy.linalg.eigh 旨在处理 Hermitian 矩阵。因为每个实对称矩阵都是厄米矩阵,所以这些模块也可以用于我的目的。
混淆行和列
@Mark Dickinson & @bnaecker:谢谢,我认为你是对的。文档也这么说,我应该更仔细地阅读它。但即使您比较的是列而不是行,您仍然会发现 4 种方法会产生 3 种不同的结果。但是如果结果包含一个只能用 7 个实基向量来描述的 7 维子空间,我仍然觉得很奇怪,一个算法会产生一个复基。
“一个错误会令人惊讶”
@bnaecker:这是真的,但确实存在令人惊讶的错误。(比如Heartbleed和其他一些。)所以,这不是一个真正的论点。
“我得到实数” - “你的样本矩阵不包含浮点数”
@Stef & @JohanC:抱歉,你没有足够仔细地阅读我的程序。在计算特征值和特征向量之前,我添加了一个1e-12toA[5,7]和A[7,5]来模拟在我的真实应用程序中不可避免地出现的微小舍入误差。(我在这里发布的只是一个很小的测试程序,大到足以证明这个问题。)
你是对的,Stef:不添加这个微小的噪音,我也得到了真实的结果。但是只有百万分之一百万分之一的微小变化会产生如此大的差异,我不明白为什么。
对 DavidB2013 回答的反应:
我尝试了您建议的工具,但得到了不同的结果。我想你也忘了添加1e-12toA[5,7]和 的小噪音A[7,5]。但是,所有结果仍然是真实的。我确实得到了这些特征值:
12.000000000000249
20
20.00000000000075
19.999999999999
20
20
20
20
Run Code Online (Sandbox Code Playgroud)
和这些特征向量:
0.3535533905932847 0.9128505045937204 0.20252576206455747 0.002673672081814904 -0.09302397289286794 -0.09302397289286794 -0.09302397289286794 -0.09302397289286794
0.3535533905932848 -0.18259457246238117 0.20444330131542393 -0.00009386949436945406 -0.20415317121194954 -0.20415317121194954 -0.20415317121194954 -0.20415317121194954
0.3535533905932848 -0.18259457246238117 0.20444330131542393 -0.00009386949436945406 -0.20415317121194954 -0.20415317121194954 -0.20415317121194954 0.9080920678356449
0.3535533905932848 -0.18259457246238117 0.20444330131542393 -0.00009386949436945406 -0.20415317121194954 0.9080920678356449 -0.20415317121194954 -0.20415317121194954
0.3535533905932848 -0.18259457246238117 0.20444330131542393 -0.00009386949436945406 0.9080920678356449 -0.20415317121194954 -0.20415317121194954 -0.20415317121194954
0.35355339059324065 -0.00011103276380548543 -0.6116010247648269 0.7060012169461334 0.0005790869815273477 0.0005790869815273477 0.0005790869815273477 0.0005790869815273477
0.3535533905932848 -0.18259457246238117 0.20444330131542393 -0.00009386949436945406 -0.20415317121194954 -0.20415317121194954 0.9080920678356449 -0.20415317121194954
0.35355339059324054 0.0002333904819895115 -0.6131412438770024 -0.7082055415560993 0.0009655029234935232 0.0009655029234935232 0.0009655029234935232 0.0009655029234935232
Run Code Online (Sandbox Code Playgroud)
只有特征值 12 的向量与您的计算具有相同的值:(在 6 维中大约有 1.1e-14 的差异,在其他两个维度中大约有 3.3e-14 的差异,但我将此视为舍入误差。)所有其他向量显着不同(最小差异的大小为 0.02)。让我感到困惑的是,输入矩阵的 2 个元素中 1e-12 的微小舍入误差会产生如此大的差异。
我用另一种方法计算了特征值(在https://www.wolframalpha.com的帮助下),当我没有添加应该模拟舍入误差的微小增量值时,我只得到两个不同的特征值,它们是 12和 20。
给定矩阵的特征多项式为:
(20 - ?)^7 * (12 - ?)
Run Code Online (Sandbox Code Playgroud)
因此,它在 ?=12 处有一个根,在 ?=20 处有 7 个根,这 8 个根是 8 个特征值。都是实数。
当我添加微小的 delta 值时,我得到了这个特征多项式:
(20 - ?)^5 * (19999999999999/1000000000000 - ?) * (1000000000000 ?^2 - 32000000000001 ? + 240000000000014)/1000000000000
Run Code Online (Sandbox Code Playgroud)
它有这些根源:
?=12.00000000000024999999999998 (rounded)
?=19.999999999999 (exact value)
?=20 (exact value)
?=20 (exact value)
?=20 (exact value)
?=20 (exact value)
?=20 (exact value)
?=20.00000000000075000000000002 (rounded)
Run Code Online (Sandbox Code Playgroud)
同样,所有 8 个特征值都是实数。
然后我计算了特征向量。不添加 1e-12 我得到这个结果:
特征值 12 的向量:
v = (1,1,1,1,1,1,1,1)
Run Code Online (Sandbox Code Playgroud)
这个向量的长度是 sqrt(8),如果你把这个向量乘以 1/sqrt(8),你就会得到其他计算的结果(每个维度为 0.35355339)。
但是特征值 20 的七个特征向量非常不同。他们是:
(-1,1,0,0,0,0,0,0)
(-1,0,1,0,0,0,0,0)
(-1,0,0,1,0,0,0,0)
(-1,0,0,0,1,0,0,0)
(-1,0,0,0,0,1,0,0)
(-1,0,0,0,0,0,1,0)
(-1,0,0,0,0,0,0,1)
Run Code Online (Sandbox Code Playgroud)
即使将它们带到长度 1,它们也与所有其他结果不同,很容易看出它们是正确的。其他结果也是正确的,但我更喜欢这些简单的结果。
我还计算了带有微小噪声的版本的特征值。所有 8 个向量都非常接近无噪声结果,甚至 Wolfram Alpha 也将它们四舍五入到与以前完全相同的值。这正是我对计算特征值和特征向量的算法所期望的行为:
小智 2
据我所知,假设1是正确的,但假设2不正确。
实数对称矩阵产生仅实数的特征值和特征向量。
然而,对于给定的特征值,关联的特征向量不一定是唯一的。
此外,对于实际上不是那么大或包含不是很小的数字的矩阵,舍入误差不应该那么大。
为了进行比较,我通过 JavaScript 版本的 RG.F(Real General,来自 EISPACK 库)运行了您的测试矩阵: 特征值和特征向量计算器
这是输出:
特征值:
20
12
20
20
20
20
20
20
Run Code Online (Sandbox Code Playgroud)
特征向量:
0.9354143466934854 0.35355339059327395 -0.021596710639534 -0.021596710639534 -0.021596710639534 -0.021596710639534 -0.021596710639533997 -0.021596710639533997
-0.1336306209562122 0.3535533905932738 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797
-0.1336306209562122 0.3535533905932738 0.9286585574999623 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797
-0.1336306209562122 0.3535533905932738 -0.15117697447673797 0.9286585574999623 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797
-0.1336306209562122 0.3535533905932738 -0.15117697447673797 -0.15117697447673797 0.9286585574999623 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797
-0.1336306209562122 0.3535533905932738 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 0.9286585574999623 -0.15117697447673797 -0.15117697447673797
-0.1336306209562122 0.3535533905932738 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 0.9286585574999622 -0.15117697447673797
-0.1336306209562122 0.3535533905932738 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 -0.15117697447673797 0.9286585574999622
Run Code Online (Sandbox Code Playgroud)
没有虚部。
为了确认或否认结果的有效性,您始终可以编写一个小程序,将结果插回到原始方程中。简单的矩阵和向量乘法。然后你就可以确定输出是否正确。或者,如果他们错了,他们离正确答案有多远。
| 归档时间: |
|
| 查看次数: |
326 次 |
| 最近记录: |