Distribution to generate a matrix

Yanjie Fu
  • Yanjie Fu
    Yanjie Fu

    Hi, All,

    Suppose a vector is drawn from a population with multivariate normal distribution, I can use

    vector ~ dmnorm(mean, tau),
    mean ~ dunif(0,100)
    tau ~ power(sigma, -2)
    sigma ~ dunif(0, 100)

    Now, suppose I have a matrix is drawn from a distribution, which distribution in jags can generate a matrix?

    If i want to estimate it, how to do that?

    Last edit: Yanjie Fu 2012-12-18
  • Luc Coffeng
    Luc Coffeng

    The Wishart distribution is the natural choice as this is the conjugate prior for a variance-covariance matrix Sigma2. In JAGS, you should use the inverse-Wishart distribution for the precision matrix Tau. From your example, I gather that you are interested in using an uninformative prior. In terms of an inverse-Wishart distribution, this translates to setting its scale to one plus the number of vectors that you want to produce with your multivariate normal distribution (e.g. in the case of a bivariate outcome, set scale to three), and set the prior expectation for the precision matrix to the indentity matrx (matrix of same dimensions of precision matrix filled with zeros and only ones across the diagonal). If you want to do things really correctly, specify a scaled inverse Wishart distribution as follows:

    vector[i,1] <- vector.raw[i,1] * xi1
    vector[i,2] <- vector.raw[i,2] * xi2
    vector.raw[i,1:2] ~ dmnorm(Mu.raw,Sigma2.inv.raw)
    # Priors for means of multivariate normal distribution
    Mu.raw[1] ~ dunif(0,100)
    Mu.raw[2] ~ dunif(0,100)
    # Scaled inverse Wishart prior
    Sigma2.inv.raw ~ dwish(R,scale) # R = 2x2 identity matrix, scale = 3
    xi1 ~ dunif(0,100) # scaling parameter for inverse wishart
    xi2 ~ dunif(0,100) # scaling parameter for inverse wishart
    # Rescale parameters to the scale of the data
    Sigma2.raw <- inverse(Sigma2.inv.raw)
    sigma1 <- pow(Sigma2.raw[1,1],0.5) * xi1 # rescaled standard deviation of vector[,1]
    sigma2 <- pow(Sigma2.raw[2,2],0.5) * xi2 # rescaled standard deviation of vector[,2]
    rho <- Sigma2.raw[1,2]/sqrt(Sigma2.raw[1,1]*Sigma2.raw[2,2]) # correlation between vector[,1] and vector[,2]
    mu[1] <- mu.raw[1] * xi1 # rescaled mean of vector[,1]
    mu[2] <- mu.raw[2] * xi2 # rescaled mean of vector[,2]

    Let me know if it works. Working with scaled inverse wishart priors is very prone to typing errors, in my experience!


    Last edit: Luc Coffeng 2013-03-06