MATLAB与Python（PyTorch/NumPy）矩阵运算及特征向量结果差异问询

问题背景 #

将MATLAB代码迁移至Python时，计算矩阵P与Q的乘积PQ的特征值和特征向量，出现了MATLAB、PyTorch、NumPy结果不一致的情况。

矩阵定义 #

P = [0.153833459882762  0.733552631578947;
     -0.911184210526316 -0.153833459882762];
Q = [0.649519052838329 3.09722222222222;
    -3.84722222222222 -0.649519052838329];

各平台计算结果 #

MATLAB结果

>> PQ = P*Q
PQ =
         -2.72222222222222                         0
     -1.11022302462516e-16         -2.72222222222222

>> [eigvec, eigval] = eig(PQ)
eigvec =
                         0         0.983547166232033
                         1         0.180651520328333
eigval =
         -2.72222222222222                         0
                         0         -2.72222222222222

PyTorch结果

import torch
P = torch.tensor([[0.153833459882762, 0.733552631578947],
     [-0.911184210526316, -0.153833459882762]],dtype = torch.complex64)
Q = torch.tensor([[0.649519052838329, 3.09722222222222],
    [-3.84722222222222, -0.649519052838329]],dtype = torch.complex64)
PQ = P@Q
print(PQ)
torch.linalg.eig(PQ)

输出：

tensor([[-2.7222e+00+0.j,  2.9802e-08+0.j],
        [ 1.1921e-07+0.j, -2.7222e+00+0.j]])
torch.return_types.linalg_eig(
eigenvalues=tensor([-2.7222+0.j, -2.7222+0.j]),
eigenvectors=tensor([[0.-0.j, 1.+0.j],
        [1.+0.j, 0.-0.j]]))

NumPy结果

import numpy as np
np.linalg.eig(P@Q)

输出：

(array([-2.7222+0.j, -2.7222+0.j], dtype=complex64),
 array([[ 0.4472-0.j, -0.4472+0.j],
        [ 0.8944+0.j,  0.8944+0.j]], dtype=complex64))

疑问 #

MATLAB与Python的矩阵乘积结果存在差异，且MATLAB、PyTorch、NumPy的特征向量完全不同。推测差异源于浮点运算的数值精度，但结果差异过大。请问：

1. PyTorch与NumPy特征向量不同的原因是什么？
2. 有无方法让Python结果与MATLAB对齐？

解答 #

1. PyTorch与NumPy特征向量不同的原因 #

• 特征向量的非唯一性：对于重特征值（此处特征值均为-2.7222，是二重根），对应的特征向量构成一个线性子空间，而非唯一向量。该子空间内的任意线性组合都是合法的特征向量。
  • NumPy返回的是该子空间的一组正交基（近似[1,2]和[-1,2]方向的归一化向量）。
  • PyTorch返回的是标准基向量[0,1]和[1,0]，这两个向量同样属于该特征子空间，满足PQ * v = λv的条件。
• 底层实现与算法差异：NumPy和PyTorch依赖不同的线性代数库：
  • NumPy默认使用LAPACK库的geev函数。
  • PyTorch根据设备（CPU/GPU）选择不同实现：CPU上可能调用MKL或OpenBLAS，GPU上使用CUDA线性代数库。不同算法处理重特征值时，会选择不同的基向量表示特征子空间。
• 浮点精度与数值稳定性：矩阵PQ的非对角元存在微小浮点误差（MATLAB中为1e-16，PyTorch中为1e-7~1e-8），这些误差会影响算法对特征子空间基向量的选择，导致最终返回的特征向量形式不同。

2. 让Python结果与MATLAB对齐的方法 #

对齐矩阵乘积结果

MATLAB默认使用双精度（double）计算，而Python代码中使用了complex64（单精度复数），这是乘积误差差异的主要原因。将张量/数组类型改为双精度：

• PyTorch：将dtype=torch.complex64改为dtype=torch.complex128（若无需复数，可使用torch.float64）。
• NumPy：使用dtype=np.complex128（或np.float64），计算后的PQ非对角元会接近MATLAB的1e-16量级。

对齐特征向量

由于特征向量的非唯一性，无法强制让Python返回与MATLAB完全相同的向量，但可通过以下方式对齐行为：

• 归一化方式匹配：MATLAB的eig函数默认返回列向量归一化（每个特征向量的2-范数为1）的结果。对Python返回的特征向量做同样处理：

  # NumPy示例：归一化特征向量列
  vals, vecs = np.linalg.eig(PQ.astype(np.complex128))
  vecs_normalized = vecs / np.linalg.norm(vecs, axis=0)
  

• 调用MATLAB引擎：在Python中直接调用MATLAB的eig函数，保证结果完全一致：

  import matlab.engine
  eng = matlab.engine.start_matlab()
  P_mat = matlab.double([[0.153833459882762, 0.733552631578947],
                         [-0.911184210526316, -0.153833459882762]])
  Q_mat = matlab.double([[0.649519052838329, 3.09722222222222],
                         [-3.84722222222222, -0.649519052838329]])
  PQ_mat = eng.mtimes(P_mat, Q_mat)
  eigvec_mat, eigval_mat = eng.eig(PQ_mat, nargout=2)
  print(eigvec_mat, eigval_mat)
  eng.quit()
  

• 手动构造特征向量：若业务逻辑依赖特定形式的特征向量，可通过线性组合Python返回的特征向量，构造与MATLAB一致的向量，验证其满足PQ * v = λv即可。

内容的提问来源于stack exchange，提问作者galib