mmVarFit: Fit Mixed Membership models using variational EM
In mixedMem: Tools for Discrete Multivariate Mixed Membership Models

Description Usage Arguments Details Value References See Also Examples

mmVarFit is the primary function of the mixedMem package. The function fits parameters φ and δ for the variational distribution of latent variables as well as psuedo-MLE estimates for the population parameters α and θ. See documentation for mixedMemModel or the package vignette, "Fitting Mixed Membership Models Using mixedMem" for a more detailed description of the variables and notation in a mixed membership model.

mmVarFit(model, printStatus = 1, printMod = 1, stepType = 3,
  maxTotalIter = 500, maxEIter = 1000, maxAlphaIter = 200,
  maxThetaIter = 1000, maxLSIter = 400, elboTol = 1e-06,
  alphaTol = 1e-06, thetaTol = 1e-10, aNaught = 1, tau = 0.899,
  bMax = 3, bNaught = 1000, bMult = 1000, vCutoff = 13,
  holdConst = c(-1))

`model`	a `mixedMemModel` object created by the `mixedMemModel` constructor.
`printStatus`	an integer 0 or 1. When `printStatus` is 1 `mmVarFit` will print status updates, when `printStatus` is 0 `mmVarFit` will not print status updates.
`printMod`	a positive integer which specifies how often to print status updates. The status will be printed at each step which is a multiple of `printMod`. The default value is 1.
`stepType`	an integer from 0-3 which specifies what steps to carry out in estimation. 0 performs a single E-Step; this can be used to find the held out ELBO. 1 performs iterates between E-steps and M-steps for estimating theta but keeps alpha constant. 2 iterates between E-steps and M-steps for estimating alpha, but keeps theta constant. 3 completes full E and M steps. The default is 3.
`maxTotalIter`	the maximum total steps before termination. A full E and M step together count as 1 step. If this maximum is ever achieved, a warning message will be printed at convergence. The default is 500.
`maxEIter`	the maximum iterations for each E-Step. If this maximum is ever achieved, a warning message will be printed at convergence. The default is 1000.
`maxAlphaIter`	the maximum iterations when fitting alpha. If this maximum is ever achieved, a warning message will be printed at convergence. The default is 200.
`maxThetaIter`	the maximum iterations when fitting theta. If this maximum is ever achieved, a warning message will be printed at convergence. The default is 1000.
`maxLSIter`	the maximum backtracking iterations in the line search for updating alpha and theta for rank data.
`elboTol`	the convergence criteria for the EM Algorithim. When the relative increase in the ELBO is less than the convergence critiera, the algorithm converges. The default is 400.
`alphaTol`	the convergence criteria for updates to alpha. When the relative increase in the ELBO is less than the convergence critiera, the update to alpha converges. The defualt is 1e-6.
`thetaTol`	the convergence criteria for updates to theta. When the relative increase in the ELBO is less than the convergence critiera, the update to theta converges. The default is 1e-10.
`aNaught`	the first step size in the backtracking line search used to update alpha and theta for rank data. The default is 1.
`tau`	the backtracking factor in the backtracking line search used to update alpha and theta for rank data. The default is .899.
`bMax`	the number of times the penalization factor is increased in the interior point method for fitting theta for rank data. The default is 3.
`bNaught`	the initial penalization factor in the interior point method for fitting theta for rank data. `bNaught` must be positive and the default value is 1000.
`bMult`	a positive number by which bNaught is multiplied by in each iteration of the interior point method for fitting theta for rank data. The default is 1000.
`vCutoff`	a positive integer cutoff for Vj at which a gradient ascent method is used instead of the Newton Raphson interior point method. This is used to avoid inverting a large matrix. The default is 13.
`holdConst`	a vector of integers specifying groups to be held fixed during the estimation procedure. The estimation algorithim will hold the theta parameters of these specific groups constant, but update all other parameters. The group numbers range from 0 to K - 1. To estimate all groups, use the default value of c(-1).

mmVarFit selects psuedo-MLE estimates for α and θ and approximates the posterior distribution for the latent variables through a mean field variational approach. The variational lower bound on the log-likelihood at the data given model parameters is given by Jensen's inequality:

E_Q{\log[p(X,Z, Λ)]} - E_Q{\log[Q(Z, Λ|φ, δ)]} ≤ \log P(obs |α, θ) where

Q = ∏_i [Dirichlet(λ_i|φ_i) ∏_j^J ∏_r^{R_j} ∏_n^{N_{i,j,r}}Multinomial(Z_{i,j,r,n}|δ_{i,j,r,n})]

The left hand side of the above equation is often referred to as the Evidence Lower Bound (ELBO). The global parameters α and θ are selected to maximize this bound as a psuedo-MLE procedure. In addition, it can be shown that maximizing the ELBO with respect to the variational parameters φ and δ is equivalent to minimizing the KL divergence between the tractable variational distribution and the true posterior.

The method uses a variational EM approach. The E-step considers the global parameters α and θ fixed and picks appropriate variational parameters φ and δ to minimize the KL divergence between the true posterior and the variational distribution (or maximize the ELBO). The M-step considers the variational parameters fixed and global parameters θ and α are selected which maximize the lower bound on the log-likelihood (ELBO).

mmVarFit returns a mixedMemModel object containing updated variational parameters and global parameters.

Beal, Matthew James. Variational algorithms for approximate Bayesian inference. Diss. University of London, 2003.

mixedMemModel

## The following example generates multivariate observations (both variables are multinomial)
## from a mixed membership model. The observed data is then used to fit a
## mixed membership model using mmVarFit

## Generate Data
Total <- 30 #30 Individuals
J <- 2 # 2 variables
dist <- rep("multinomial",J) # both variables are multinomial
Rj <- rep(100,J) # 100 repititions for each variable
# Nijr will always be 1 for multinomials and bernoulli's
Nijr <- array(1, dim = c(Total, J, max(Rj)))
K <- 4 # 4 sub-populations
alpha <- rep(.5, K) # hyperparameter for dirichlet distribution
Vj <- rep(5, J) # each multinomial has 5 categories
theta <- array(0, dim = c(J, K, max(Vj)))
theta[1,,] <- gtools::rdirichlet(K, rep(.3, 5))
theta[2,,] <- gtools::rdirichlet(K, rep(.3, 5))
obs <- rmixedMem(Total, J, Rj, Nijr, K, Vj, dist, theta, alpha)$obs

## Initialize a mixedMemModel object
test_model <- mixedMemModel(Total = Total, J = J,Rj = Rj, Nijr= Nijr,
 K = K, Vj = Vj,dist = dist, obs = obs,
  alpha = alpha, theta = theta+0)

## Fit the mixed membership model
out <-mmVarFit(test_model)