Uncollapse output from optimPibbleCollapsed to full pibble Model

Source: R/RcppExports.R

See details for model. Should likely be called following optimPibbleCollapsed. Notation: N is number of samples, D is number of multinomial categories, Q is number of covariates, iter is the number of samples of eta (e.g., the parameter n_samples in the function optimPibbleCollapsed)

uncollapsePibble(
  eta,
  X,
  Theta,
  Gamma,
  Xi,
  upsilon,
  seed,
  ret_mean = FALSE,
  ncores = -1L
)

Arguments

eta

array of dimension (D-1) x N x iter (e.g., Pars output of function optimPibbleCollapsed)

X

matrix of covariates of dimension Q x N

Theta

matrix of prior mean of dimension (D-1) x Q

Gamma

covariance matrix of dimension Q x Q

Xi

covariance matrix of dimension (D-1) x (D-1)

upsilon

scalar (must be > D) degrees of freedom for InvWishart prior

seed

seed to use for random number generation

ret_mean

if true then uses posterior mean of Lambda and Sigma corresponding to each sample of eta rather than sampling from posterior of Lambda and Sigma (useful if Laplace approximation is not used (or fails) in optimPibbleCollapsed)

ncores

(default:-1) number of cores to use, if ncores==-1 then uses default from OpenMP typically to use all available cores.

Value

List with components

  1. Lambda Array of dimension (D-1) x Q x iter (posterior samples)

  2. Sigma Array of dimension (D-1) x (D-1) x iter (posterior samples)

  3. The number of cores used

  4. Timer

Details

Notation: Let Z_j denote the J-th row of a matrix Z. While the collapsed model is given by: $$Y_j ~ Multinomial(Pi_j)$$ $$Pi_j = Phi^{-1}(Eta_j)$$ $$Eta ~ T_{D-1, N}(upsilon, Theta*X, K, A)$$ Where A = I_N + X * Gamma * X', K = Xi is a (D-1)x(D-1) covariance matrix, Gamma is a Q x Q covariance matrix, and Phi^-1 is ALRInv_D transform.

The uncollapsed model (Full pibble model) is given by: $$Y_j ~ Multinomial(Pi_j)$$ $$Pi_j = Phi^{-1}(Eta_j)$$ $$Eta ~ MN_{D-1 x N}(Lambda*X, Sigma, I_N)$$ $$Lambda ~ MN_{D-1 x Q}(Theta, Sigma, Gamma)$$ $$Sigma ~ InvWish(upsilon, Xi)$$ This function provides a means of sampling from the posterior distribution of Lambda and Sigma given posterior samples of Eta from the collapsed model.

References

JD Silverman K Roche, ZC Holmes, LA David, S Mukherjee. Bayesian Multinomial Logistic Normal Models through Marginally Latent Matrix-T Processes. 2019, arXiv e-prints, arXiv:1903.11695

Examples

sim <- pibble_sim()

# Fit model for eta
fit <- optimPibbleCollapsed(sim$Y, sim$upsilon, sim$Theta%*%sim$X, sim$KInv, 
                             sim$AInv, random_pibble_init(sim$Y))  

# Finally obtain samples from Lambda and Sigma
fit2 <- uncollapsePibble(fit$Samples, sim$X, sim$Theta, 
                                   sim$Gamma, sim$Xi, sim$upsilon, 
                                   seed=2849)