A Bayesian EAP-Based Nonlinear Extension of Croon and Van Veldhoven’s Model for Analyzing Data from Micro–Macro Multilevel Designs

: Croon and van Veldhoven discussed a model for analyzing micro–macro multilevel designs in which a variable measured at the upper level is predicted by an explanatory variable that is measured at the lower level. Additionally, the authors proposed an approach for estimating this model. In their approach, estimation is carried out by running a regression analysis on Bayesian Expected a Posterior (EAP) estimates. In this article, we present an extension of this approach to interaction and quadratic effects of explanatory variables. Speciﬁcally, we deﬁne the Bayesian EAPs, discuss a way for estimating them, and we show how their estimates can be used to obtain the interaction and the quadratic effects. We present the results of a “proof of concept” via Monte Carlo simulation, which we conducted to validate our approach and to compare two resampling procedures for obtaining standard errors. Finally, we discuss limitations of our proposed extended Bayesian EAP-based approach.


Introduction
In organizational research, one may be interested in which factors determine an organization or team's outcome, such as the team's productivity in terms of the number of sales. For a general framework, see Schneider et al. [1]. The research conducted for the purpose of studying these factors typically face a data structure in which employees are nested in organizations/teams. To analyze such data, multilevel analysis is often used. The term "multilevel analysis" subsumes various approaches, some of which were developed in parallel in very different disciplines, such as economics, psychology, and education sciences. The multilevel analysis involves the use of hierarchical models, models with mixed effects, or multilevel models, terms that are largely used synonymously in the statistical literature. In the following, we will use the term "multilevel model". Over the past 20 years, these models have established themselves as the standard for analyzing multilevel data, mainly because of their usefulness to examine cross-level interactions and the availability of free software that has been developed and complemented commercial software for performing multilevel analyses. Moreover, a number of influential textbooks have been published (e.g., [2][3][4][5][6]) that make it easier to get started with multilevel analysis.
A research design can be considered a multilevel design if at least one of the variables has a multilevel structure. Following this view and focusing on one level of clustering (i.e., two-level designs), designs can then roughly be differentiated according to whether the dependent variable is measured at the lower or the upper level. Among the designs in which the dependent variable is measured at the upper level, the micro-macro design [6] stands out because it is of great relevance in organizational research. In this design, the dependent variable is predicted by a variable measured at the lower level (e.g., [7]; see also [8]). Croon and van Veldhoven [9] proposed an approach for analyzing data from micro-macro designs. An important feature of their approach is that estimation is carried out by an ordinary least squares regression analysis with Bayesian Expected a Posterior (EAP) estimates as input. Almost one decade after Croon and van Veldhoven [9] published their article, there has been renewed interest in the use of this type of estimation. Despite the criticism of its lack of flexibility [10], the method offers simple computations for an otherwise computationally very demanding class of models. For example, Zitzmann [11] emphasized the method's favorable statistical properties and showed how the method can be used to estimate multiple indicator two-level measurement models. Zitzmann and Helm [12] extended the method in such a way that it can be used for performing two-level mediation and moderation analyses. Aydin et al. [13] and Aydin et al. [14] discussed the necessity of an extension to three-level models, and Aydin and Algina [15] proposed such an extension.
In this article, we present an extension of Croon and van Veldhoven's approach for analyzing data from micro-macro designs to interaction and quadratic effects of explanatory variables. Specifically, (a) we define the Bayesian EAPs, (b) present an efficient way for estimating them, and (c) we explain how their estimates can be used to obtain the effects in the model. Moreover, (d) we conduct a "proof of concept" Monte Carlo study in which we also compare the performance of a standard bootstrap procedure for obtaining standard errors with the jackknife proposed by Zitzmann [11] jackknife. Finally, (e) we discuss one possibility to extend the model and a limitation of our proposed method.

Nonlinear Extension of Croon and Van Veldhoven's Model
We begin with the dependent variable, which is measured at the upper level. This may be a measure of productivity, such as a team's productivity Y. The explanatory variable, on the other hand, is measured at the lower level. For example, it is possible to assess a characteristic of the team leader such as leadership behavior by his or her employees' perceptions X 1 of the team leader's capacity to set a vision and goals and to support the team to achieve them. To this end, questionnaires or interviews are typically used [9]. One efficient way of assessing the leadership behavior of their team leader is ask employees to rate this characteristic and then average the ratings across the employees in the team [16]. This mean is assumed to reflect the employees' shared perception of their team leader and should ideally be formed by a latent instead of the common manifest aggregation procedure (see [17] for a detailed discussion). The latent mean ξ 1 varies between teams but not across employees in a team, whereas the individual deviations δ 1 from this mean vary across the employees within a team. Thus, the measurement model for a team leader's behavior in terms of his or her employees' shared perception can be formulated as: for an employee i = 1, . . . , n in a team j = 1, . . . , J. A central assumption is that employees from the same team differ less in their perceptions from each other than employees from different teams, which implies that the perceptions within a team are not independent of one another. This dependency can be quantified by the Intraclass Correlation (ICC). Formally, the ICC quantifies the proportion of the total variance that is located between the teams, and it can take on values between 0 and 1.
To investigate how leadership behavior ξ 1 is related to team productivity Y (controlling for a covariate Z), the following regression model can be used, which was proposed by Croon and van Veldhoven [9]: where β 1 denotes the effect of the explanatory variable of interest, which describes the effect of the leadership behavior on team productivity. γ is the coefficient of the covariate. The effect of the leadership behavior indicates by how many units the productivity of the team will increase on average when the leadership behavior is increased by one unit, while the covariate is kept constant. The covariate is a characteristic that is measured directly at the upper level. One may consider the personality of the team leader as a covariate, for example. The team leader's personality can be measured by asking the team leader rather than his or her employees. To extend the model, we consider an additional explanatory variable, which may be a team's climate. Similar to leadership behavior, this variable can be assessed by employees' perceptions X 2 . The latent mean ξ 2 across the employees in the team reflects the employees' shared perception of the team, and thus, the measurement model for a team's climate reads similar to the measurement model in Equation (1). Both explanatory variables may be assumed to interact in such a way that a positive team climate enlarges the effect of leadership behavior. To complete the nonlinear specification, it is further assumed that the explanatory variables have quadratic effects. The resulting model takes the form of the prototype for simultaneously estimating interaction and quadratic effects [18]. The model is expressed by the following regression. For better readability, we omit the employee and team indices: The coefficient β 3 denotes the interaction effect of the two explanatory variables, and β 4 and β 5 are their quadratic effects. The interaction effect indicates that the strength and perhaps even the direction of the relation between team productivity and leadership behavior is moderated by the team climate. In a similar vein, the quadratic effects of leadership behavior and team climate describe how the relations of the team productivity with leadership behavior and team climate depend on the levels of leadership behavior and team climate, respectively. Next, we give a detailed description of how this extended nonlinear model can be estimated with the help of Bayesian EAP estimates.

Bayesian EAP-Based Model Estimation
Technically speaking, the measurement model in Equation (1) can also be interpreted as decomposing each employee's perception into two orthogonal or uncorrelated parts: the latent mean across the perceptions of the employees of a team regarding their team leader's behavior and an individual deviation from this mean. Thus, the model allows the quantification of the differences in employees' perceptions between and within the teams by two variance components, which we hereafter call the between variance and the within variance.

Definition of (Adjusted) EAPs
To formally define a Bayesian EAP, we first make some distributional assumptions. Specifically, we assume the employees' perceptions X 1 to have a grand mean of zero (centering) to facilitate the presentation and later estimation. Moreover, we assume that the latent mean ξ 1 across the perceptions of the employees of the team follows a normal distribution: where ∼ reads "is distributed as", N(·, ·) stands for the Gaussian or (univariate) normal or distribution, and σ 2 ξ 1 is the between variance of the employee's perceptions. The individual deviations δ 1 from each latent mean are normally distributed around 0 with a within variance of σ 2 δ 1 . As the deviations' variance is the same for all teams, homoscedasticity is inherently assumed.
The latent mean across the perceptions of the employees of the jth team corresponds with the observed meanX 1•j = ∑ n i=1 X 1ij /n, where n is the number of employees surveyed in the team. Note that n does not vary across teams; that is, the same number of employees is surveyed in each team (balanced data). Thus, the team's likelihood function is: From a Bayesian perspective, Equation (4) can be viewed as a Bayesian prior distribution [19]. As mentioned in other articles of this Special Issue, Bayesian estimates can be obtained from the posterior distribution, which results from combining the prior with the likelihood (i.e., from Bayes' theorem; [20]). Now, if we combine the prior in Equation (4) with the likelihood in Equation (5), we yield the following (team-specific) normal posterior: where var(·) denotes the variance (thus, var X 1•j denotes the variance of the observed mean across the employees' perceptions). For detailed information about the derivation of this equation, see Appendix A. The mean of this posterior is called the EAP of the latent mean. ω is the weighting factor, which can be interpreted as the reliability of the observed counterpart of the latent mean across the perceptions of the employees of a team. The fact that the EAP is the mean of the posterior can be further substantiated by noticing that the setup is similar to the simple normal model in which Equation (4) is the prior for the variable's mean, and Equation (5) is the likelihood (see also [19] p. 134). Although we made explicit use of the Bayes' theorem in the definition of the EAP, we will skip the word "Bayesian" for the sake of simplicity. It is also interesting to note that from a non-Bayesian perspective, as indicated by Equation (7), the EAP is simply the prediction from a regression of the latent on the observed mean [21], and this is why the weight can be analytically derived via the ordinary least squares principle.
In models such as Croon and van Veldhoven's model and our extended nonlinear model, effects can be biased when EAP estimates are used as input in ordinary least squares regression analyses [22]. Therefore, it is necessary to adjust EAPs for other variables in the model. The model in Equation (2) contains also the team leader's personality Z. Thus, an adjustment needs to be made to the EAP of the latent mean across the perceptions of the employees of a team. The adjusted EAP can be expressed as follows: where β is the result of a regression ofX 1•j on Z. In this equation, ω is the conditional reliability (i.e., conditioned on the covariate; [23]) (see Appendix B for detailed information about the derivation). Because the nonlinear model in Equation (3) contains the team climate as another explanatory variable, this variable also needs to be taken into account. Thus, one way to express the adjusted EAP is: whereX 2•j = ∑ n i=1 X 2ij /n is the observed mean across the perceptions of the employees of a team regarding the team's climate, σ ξ 1 ξ 2 is the between covariance of the employee's perceptions X 1 and X 2 , and σ ξ 1 Z is the between covariance of X 1 and Z. var −1 (·) gives the inverted covariance matrix, T is the transpose operator, and ω is a vector containing the weights ofX 1•j ,X 2•j , and Z j . The equation indicates once more that the EAP is the prediction from a regression and that the weights are the coefficients in this regression and can thus be derived by ordinary least squares. To adjust the EAP of the latent mean across the perceptions of the employees of a team regarding the team's climate, an analogous adjustment is made.
What is more interesting is how an adjusted EAP can be defined for the product of the latent means in Equation (3). To tackle this problem, it is instructive to consider one of its solutions in single-level models. Here, products of indicators from two explanatory variables are formed and used as derived "product indicators" of a latent variable whose effect is the interaction effect (i.e., the product indicator approach; e.g., [24][25][26]). Moreover, scholars have pointed out that the measurement models in our extended nonlinear model (see Equation (1)) can be viewed as describing the relation between a latent variable (i.e., the latent mean) and its indicators (i.e., the perceptions of the employees of a team; [27]). This similarity suggests that the EAP of the product of the two latent means can be defined as the prediction from a regression on the observed mean across the products of the employees' perceptions X 1 of their team leaders and X 2 of their team, which is given as (3) includes also the quadratic effects, an adjustment for the squares of the latent means needs to be made. The adjusted EAP is: where X 2 1 •j = ∑ n i=1 X 2 1ij /n and X 2 2 •j = ∑ n i=1 X 2 2ij /n. In this equation, the variance of the product of the latent means and the covariances between this product and the squares of the latent means are expressed in terms of the variances and covariances of the latent means applying Theorem 13 of Bohrnstedt and Goldberger [28]. It is important to note that although the nonlinear model contains the latent means and the covariate as variables, these variables need not to be taken into account in the adjustment. This is because the product of the latent means and their squares are uncorrelated with these other variables due to centering. The adjusted EAP of the square of the latent mean across the perceptions of the employees of a team regarding their team leader's behavior is: The adjusted EAP of the square of the latent mean across the perceptions of the employees of a team regarding their team's climate can analogously be expressed. So far, we have only defined the EAPs. In the next section, we will discuss how these EAPs can efficiently be estimated in order to help applied researchers implement the extended nonlinear model on their own.

Estimating the Adjusted EAPs
To estimate the adjusted EAPs, we first need to estimate the variances and covariances in the equations. Estimates will be indicated by a hat (ˆ) symbol. The estimates of the variances and covariances can then be plugged into the equations in order to obtain the EAP estimates. A straightforward way of estimating the variances and covariances is using estimators from the Analysis of Variance (ANOVA) literature (see [9,15]). We begin with the estimation of the variances and covariances of the observed quantitiesX 1•j ,X 2•j , and Z j . In the case of balanced data (i.e., equal numbers of employees across teams), the estimate of the variance of the observed mean across the perceptions of the employees of a team regarding their team leader's behavior is given as the mean sum of squares: /n is the grand mean across the employees' perceptions. Note that Equation (12) is simply the sample variance ofX 1•j . For the variance of the observed mean across the perceptions of the employees of a team regarding the team's climate, an analogous estimate is formulated. To estimate the variance of the personality Z of the team leader, one can use: The estimate of the covariance betweenX 1•j and Z is the mean sum of cross products: An analogous estimate is used for the covariance ofX 2•j and Z. The covariance betweenX 1•j andX 2•j is estimated by: One efficient way for computing all estimates at once is computing the estimate of a covariance matrix: In the estimation of the adjusted EAP of the latent mean across the perceptions of the employees of a team, this estimate can be used in place of the covariance matrix of the observed quantities in Equation (9). However, to estimate the adjusted EAP, one also needs to estimate the variance of the latent mean and its covariances. Estimating these between variances and covariances is less straightforward but still not very difficult. What helps to find estimates is that scholars have noticed that the measurement model of a latent mean is a random-effects model [27], with the latent mean being the team-specific random effect. Due to this equivalence between a measurement model and random-effects model, it is possible to employ the ANOVA estimator for the variance of the random effect to estimate the between variance (see, e.g., [29] for illustrations of how these estimates can be drived using the ANOVA method). The estimate of the between variance in the measurement model for the latent mean across the perceptions of the employees of a team regarding their team leader's behavior is given by: As the measurement model for the latent mean across the perceptions of the employees of a team regarding the team's climate is structurally equivalent, the estimate of the between variance σ 2 ξ 2 looks similar. The estimates of the between covariances are: If we plug these estimates of the between variance and covariances together with the estimate of covariance matrix into Equation (9), then the estimate of the adjusted EAP is: From a computational perspective, it makes sense to compute this estimate and the estimate of the adjusted EAP of the latent mean across the perceptions of the employees of a team regarding the team's climate in one step in order to avoid the need for a second matrix inversion and thus to save computing time. To this end, we compute: Note thatΩ is the estimate of a weight matrix rather than of a weight vector in this equation. Alternatively, the result of the matrix inversion could first be stored, and then, both EAP estimates could be computed separately using this stored result.
To estimate the adjusted EAPs of the product of the two latent means and their squares, we only need to find estimates of the variances and covariances of the observed counterparts in the equations for these EAPs. The estimation equations for the between variances and covariances have already been discussed. The efficient method for computing estimates of the observed variances and covariances is computing the following covariance matrix: Together with the estimates of the between variances and covariances, this covariance matrix estimate can be plugged into the equations in order to estimate the adjusted EAPs of the product of the two latent means and their squares. As mentioned, to avoid multiple matrix inversions, it makes sense to estimate these EAPs simultaneously by: Once all adjusted EAP estimates are computed, these estimates can be used to obtain the coefficients in our extended nonlinear model.

Using Adjusted EAP Estimates to Estimate the Model
The nonlinear model in Equation (3) is estimated by running an ordinary least squares regression analysis on the estimates of the adjusted EAPs of the latent means, their product, their squares, and the covariate. To this end, the following regression model is specified and estimated: where the estimates of the coefficients are the usual ordinary least squares estimates, which are assumed to be unbiased for the corresponding coefficients in Equation (3). In other words, when certain conditions are met (e.g., a sufficiently large sample to overcome smallsample bias), the regression analysis will, on average, yield the actual coefficients-an assumption that we will test empirically in the next section of this article by using computer simulations. Apart from this favorable feature, there is also a significant drawback. It has been emphasized that regression analysis can provide incorrect standard errors, particularly when the data are very unbalanced or homoscedasticity is violated [9]. Some suggestions have been made in order to overcome this limitation. For example, Croon and van Veldhoven [9] suggested the method developed by Davidson and MacKinnon [30] be used in this case. Another suggestion is using resampling procedures such as Zitzmann's jackknife. For example, for the interaction effect β 3 in the model, the procedure first computes estimates of this effect from R subsamples, each omitting d teams. d is typically much smaller than the number J of teams. In this specific variant of the jackknife, the subsamples are obtained by dividing the indices (1, . . . , J) into R = J/d non-overlapping subsets and then using these subsets to create the subsamples. The standard error is then estimated on the basis of the effect estimates by: j=1β 3j /R is the mean across the estimates. This jackknife is computationally very efficient because it can perform well with only 20 subsamples (see [11]). However, some possible problems have been encountered with the use of this procedure (e.g., a standard error that is too large due to a few extreme estimates), and this is why a bootstrap (e.g., [31]) could be an alternative. The bootstrap repeatedly estimates the effect from R * subsamples, which are random samples from the original sample (drawn with replacement). As the size of these samples, one could choose J − d (i.e., the subsample size in the jackknife). Using the estimates, the standard error is: R * is typically a large number, and Preacher and Hayes [32] recommended at least 1000 subsamples, for example. In the next section, we will present the results of a simulation study, which we conducted to validate our extension of the EAP-based approach to interaction and quadratic effects and to compare the two resampling procedures for obtaining standard errors.

Monte Carlo Study
We conducted a Monte Carlo simulation study to find out whether the EAP-based approach performs as expected, meaning that whether it yields unbiased results for the coefficients in the nonlinear model when the sample of teams is sufficiently large. To this end, we assumed the dependent variable Y (i.e., the team's productivity), the covariate Z (i.e., the team leader's self-assessed personality), and the two variables X 1 and X 2 (i.e., the employees' perceptions of their team leader and their team, respectively) to be standardnormally distributed with a variance of 1. All variables except Y were grand-mean centered, and X 1 and X 2 both had a multilevel structure (i.e., the employees' perceptions split into latent team means (i.e., shared perceptions) and individual deviations from these shared perceptions). The between correlation of X 1 and X 2 (i.e., the correlation between ξ 1 and ξ 2 ) was set to 0.3 as were the between correlations between the two X variables and Z. We selected this value because 0.3 lies in the middle of the typical range of empirically observed correlations.
We set the number of teams to 60, 100, or 200 teams. Whereas 60 and 100 teams served to study the behavior of our method in small-sample, 200 teams served as the test case for the absence of bias in large samples and thus for the validity of our approach (asymptotic unbiasedness). The number of employees surveyed in a group was either 10 or 30 employees, and the ICCs of X 1 and X 2 were both set to 0.1 or 0.3, implying reliabilities of the observed meansX 1•j andX 2•j between 0.53 and 0.93. For each of the 3 × 2 × 2 = 12 data constellations, we generated 1000 data sets from the following nonlinear model: and we estimated this model as described above in the popular statistical computing environment R, which is also known as "the R Project for Statistical Computing" [33]. We studied the bias of our approach in recovering the main, interaction, and quadratic effects of the explanatory variables and the coefficient of the covariate in the model. The bias is defined as the deviation of the expected value of an estimate from the true value in the data-generating model. We divided the bias by the true value (e.g., [34]) and multiplied it by 100% in order to obtain the percentage of bias-a measure that helped to judge whether the bias was negligible or not. An absolute value of less than 10% was considered negligible [35]. In addition, we studied the performance of the bootstrap for obtaining the standard errors of the effects with 1000 subsamples. As the criterion for the evaluation of the bootstrap, we used the coverage rate, which is the probability that the 95% Confidence Interval (CI), which is, for example,β 3 − z 1−α/2σβ 3 ≤ x ≤β 3 + z 1−α/2σβ 3 for the interaction effect in the model, captures the true value. A coverage rate of less than 91% (more than 98%) was considered to indicate that the standard errors were underestimated (overestimated; [35]). In addition, we compared the bootstrap with Zitzmann's jackknife. The CI from the jackknife differs only in that Equation (25) is used in place of the bootstrapped standard error. Table 1 shows the percentages of bias in estimating the coefficients in the model. As the biases were similar for both main effects and both quadratic effects, the table shows only the biases for the first main effect and the first quadratic effect. In addition, it shows the coverage rates for the bootstrap and Zitzmann's jackknife. Despite the match between the nonlinear data-generating model and the analysis model, our approach for estimating this model provided substantially biased estimates of the main, interaction, and quadratic effects (i.e., it exhibited biases of more than 10%) but not for the coefficient of the covariate. The biases were particularly pronounced when the number of teams was rather small and the ICCs of the two Xs was low. However, these large small-sample biases could have been expected. What is more important for the validation of the approach is that the approach provided unbiased estimates in large samples, which suggests that the approach is indeed an asymptotically unbiased one.
A similar picture emerged with regard to the coverage rates. Some of the coverage rates were too high (i.e., >98%) when the ICCs were low, particularly those of the main, interaction, and quadratic effects. It is interesting to note that the coverage rates for the bootstrap were slightly less accurate than those for Zitzmann's jackknife. However, they both tended to become close to the nominal level of 95% when the ICCs increased, indicating that the procedures provide correct standard errors.

Discussion
Multilevel data with employees nested in organizations/teams are often collected to examine which factors determine an organization/team's performance. Of great importance in organizational research are micro-macro designs, in which an upper level variable is predicted by a variable measured at the lower level. To analyze the data from such designs, Croon and van Veldhoven [9] proposed an approach in which an ordinary least squares regression analysis is carried out with Bayesian EAP estimates (which we also called EAPs for the sake of simplicity) as input. In this article, we developed and showed an extension of this approach to interaction and quadratic effects of explanatory variables measured at the lower level. After defining the EAPs, we discussed how these EAPs can be estimated and used to obtain interaction and quadratic effects. We conducted a simulation study to validate our extended approach, and we compared two procedures for obtaining standard errors with each other: a standard bootstrap procedure and Zitzmann's jackknife. The main findings were the following. First and foremost, the results with regard to the bias showed no asymptotic bias, indicating the validity of our approach. Moreover, the bootstrap and the jackknife performed very similarly overall, with a small disadvantage of the bootstrap under challenging conditions. However, because the bootstrap uses much more subsamples than the jackknife (1000 vs. only 20 subsamples), the different procedures differ in terms of the computing effort. As often, a series of different variants of a model or even different variants of different models are estimated, computing effort adds up, which speaks for the use of the jackknife in research practice.
Regarding possible future extensions, it would be interesting to add further dependent variables to the model. However, such a multivariate model is not estimable with ordinary least squares regression analysis of EAP estimates. Thus, multiple univariate models need to be estimated (one per dependent variable), which highlights once more that our approach is a stepwise and thus only a limited information method [9]. However, the fact that it does not estimate the model at once but divides it into simpler submodels can also be seen as a feature ( [36,37]; see also [38]). Studying the potential of this type of estimation for estimating complex models is an interesting subject for future research.
A limitation of the use of EAP estimates should nevertheless be mentioned. The approach places relatively high demands on the data. Specifically, the number of teams should be rather large, and the ICCs should not be too low. It is interesting to note that alternatively, the analysis can be carried out with Structural Equation Modeling (SEM) software. However, the demands tend to be even higher with this software, at least when maximum likelihood methods are used [12], which are the default in commercial software such as Mplus [39]. Prominent examples of maximum likelihood methods for nonlinear models are latent moderated structural equations (LMS; [40]) and quasimaximum likelihood (QML; [41]). These methods determine the estimates in such a way that the probability of the data is maximized under the model. They give unbiased results and small standard errors only when the data provide a large amount of information (i.e., many teams, high ICCs). One possible alternative to our approach and the maximum likelihood methods is the Bayesian Markov Chain Monte Carlo (MCMC) method, which iteratively samples from conditional distributions, thereby creating an MCMC chain from which the estimates can be computed. Similar to our approach, this method places less demands on the data than maximum likelihood ( [42,43]; see also [44,45]). Furthermore, it enables more flexibility in specifying models [46], and it leads less frequently to estimation problems [47,48]. However, the Bayesian MCMC method can be problematic despite the mentioned advantages. For example, this method comes at the cost of long computing times because a single iteration of the algorithm is slow, and the method requires many iterations to reach convergence (see [49] for a discussion; see also [50]), particularly when the model is empirically hardly identified. See Hecht et al. [51], Merkle et al. [52], Xu and Liao [53] and Yi and Tang [54] for ways to speed up the method. See Hecht et al. [55] and Hecht and Zitzmann [56] for applications to the analysis of longitudinal data. Whether and under which conditions the Bayesian MCMC method outperforms our approach is however an open question, which can best be addressed in extensive simulation work.
To conclude, we showed that Croon and van Veldhoven's EAP-based approach for analyzing micro-macro multilevel designs can be extended to interaction and quadratic effects, and we hope that the article will contribute to the use of this approach in organizational research. The approach can easily be implemented in any statistics software, such as SPSS, SAS or R. Finally, we would like to stress that the application of this approach is not limited to the organizational context and that it may also be an interesting option in other areas of multilevel research in which latent means are used (e.g., education sciences, psychology; [14,57]).

Conflicts of Interest:
The authors declare no conflict of interest. (6) To derive Equation (6), consider the simple normal model first:

Appendix A. Derivation of Equation
where µ is the variable's mean, and ε i (i = 1, . . . , N) are normally distributed residuals with variance σ 2 . If we select the normal prior µ ∼ N 0, σ 2 µ (A2) and combine this prior with the likelihood X • ∼ N µ, σ 2 /N (A3) whereX • = ∑ N i=1 X i /N is the observed counterpart of µ, we obtain the following posterior: The mean of this distribution (i.e., the EAP of µ) is Now, consider the above mentioned measurement model (Equation (1)) for a team leader's behavior in terms of his or her employees' shared perception: (see also Equation (1) in the main body of the text). We choose the prior and combine this prior with the the team's likelihood in order to obtain the posterior. Because the setup (i.e., normal prior, likelihood) is similar to the simple normal model, the posterior has to look similar; that is, it has to be of the same form. Hence, the (team-specific) posterior reads: Appendix B. Derivation of Equation (8) To derive Equation (8), recall that when the model contains an additional covariate, the adjusted EAP is the prediction from a regression of the latent mean on the observed mean and the covariate. More formally, this means: In order to obtain ω, we first compute the inverted covariance matrix: By using Equation (A11) and the term β = σ ξ 1 Z var(Z j ) , we yield ω as: +var(X 1•j )σξ 1 Z var(X 1•j )var(Zj)−σ 2 −β 2 var(Z j ) var(X 1•j )−β 2 var(Z j ) var(X 1•j )−σ 2 ξ 1 −β 2 var(Z j )+β 2 var(Z j ) var(X 1•j )−β 2 var(Z j ) β Finally, substituting Equation (A12) into Equation (A10) yields: