I have repeated observations (count data) conducted in different plots (say, 5 times each plot). I would like to conduct a regression analysis of the observation data at the plot level with plot-level independent variables. I'd like to consider approaches besides putting the plot-level independent variables into an observation-level model because I plan to build a structural equation model that will relate to other datasets at the plot level.

I could average the repeated observations for each plot, but each observation has individual sources of variability (e.g. weather conditions, time of day, surveyor ID, etc.) that could be systematically affecting the data, so I would like to factor them out.

Would it be possible to remove these nuisance sources of variability by either:

  1. Making an observation-level mixed-effects model with plot as a random effect grouping and running the plot-level analysis on the fitted random effects coefficients.


  1. Fitting an observation-level model without mixed effects and averaging the residuals of this model by plot groupings and running the plot-level analysis on the averaged residuals.

Are either of these an appropriate approach? Or is there a more recommended or standard approach that I'm naively not aware of?


1 Answer 1


The multilevel model allows you to simultaneously model within- and between-plot variation. Consider a simplified multilevel model, with the Xs representing within-plot variables and W representing a between-plot variable.

The within-plot model: $y_{ij}$ = $\beta0_j$ + $\beta1_jX_{1ij}$ + $\beta2X_{2ij}...$ + $e_{ij}$

Here you can include as many X variables as are necessary to explain within-plot variation in your outcome. Recognize that X variables, because they are measured within-plots, also have between-plot variation, and thus you might also consider including their plot mean in the model for the between-plot mean outcome (intercept):

$\beta0_j$ = $\gamma_{00}$ + $\gamma_{01}W_j$ + $u_{0j}$

Where W could include plot means of all X variables as well as those plot-level variables you are interested in using as predictors of between-plot variance.

This model could be further extended by modeling between-plot variation in a within-plot slope. For example, the j in $\beta_{1j}$:

$\beta1_j$ (cluster slope) = $\gamma_{10}$ + $\gamma_{11}W_j$ + $u_{1j}$

Note that for both $\beta0_j$ and $\beta1_j$, the model assumes that these latent variables have a Normal distribution, with a mean of 0 (relative to the fixed effect parameter estimate) and an estimated variance ($\sigma^2_0$ and $\sigma^2_1$ for the intercept and slope, respectively).

Unless I am missing something, it seems that the multilevel model can accommodate both your desire to control for within-plot differences that might affect between-plot variance in the outcome and allow you to model between-plot variance in the mean level of the outcome. No need for a two-step approach where you mess with residuals from the first step.

  • $\begingroup$ I might not have been clear in my question. I understand your suggestion as adding plot-level variables to the multilevel model. What I am more interested in is whether there is a way to summarize the response variable (count data in groups of 5 observations) by groups (plot) that takes into account the "nuisance" observation-level variations. Like a mean, but "better". This is because I would then like to use this plot-level summary statistic to relate to other data in a plot-level structural equation model, some of which may be similarly summarized over their own levels of replication. $\endgroup$
    – K Li
    Jan 14, 2020 at 9:03
  • $\begingroup$ You can run a multilevel structural equation model that extends the multilevel regression framework I described above to SEM. Mplus is a really good program for this, but it is also possible to do it in R with lavaan. See lavaan.ugent.be/tutorial/multilevel.html $\endgroup$
    – Erik Ruzek
    Jan 14, 2020 at 20:21
  • $\begingroup$ Thanks for the suggestion to look at lavaan (I use R). I had been trying to use piecewiseSEM because it takes other distributions (as I understand it, lavaan assumes data are normally distributed) but piecewiseSEM doesn't have the capability to make multilevel SEMs, yet. $\endgroup$
    – K Li
    Feb 11, 2020 at 8:55
  • $\begingroup$ Just FYI, Mplus can deal with non-normal outcomes - binary, ordinal, count, zero-inflation, etc. Unfortunately it costs 895 US dollars and there is a 100ish dollars yearly fee to maintain the license. However you can still use the software after your license expires. $\endgroup$
    – Erik Ruzek
    Feb 11, 2020 at 12:45

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.