如何通过pymc.LKJCorr正确构造相关矩阵？

你的代码错误在于生成的矩阵不满足相关矩阵的对称要求，pm.expand_packed_triangular默认仅填充指定三角区域的元素，其余位置为0，手动加单位矩阵后无法得到对称的合法相关矩阵。

以下是两种正确的实现方式：

方法一：通过三角矩阵转置生成对称相关矩阵 #

import arviz as az
import matplotlib.pyplot as plt
import numpy as np
import pymc as pm

with pm.Model() as model:
    # 获取n=2维度的上三角非对角元素（LKJCorr默认返回上三角packed形式）
    packed_corr = pm.LKJCorr('packed_corr', n=2, eta=2)
    # 展开为上三角矩阵，指定对角元素为1
    upper_tri = pm.expand_packed_triangular(
        n=2, 
        packed=packed_corr, 
        lower=False, 
        diag=np.ones(2)
    )
    # 转置后相加，减去重复的对角元素，得到对称相关矩阵
    corr = pm.Deterministic('corr', upper_tri + upper_tri.T - np.diag(np.diag(upper_tri)))

with model:
    idata = pm.sample()
    az.plot_trace(idata, combined=True)
    plt.tight_layout()
    plt.show()

方法二：利用Cholesky分解生成正定相关矩阵 #

这种方式更严谨，可确保生成的矩阵满足对称、对角为1、正定的相关矩阵核心要求：

import arviz as az
import matplotlib.pyplot as plt
import numpy as np
import pymc as pm

with pm.Model() as model:
    # 获取下三角的packed非对角元素（指定lower=True）
    packed_corr = pm.LKJCorr('packed_corr', n=2, eta=2, lower=True)
    # 展开为下三角矩阵，对角设为1
    lower_tri = pm.expand_packed_triangular(
        n=2, 
        packed=packed_corr, 
        lower=True, 
        diag=np.ones(2)
    )
    # 下三角矩阵乘自身转置，得到正定对称的相关矩阵
    corr = pm.Deterministic('corr', lower_tri @ lower_tri.T)

with model:
    idata = pm.sample()
    az.plot_trace(idata, combined=True)
    plt.tight_layout()
    plt.show()

错误原因说明 #

原代码中pm.expand_packed_triangular默认生成的矩阵对角为0，非指定三角区域的元素也为0，加单位矩阵后得到的是[[1, R_upper], [0, 1]]，这不是对称的合法相关矩阵，因此出现不符合预期的元素分布。

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