我可以使用dask创建multivariate_normal矩阵吗？

问题描述

与this post有关，我正在尝试在multivariate_normal中复制dask： 使用numpy，我可以使用以下方法创建具有指定协方差的多元正态矩阵：

import numpy as np
n_dim = 5
size = 300
A = np.random.randn(n_dim,n_dim) # a matrix
covm = A.dot(A.T) # A*A^T is positive semi-definite,as a covariance matrix
x = np.random.multivariate_normal(size=300,mean=np.zeros(len(covm)),cov=covm) # generate data

但是，我需要一个n_dim = 4_500_000和size = 100000的大矩阵。计算cpu和内存的成本将非常昂贵。幸运的是，我可以访问Cloudera DataScience Workbench集群，并正在尝试使用dask解决此问题：

import dask.array as da
n_dim = 4_500_000
size = 100000
A = da.random.standard_normal((n_dim,n_dim))  
covm = A.dot(A.T)
#x = da.random.multivariate_normal(size=300,cov=covm) # generate data

在documentation中，我找不到任何似乎可以满足我需要的功能。是否有人知道解决方案/工作环境，可能使用xarray或在集群上运行的任何其他模块？

解决方法

目前的一项工作是使用cholesky分解。注意，任何协方差矩阵C都可以表示为C = G * G'。然后，如果y为标准正态，则x = G'* y如C中指定的那样相关（请参见StackExchange Mathematic上的excellent post）。在代码中：

脾气暴躁

n_dim =4
size = 100000
A = np.random.randn(n_dim,n_dim)
covm = A.dot(A.T)

x=  np.random.multivariate_normal(size=size,mean=np.zeros(len(covm)),cov=covm)
## verify numpys covariance is correct
np.cov(x,rowvar=False)
covm

黄昏

## create covariance matrix
A = da.random.standard_normal(size=(n_dim,n_dim),chunks=(2,2))
covm = A.dot(A.T)

## get cholesky decomp
L = da.linalg.cholesky(covm,lower=True)

## drawn standard normal 
sn= da.random.standard_normal(size=(size,chunks=(100,100))

## correct for correlation
x =L.dot(sn.T)
x.shape

## verify
covm.compute()
da.cov(x,rowvar=True).compute()