Python/Numpy：来自多变量分布的条件模拟答案

【问题标题】：Python/Numpy: Conditional simulation from a multivatiate distributionPython/Numpy：来自多变量分布的条件模拟
【发布时间】：2016-12-07 09:55:12
【问题描述】：

使用 numpy 我可以无条件地从多元正态分布中模拟

mean = [0, 0]
cov = [[1, 0], [0, 100]]  # diagonal covariance
x, y = np.random.multivariate_normal(mean, cov, 5000).T

如果我有 5000 个 x 的实现，我如何从同一个分布中模拟 y？我正在寻找可以扩展到任意维度的通用解决方案。

【问题讨论】：

如果这些是独立的（协方差表明是这样）你不能只从单变量正态分布（均值=0，方差=100）生成吗？
确实如此。这只是一个写得不好的例子——我希望能够处理满秩协方差矩阵的情况

标签： python numpy simulation

【解决方案1】：

在伊顿仰望，莫里斯 L. (1983)。多元统计：向量空间方法，我收集了以下示例解决方案，用于 4 变量系统，具有 2 个因变量（前两个）和 2 个自变量（后两个）

import numpy as np

mean = np.array([1, 2, 3, 4])
cov = np.array(
    [[ 1.0,  0.5,  0.3, -0.1], 
     [ 0.5,  1.0,  0.1, -0.2], 
     [ 0.3,  0.1,  1.0, -0.3], 
     [-0.1, -0.2, -0.3,  0.1]])  # diagonal covariance

c11 = cov[0:2, 0:2] # Covariance matrix of the dependent variables
c12 = cov[0:2, 2:4] # Custom array only containing covariances, not variances
c21 = cov[2:4, 0:2] # Same as above
c22 = cov[2:4, 2:4] # Covariance matrix of independent variables

m1 = mean[0:2].T # Mu of dependent variables
m2 = mean[2:4].T # Mu of independent variables

conditional_data = np.random.multivariate_normal(m2, c22, 1000)

conditional_mu = m2 + c12.dot(np.linalg.inv(c22)).dot((conditional_data - m2).T).T
conditional_cov = np.linalg.inv(np.linalg.inv(cov)[0:2, 0:2])

dependent_data = np.array([np.random.multivariate_normal(c_mu, conditional_cov, 1)[0] for c_mu in conditional_mu])

print np.cov(dependent_data.T, conditional_data.T)
>> [[ 1.0012233   0.49592165  0.28053086 -0.08822537]
    [ 0.49592165  0.98853341  0.11168755 -0.22584691]
    [ 0.28053086  0.11168755  0.91688239 -0.27867207]
    [-0.08822537 -0.22584691 -0.27867207  0.94908911]]

可以接受地接近预定义的协方差矩阵。解决方案在Wikipedia上也有简述

【讨论】：

小笨蛋：solve 几乎总是比inv 更受欢迎。
真的吗？谢谢！似乎不合逻辑
耸耸肩。关键词是稳定性。
小错误：conditional_mu = m2 + c12.dot(np.linalg.inv(c22)).dot((conditional_data - m2).T).T 应改为：conditional_mu = m1 + c12.dot(np.linalg.inv(c22)).dot((conditional_data - m2).T).T