我需要从资产类别返回中建模并估算一个差异 - 可达矩阵，所以我正在查看https://github.com/camdavidsonpilon/的第6章中给出的股票返回示例概率编程和bayesian-methods-for-backers
这是我的简单实现，在该实现中，我使用带有已知均值和方差互动矩阵的多元正常人开始样本。然后，我尝试使用非信息性冠军对其进行估算。
估计与已知的先验不同，因此我不确定我的实现是否正确。如果有人可以指出我在做什么错，我会很感激？
import numpy as np
import pandas as pd
import pymc as pm
p=3
mu=[.03,.05,-.02]
cov_matrix= [[.025,0.0075, 0.00175],[0.0075,.007,0.00135],[0.00175,0.00135,.00043]]
n_obs=10000
x=np.random.multivariate_normal(mu,cov_matrix,n_obs)
prior_mu=np.ones(p)
prior_sigma = np.eye(p)
post_mu = pm.Normal("returns",prior_mu,1,size=p)
post_cov_matrix_inv = pm.Wishart("cov_matrix_inv",n_obs,np.linalg.inv(cov_matrix))
obs = pm.MvNormal( "observed returns", post_mu, post_cov_matrix_inv, observed = True, value = x )
model = pm.Model( [obs, post_mu, post_cov_matrix_inv] )
mcmc = pm.MCMC()
mcmc.sample( 5000, 2000, 3 )
mu_samples = mcmc.trace("returns")[:]
mu_samples.mean(axis=0)
cov_inv_samples = mcmc.trace("cov_matrix_inv")[:]
mean_covariance_matrix = np.linalg.inv( cov_inv_samples.mean(axis=0) )
这是我提出的一些建议，可以改善代码 推理：
-
我会
pm.Wishart("cov_matrix_inv",n_obs,np.linalg.inv(cov_matrix))pm.Wishart("cov_matrix_inv",n_obs,np.eye(3) )是更客观的（并且有10000个数据点，您的先验都不会很重要） -
mcmc = pm.MCMC()应该是mcmc = pm.MCMC(model) -
mcmc.sample( 5000, 2000, 3 )这里有很少的样本。MCMC的后半部分是蒙特卡洛（Monte Carlo），当有很多样本时最强：我的意思是成千上万。在这里，您只有1000，因此由蒙特卡洛引起的误差将很高（随着样本量的增加而误差降低）。此外，在2000年燃烧样品后，MCMC可能没有融合。您可以在pymc.Matplot中检查与plot的收敛性并调用plot(mcmc)。我使用了mcmc.sample( 25000, 15000, 1 )，并取得了更好的结果。
我想您使用如此低样本的原因是性能。其中大部分是由大量样本引起的：您有10000个观测值。在实践中，这可能很高。
记住，贝叶斯推论的大部分价值都得到了后验样本：以这些样本的平均值似乎是浪费 - 考虑在损失功能中使用样品（请参阅书中的第5章）。
请注意，如果您想使用信息性之前，则不应将np.linalg.inv(cov_matrix)用于Wishart，而应仅使用cov_matrix。确切地说，您应该使用cov_matrix * n_obs来正确缩放