具有分离均值、维度不匹配错误的多元正态模型
Model multivariate normal with separate means, dimension mismatch error
我想用单独的均值和协方差矩阵对多变量(示例中的双变量)法线建模,但不知何故我无法获得匹配的维度。模型比较简单,下面是数据生成过程:
import numpy as np
import pymc3 as pm
import theano.tensor as tt
np.random.seed(123)
N = 200
# x1 and x2 are independent variables
x1 = np.random.randn(N)
x2 = np.random.randn(N)
α1, β1, α2, β2 = np.random.randn(4)
# the means of the joint distribution
μ1 = α1 + β1 * x1
μ2 = α2 + β2 * x2
# the covariance of the joint distribution
σ = 1, 2
Ω = np.array([[1, 0.5], [0.5, 1]])
Σ = np.diag(σ) @ Ω @ np.diag(σ)
# the error term
ϵ = np.random.multivariate_normal(mean=(0, 0), cov=Σ, size=N)
# the dependent variable as bivariate normal
y = np.stack((μ1, μ2)).T + ϵ
这是我的模型:
with pm.Model() as joint_model:
a1 = pm.Normal('a1', 0, 1)
b1 = pm.Normal('b1', 0, 1)
a2 = pm.Normal('a2', 0, 1)
b2 = pm.Normal('b2', 0, 1)
mu1 = a1 + b1 * x1
mu2 = a2 + b2 * x2
sigma = pm.HalfCauchy('sigma', 1, shape=2)
# build the covariance matrix
C_triu = pm.LKJCorr('omega', n=2, p=2)
C = tt.fill_diagonal(C_triu[np.zeros((2, 2), dtype=np.int64)], 1)
sigma_diag = tt.nlinalg.diag(sigma)
cov = tt.nlinalg.matrix_dot(sigma_diag, C, sigma_diag)
joint_obs = pm.MvNormal('joint', mu=(mu1, mu2), cov=cov, observed=y)
错误信息是ValueError: Input dimension mis-match. (input[0].shape[0] = 200, input[1].shape[0] = 2)
。
很明显,我在指定多元正态分布时犯了一个错误,mu
、cov
和观察到的数据之间存在(可能不止一个)维度不匹配,我只是想不通在哪里以及如何修复它。这是(很长的)错误消息,上面的代码应该能够复制错误:
Traceback (most recent call last):
File "<ipython-input-4-b3658e712eaf>", line 18, in <module>
joint_obs = pm.MvNormal('joint', mu=(mu1, mu2), cov=cov, observed=y)
File "/Users/oma/anaconda/lib/python3.6/site-packages/pymc3/distributions/distribution.py", line 39, in __new__
return model.Var(name, dist, data, total_size)
File "/Users/oma/anaconda/lib/python3.6/site-packages/pymc3/model.py", line 545, in Var
total_size=total_size, model=self)
File "/Users/oma/anaconda/lib/python3.6/site-packages/pymc3/model.py", line 970, in __init__
self.logp_elemwiset = distribution.logp(data)
File "/Users/oma/anaconda/lib/python3.6/site-packages/pymc3/distributions/multivariate.py", line 200, in logp
delta = value - mu
File "/Users/oma/anaconda/lib/python3.6/site-packages/theano/tensor/var.py", line 147, in __sub__
return theano.tensor.basic.sub(self, other)
File "/Users/oma/anaconda/lib/python3.6/site-packages/theano/gof/op.py", line 674, in __call__
required = thunk()
File "/Users/oma/anaconda/lib/python3.6/site-packages/theano/gof/op.py", line 843, in rval
fill_storage()
File "/Users/oma/anaconda/lib/python3.6/site-packages/theano/gof/cc.py", line 1698, in __call__
reraise(exc_type, exc_value, exc_trace)
File "/Users/oma/anaconda/lib/python3.6/site-packages/six.py", line 686, in reraise
raise value
ValueError: Input dimension mis-match. (input[0].shape[0] = 200, input[1].shape[0] = 2)
问题已在 discourse.pymc.io 上得到解答。
pm.MvNormal
中 mu
的规格应为 mu = tt.stack([mu1, mu2]).T
。
我想用单独的均值和协方差矩阵对多变量(示例中的双变量)法线建模,但不知何故我无法获得匹配的维度。模型比较简单,下面是数据生成过程:
import numpy as np
import pymc3 as pm
import theano.tensor as tt
np.random.seed(123)
N = 200
# x1 and x2 are independent variables
x1 = np.random.randn(N)
x2 = np.random.randn(N)
α1, β1, α2, β2 = np.random.randn(4)
# the means of the joint distribution
μ1 = α1 + β1 * x1
μ2 = α2 + β2 * x2
# the covariance of the joint distribution
σ = 1, 2
Ω = np.array([[1, 0.5], [0.5, 1]])
Σ = np.diag(σ) @ Ω @ np.diag(σ)
# the error term
ϵ = np.random.multivariate_normal(mean=(0, 0), cov=Σ, size=N)
# the dependent variable as bivariate normal
y = np.stack((μ1, μ2)).T + ϵ
这是我的模型:
with pm.Model() as joint_model:
a1 = pm.Normal('a1', 0, 1)
b1 = pm.Normal('b1', 0, 1)
a2 = pm.Normal('a2', 0, 1)
b2 = pm.Normal('b2', 0, 1)
mu1 = a1 + b1 * x1
mu2 = a2 + b2 * x2
sigma = pm.HalfCauchy('sigma', 1, shape=2)
# build the covariance matrix
C_triu = pm.LKJCorr('omega', n=2, p=2)
C = tt.fill_diagonal(C_triu[np.zeros((2, 2), dtype=np.int64)], 1)
sigma_diag = tt.nlinalg.diag(sigma)
cov = tt.nlinalg.matrix_dot(sigma_diag, C, sigma_diag)
joint_obs = pm.MvNormal('joint', mu=(mu1, mu2), cov=cov, observed=y)
错误信息是ValueError: Input dimension mis-match. (input[0].shape[0] = 200, input[1].shape[0] = 2)
。
很明显,我在指定多元正态分布时犯了一个错误,mu
、cov
和观察到的数据之间存在(可能不止一个)维度不匹配,我只是想不通在哪里以及如何修复它。这是(很长的)错误消息,上面的代码应该能够复制错误:
Traceback (most recent call last):
File "<ipython-input-4-b3658e712eaf>", line 18, in <module>
joint_obs = pm.MvNormal('joint', mu=(mu1, mu2), cov=cov, observed=y)
File "/Users/oma/anaconda/lib/python3.6/site-packages/pymc3/distributions/distribution.py", line 39, in __new__
return model.Var(name, dist, data, total_size)
File "/Users/oma/anaconda/lib/python3.6/site-packages/pymc3/model.py", line 545, in Var
total_size=total_size, model=self)
File "/Users/oma/anaconda/lib/python3.6/site-packages/pymc3/model.py", line 970, in __init__
self.logp_elemwiset = distribution.logp(data)
File "/Users/oma/anaconda/lib/python3.6/site-packages/pymc3/distributions/multivariate.py", line 200, in logp
delta = value - mu
File "/Users/oma/anaconda/lib/python3.6/site-packages/theano/tensor/var.py", line 147, in __sub__
return theano.tensor.basic.sub(self, other)
File "/Users/oma/anaconda/lib/python3.6/site-packages/theano/gof/op.py", line 674, in __call__
required = thunk()
File "/Users/oma/anaconda/lib/python3.6/site-packages/theano/gof/op.py", line 843, in rval
fill_storage()
File "/Users/oma/anaconda/lib/python3.6/site-packages/theano/gof/cc.py", line 1698, in __call__
reraise(exc_type, exc_value, exc_trace)
File "/Users/oma/anaconda/lib/python3.6/site-packages/six.py", line 686, in reraise
raise value
ValueError: Input dimension mis-match. (input[0].shape[0] = 200, input[1].shape[0] = 2)
问题已在 discourse.pymc.io 上得到解答。
pm.MvNormal
中 mu
的规格应为 mu = tt.stack([mu1, mu2]).T
。