In this article, we develop a class of two-stage models to accommodate three common characteristics of behavioral data. First, behavior is invariably multivariate in its conceptualization and communication. Separate univariate analyses of related outcome variables are fraught with potential interpretive blind spots for the researcher. This practice also suffers, from an inferential standpoint, because it fails to take advantage of any redundant information in the outcomes. Second, studies of behavior, especially in experimental research, employ smaller samples. This situation raises issues of robustness of inference with respect to outlying individuals. Third, the outcome variable may have observations missing because of accidents or by design. The model permits the estimation of the full spectrum of plausible measurement error structures while using all the available information. Maximum likelihood estimates are obtained for various members of a multivariate hierarchical linear model (MHLM), and, in the context of several illustrative examples, these estimates match closely the results from a Bayesian approach to the normal-normal MHLM and to the normal- t MHLM.