(Restricted) Maximum Likelihood Estimation with Prior Distributions and Penalty Functions under Multivariate Normality

The mlnormal estimates statistical model for multivariate normally distributed outcomes with specified mean structure and covariance structure (see Details and Examples). Model classes include multilevel models, factor analysis, structural equation models, multilevel structural equation models, social relations model and perhaps more.

The estimation can be conducted under maximum likelihood, restricted maximum likelihood and maximum posterior estimation with prior distribution. Regularization (i.e. LASSO penalties) is also accomodated.

mlnormal(y, X, id, Z_list, Z_index, beta=NULL, theta, method="ML", prior=NULL,
    lambda_beta=NULL, weights_beta=NULL, lambda_theta=NULL, weights_theta=NULL,
    beta_lower=NULL, beta_upper=NULL,    theta_lower=NULL, theta_upper=NULL,
    maxit=800, globconv=1e-05, conv=1e-06, verbose=TRUE, REML_shortcut=NULL,
    use_ginverse=FALSE, vcov=TRUE, variance_shortcut=TRUE, use_Rcpp=TRUE,
    level=0.95, numdiff.parm=1e-04, control_beta=NULL, control_theta=NULL)

# S3 method for mlnormal
summary(object, digits=4, file=NULL, ...)

# S3 method for mlnormal
print(x, digits=4, ...)

# S3 method for mlnormal
coef(object, ...)

# S3 method for mlnormal
logLik(object, ...)

# S3 method for mlnormal
vcov(object, ...)

# S3 method for mlnormal
confint(object, parm, level=.95, ... )

Arguments

y: Vector of outcomes
X: Matrix of covariates
id: Vector of identifiers (subjects or clusters, see Details)
Z_list: List of design matrices for covariance matrix (see Details)
Z_index: Array containing loadings of design matrices (see Details). The dimensions are units $\times$ matrices $\times$ parameters.
beta: Initial vector for $\bold{\beta}$
theta: Initial vector for $\bold{\theta}$
method: Estimation method. Can be either "ML" or "REML".
prior: Prior distributions. Can be conveniently specified in a string which is processed by the function prior_model_parse. Only univariate prior distributions can be specified.
lambda_beta: Parameter $\lambda_{\bold{\beta}}$ for penalty function $P( \bold{\beta} )=\lambda_{\bold{\beta}} \sum_h w_{\bold{\beta}h} | \beta _h |$
weights_beta: Parameter vector $\bold{w}_{\bold{\beta}}$ for penalty function $P( \bold{\beta} )=\lambda_{\bold{\beta}} \sum_h w_{\bold{\beta}h} | \beta _h |$
lambda_theta: Parameter $\lambda_{\bold{\theta}}$ for penalty function $P( \bold{\theta} )=\lambda_{\bold{\theta}} \sum_h w_{\bold{\theta}h} | \theta _h | $
weights_theta: Parameter vector $\bold{w}_{\bold{\theta}}$ for penalty function $P( \bold{\theta} )=\lambda_{\bold{\theta}} \sum_h w_{\bold{\theta}h} | \theta _h | $
beta_lower: Vector containing lower bounds for $\bold{\beta}$ parameter
beta_upper: Vector containing upper bounds for $\bold{\beta}$ parameter
theta_lower: Vector containing lower bounds for $\bold{\theta}$ parameter
theta_upper: Vector containing upper bounds for $\bold{\theta}$ parameter
maxit: Maximum number of iterations
globconv: Convergence criterion deviance
conv: Maximum parameter change
verbose: Print progress?
REML_shortcut: Logical indicating whether computational shortcuts should be used for REML estimation
use_ginverse: Logical indicating whether a generalized inverse should be used
vcov: Logical indicating whether a covariance matrix of $\bold{\theta}$ parameter estimates should be computed in case of REML (which is computationally demanding)
variance_shortcut: Logical indicating whether computational shortcuts for calculating covariance matrices should be used
use_Rcpp: Logical indicating whether the Rcpp package should be used
level: Confidence level
numdiff.parm: Numerical differentiation parameter
control_beta: List with control arguments for $\bold{\beta}$ estimation. The default is
list( maxiter=10, conv=1E-4, ridge=1E-6).
control_theta: List with control arguments for $\bold{\theta}$ estimation. The default is
list( maxiter=10, conv=1E-4, ridge=1E-6).
object: Object of class mlnormal
digits: Number of digits used for rounding
file: File name
parm: Parameter to be selected for confint method
...: Further arguments to be passed
x: Object of class mlnormal

Details

The data consists of outcomes $\bold{y}_i$ and covariates $\bold{X}_i$ for unit $i$. The unit can be subjects, clusters (like schools) or the full outcome vector. It is assumed that $\bold{y}_i$ is normally distributed as $N( \bold{\mu}_i, \bold{V}_i )$ where the mean structure is modelled as $$ \bold{\mu}_i=\bold{X}_i \bold{\beta} $$ and the covariance structure $ \bold{V}_i$ depends on a parameter vector $\bold{\theta}$. More specifically, the covariance matrix $ \bold{V}_i$ is modelled as a sum of functions of the parameter $\bold{\theta}$ and known design matrices $\bold{Z}_{im}$ for unit $i$ ($m=1,\ldots,M$). The model is $$\bold{V}_i=\sum_{m=1}^M \bold{Z}_{im} \gamma_{im} \qquad \mathrm{with} \qquad \gamma_{im}=\prod_{h=1}^H \theta_h^{q_{imh}} $$ where $q_{imh}$ are non-negative known integers specified in Z_index and $\bold{Z}_{im}$ are design matrices specified in Z_list.

The estimation follows Fisher scoring (Jiang, 2007; for applications see also Longford, 1987; Lee, 1990; Gill & Swartz, 2001) and the regularization approach is as described in Lin, Pang and Jiang (2013) (see also Krishnapuram, Carin, Figueiredo, & Hartemink, 2005).

Value

List with entries

theta: Estimated $\bold{\theta}$ parameter
beta: Estimated $\bold{\beta}$ parameter
theta_summary: Summary of $\bold{\theta}$ parameters
beta_summary: Summary of $\bold{\beta}$ parameters
coef: Estimated parameters
vcov: Covariance matrix of estimated parameters
ic: Information criteria
V_list: List with fitted covariance matrices $\bold{V}_i$
V1_list: List with inverses of fitted covariance matrices $\bold{V}_i$
prior_args: Some arguments in case of prior distributions
...: More values

References

Gill, P. S., & Swartz, T. B. (2001). Statistical analyses for round robin interaction data. Canadian Journal of Statistics, 29, 321-331. doi:10.2307/3316080

Jiang, J. (2007). Linear and generalized linear mixed models and their applications. New York: Springer.

Krishnapuram, B., Carin, L., Figueiredo, M. A., & Hartemink, A. J. (2005). Sparse multinomial logistic regression: Fast algorithms and generalization bounds. IEEE Transactions on Pattern Analysis and Machine Intelligence, 27, 957-968. doi:10.1109/TPAMI.2005.127

Lee, S. Y. (1990). Multilevel analysis of structural equation models. Biometrika, 77, 763-772. doi:10.1093/biomet/77.4.763

Lin, B., Pang, Z., & Jiang, J. (2013). Fixed and random effects selection by REML and pathwise coordinate optimization. Journal of Computational and Graphical Statistics, 22, 341-355. doi:10.1080/10618600.2012.681219

Longford, N. T. (1987). A fast scoring algorithm for maximum likelihood estimation in unbalanced mixed models with nested random effects. Biometrika, 74, 817-827. doi:10.1093/biomet/74.4.817

Examples