从多变量正态分布中采样多个元素答案

【问题标题】：Sample more than one element from multivariable normal distribution从多变量正态分布中采样多个元素
【发布时间】：2021-01-26 11:30:46
【问题描述】：

我有一个大小为n*m 的二维均值矩阵，其中n 是样本数，m 是数据的维度。我也有nm*m 的矩阵，即 sigma 是我的形状为n*m*m 的方差矩阵。我希望从上述分布中抽取n 样本，例如x_i~N(mean[i], sigma[i])。在numpy 或任何其他不使用 for 循环运行的标准库中有什么方法可以做到这一点？

我认为唯一的选择是使用np.random.multivariate_normal()，将均值矩阵展平为一个向量，并将 3D sigma 展平为 2D 块对角矩阵。当然，之后会重塑。但这意味着我们将使用形状为 (n*m)*(n*m) 的 sigma 的样本，这很容易变得非常大，并且仅计算和分配该矩阵（如果可能）可能比在 for 循环中运行花费更长的时间。

在我的具体任务中，现在 Sigma 是所有样本的相同矩阵，这意味着我可以用m*m 表示 Sigma，并且所有 n 个点都是相同的。但我对通用解决方案感兴趣。

感谢您的帮助。

【问题讨论】：

标签： python numpy matrix distribution normal-distribution

【解决方案1】：

没有可测试的代码很难判断，但这应该很接近：

A = numpy.linalg.cholesky(sigma)                # => shape (n, m, m), same as sigma
Z = np.random.normal(size = (n, m))             # shape (n, m)
X = np.einsum('ijk, ik -> ij', A, Z) + mean     # shape (n, m)

发生了什么：

我们根据here 概述的标准 Cholesky 分解方法手动采样多元正态分布。 A 的构建使得 A@A.T = sigma。那么X（多元法线）可以由A和一元法线N(0, 1)向量Z的点积加上mean形成。 p>

在第一个（索引 = 0，'i'）轴中的第一个（索引 = 0，'i'）中保持无关维度，同时收缩最后一个（'k'）轴，形成点积。

【讨论】：

谢谢！我想知道为什么 numpy 不在 np.random.multivariate_random() 的情况下这样做 dim(mean) != 1 而不是抛出错误。
因为它在内部构建为X = A@Z + mean，除非一切都保持在 1 和 2 空间中，否则它不起作用。我讨厌作为在高维状态空间工作的人，但是由于numpy 的绝大多数用户群不知道任何高于 2 空间（矩阵）的张量数学，所以一切都被简化了到 2 空间并做任何更高-d 的事情留给用户使用 einsum 手动处理
顺便说一句，我有一个类似的问题（也许应该打开一个新问题），如果我有自己的自定义分布（一维数组不是负数，总和为 1）并且我想从中采样，不np.random.choice() 的问题。如果我有一个分布矩阵，即每一列都是一个分布向量，并且我想从每个向量中进行采样，那么没有 for 循环的任何方法吗？
请再写一个问题，这次尝试包括minimal reproducible example。我不是特别擅长读心术。