计算两个通道中复杂数据的协方差矩阵（无复杂数据类型）

Question

我在矩阵的 2 个通道中给出了复值数据（一个是实部，一个是虚部，所以矩阵维度是 (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-d 数组，并希望它能够工作。

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

Answer 1

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

---

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 是复数值，即由H 行变量W 复数观察组成。 此外，令A 和B 为实值矩阵，使得M = A + i.B。

我不会进行数学演示，感谢@zimzam，您可以找到here，但在这种情况下，cov(M) 可以分解为：

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

鉴于M 的实部和虚部（A 和B），这使得单独计算cov(M) 的实部和虚部变得简单。

---

实施与验证

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

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