如何用 MCMC (MHadaptive) 拟合多元正态分布

问题描述

我想通过 MCMC 采样来估计多元正态分布的参数。棘手的部分是协方差矩阵。我知道传统的先验选择是逆愿望先验（或精度矩阵的愿望先验）。但是，我不知道该给谁使用。当我尝试估计多元正态分布的参数时，MHadaptive 会崩溃，因为 sigma 不是正定的。我也用不同的采样器尝试过，但最终还是遇到了同样的问题。以下是我迄今为止尝试过的其他一些事情：

使用带有 cholesky 参数化的逆心愿图并估计 U 的上三角。然后sigma=t(U)%*%U。这不会导致正定问题，但它也不起作用。

改用精度矩阵，这也没有解决问题。

当然，我也尝试过只估计Sigma的上三元并重建Sigma，使其对称。然而，它仍然不是正定的。

这是我的代码：

library(mvtnorm)
library(MHadaptive)

#data
x <- rmvnorm(100,c(20,65,-44),diag(3))

#define the log prior
logprior <- function(parvect){
    logprior <- LaplacesDemon::dinvwishart(matrix(parvect[1:9],3,byrow=TRUE),diag(3),log=TRUE) +
      dnorm(parvect[10],100000,log=TRUE) + 
      dnorm(parvect[11],log=TRUE) + 
      dnorm(parvect[12],log=TRUE) 
    return(logprior)
}

#define the log likelihood
LL <- function(parvect,data){
  LL <- sum(mvtnorm::dmvnorm(data,parvect[10:12],matrix(parvect[1:9],log=TRUE))

  return(LL) 
}

LL_reg <- function(parvect,data){
  LL <- LL(parvect = parvect,data=data)
  logPrior <- logprior(parvect = parvect)
  LL_reg <- LL+logPrior
  return(LL_reg)
}

#inits etc.
df_params <- data.frame(name = c( "covmat[1,1]","covmat[1,2]",3]","covmat[2,"covmat[3,"mean[1]","mean[2]","mean[3]"),min = c(rep(-Inf,12)),max = c(rep(Inf,init=c(as.vector(t(cov(x))),colMeans(x)))


Metro_Hastings(li_func = LL_reg,pars = df_params$init,par_names = df_params$name,data=x)

解决方法

既然你知道似然和先验，为什么不试试 Rstan 来得到估计呢？定义 cov_matrix 或直接将方差-协方差矩阵放在多正态函数的方差部分会限制为 p.d.

bayesian mcmc