python - 计算两个 channel 中复杂数据的协方差矩阵(无复杂数据类型)

Question

我在矩阵的 2 个 channel 中给出了复数值数据(一个是实部，一个是虚部，所以矩阵维度是 (height, width, 2)，因为 Pytorch 确实没有原生的复杂数据类型。我现在想计算协方差矩阵。适用于 Pytorch 的精简 numpy 计算是这样的:

def cov(m, y=None):

    if m.ndimension() > 2:
        raise ValueError("m has more than 2 dimensions")

    if y.ndimension() > 2:
        raise ValueError('y has more than 2 dimensions')

    X = m
    if X.shape[0] == 0:
        return torch.tensor([]).reshape(0, 0)
    if y is not None:
        X = torch.cat((X, y), dim=0)

    ddof = 1

    avg = torch.mean(X, dim=1)

    fact = X.shape[1] - ddof

    if fact <= 0:
        import warnings
        warnings.warn("Degrees of freedom <= 0 for slice",
                      RuntimeWarning, stacklevel=2)
        fact = 0.0

    X -= avg[:, None]
    X_T = X.t()
    c = dot(X, X_T)
    c *= 1. / fact
    return c.squeeze()

现在在 numpy 中，这将透明地处理复数，但我不能简单地提供最后一个维度为 (real, imag) 的 3 维数组并希望它会起作用。

如何调整计算以获得具有实部和虚部 channel 的复数协方差矩阵？

Answer 1

[对于复杂矩阵的 cov() 的 PyTorch 实现，跳过解释并转到最后一个片段]

---

cov()操作说明

设M是一个HxW矩阵，其中H行的每一行对应一个W的变量复杂的观察。

现在令 cov(M) 为 M 的 H 变量的 HxH 协方差矩阵(定义的 numpy.cov() )。它可以按如下方式计算(忽略边缘情况):

cov(M) = 1 / (W - 1) . M * M.T

使用*矩阵乘法运算符，以及M的M.T转置。

注意:为了阐明下一个等式，令 cov_prod(X, Y) = 1/(W - 1) 。 X * Y.T，具有 X, Y HxW 矩阵。因此我们有 cov(M) = cov_prod(M, M)。

到目前为止，没有什么新内容，这与您编写的代码相对应(减去 y 权重和边缘情况的数据检查)。让我们仔细检查一下这个公式的 Pytorch 实现是否对应于 Numpy 的实值数据:

import torch
import numpy as np

def cov(m, y=None):
    if y is not None:
        m = torch.cat((m, y), dim=0)
    m_exp = torch.mean(m, dim=1)
    x = m - m_exp[:, None]
    cov = 1 / (x.size(1) - 1) * x.mm(x.t())
    return cov

# Real-valued matrix:
M_np = np.random.rand(3, 2)  
# Same matrix as torch.Tensor:
M = torch.from_numpy(M_np)

cov_real_np = np.cov(M_np)
cov_real = cov(M)
eq = np.allclose(cov_real_np, cov_real.numpy())
print("Numpy & Torch real covariance results equal? > {}".format(eq))
# Numpy & PyTorch real covariance results equal? > True

---

复杂矩阵的扩展

现在，这对复杂矩阵有何作用？从这里开始，令 M 为复值，即由 W 复杂观察的 H 行变量组成。此外，设 A 和 B 为实值矩阵，满足 M = A + i.B。

我不会进入数学演示，你可以找到here感谢@zimzam，但在那种情况下 cov(M) 可以分解为:

cov(M) = [cov_prod(A, A) + cov_prod(B, B)] + i.[-cov_prod(A, B) + cov_prod(B, A)]

给定 M 的实部和虚部(A 和 B)。

---

实现与验证

在下面找到一个优化的实现:

import torch
import numpy as np

def cov_complex(m_comp):
    # (adding further parameters such as `y` is left for exercise)
    # Supposing real and img are stored separately in the last dim:
    real, img = m_comp[..., 0], m_comp[..., 1] 
    x_real = real - torch.mean(real, dim=1)[:, None]
    x_img = img - torch.mean(img, dim=1)[:, None]
    x_real_T = x_real.t()
    x_img_T = x_img.t()
    frac = 1 / (x_real.size(1) - 1)

    cov_real = frac * (x_real.mm(x_real_T) + x_img.mm(x_img_T))
    cov_img = frac * (-x_real.mm(x_img_T) + x_img.mm(x_real_T))
    return torch.stack((cov_real, cov_img), dim=-1)

# Matrix with real/img values stored separately in last dimension:
M_np = np.random.rand(3, 2, 2)  
# Same matrix converted to np.complex format:                  
M_comp_np = M_np.view(dtype=np.complex128)[...,0]
# Same matrix as torch.Tensor:
M = torch.from_numpy(M_np)

cov_com_np = np.cov(M_comp_np)
cov_com = cov_complex(M)
eq = np.allclose(cov_com_np, cov_com.numpy().view(dtype=np.complex128)[...,0])
print("Numpy & Torch complex covariance results equal? > {}".format(eq))
# Numpy & PyTorch complex covariance results equal? > True