# How can I get a best model? An exploratory LMM

I'd like to inquire about the linear mixed model and its application to my dataset. The dataset comprises a dependent variable (DV) denoted as V, alongside three associated demographic variables (e.g., age and sex labeled as d1-d3), a single within-subject variable identified as wi, and a single between-subject variable as bw. The variable id represents unique subject identifiers within the dataset.

Now, I set the initial full model as:

V ~ d1 * d2 * d3 * wi * bw + (1 + wi | id)


Then I used step() to identify the model that best explains the variation in the "V", which resulted in finding a complex model such as:

V ~ d1 + d2 + d3 + wi + bw + (1 + wi | id) + d1:d2 + d1:d3 + ...


Along with 13 combinations with each variable. Then I evaluated the model fit with this complex model. It has a lot of fixed effects but only some displayed statistical significance (let's say d1, d1:wi, and d2:wi:bw).

Can I trust this model and conclude that the significant fixed effects (d1, d1:wi, and d2:wi:bw) are important factors for the DV? Are there better ways to find which variables or their interactions were important for the DV? Note that the dataset has 65 subjects so I wonder if the model found by step() has too many explanatory variables.

• Stepwise is not a good method. This has been discussed here many times. See this post stats.stackexchange.com/questions/20836/… in particular. Commented Dec 27, 2023 at 17:15
• What is your research question ? Are you interested in prediction or inference ? Commented Dec 30, 2023 at 18:59
• Please edit to give your question a title that is an informative summary of what you would like to know.
– Sycorax
Commented Dec 30, 2023 at 19:48
• Thank you everyone! So, is there no way to find an essential "best model" that matches the given data exploratory and automatically, like the stepwise method (I thought it was a good way to find a better model automatically)?
– TKw
Commented Jan 2 at 17:59

With only 65 subjects, the model you "found":

V ~ d1 + d2 + d3 + wi + bw + (1 + wi | id) + d1:d2 + d1:d3 + ... + 13 combinations


is likely to be over-fitted. As mentioned by Peter Flom in the comments to the question, a stepwise approach is very bad and should be avoided.

You started with a model that included a 5-way interaction and all lower order interactions. Such a model would be almost impossible to interpret. Just playing devil's advocate for a minute, and with tongue firmly in cheek, you could have also added nonlinear terms (eg quadratic, cubic, orthogonal splines etc) for all the variables too.

The best way to proceed in my opinion is to draw a causal diagram (a DAG), using expert knowledge of the data, which will encapsulate the causal relationships between the variables. The DAG will then inform the model. Note that the DAG is non-parametric so it will not tell you anything about interactions or nonlinear associations. Again, these should be informed by expert knowledge. See the following thread for more details about using DAGs to build a model

How do DAGs help to reduce bias in causal inference?

• That link is an excellent read. +1 to both answers. Commented Dec 31, 2023 at 5:15

#### Trusting Models

Your primary question from what I gather is this one:

Can I trust this model and conclude that the significant fixed effects are important factors for the DV?

It is probably obvious that no model should be trusted on statistical significance alone. There are many reasons for this. Statistical significance is not practical significance (see simulations on this point for the univariate case and the bivariate case). Statistical significance can also be completely spurious if the model has a poor fit (such as overfitting or omitted variable bias), or has theoretically poor grounds for inclusion of effects. Fitting a bunch of terms and somehow getting statistical significance (issues I note further below) are not the wisest way to determine which effects are the most influential.

Is there a better way to find which variables or their interactions were important for the DV?

The best ways are:

• Understanding the previous literature on these effects and making assumptions based on this prior information.
• Exploring the data to determine what is driving certain parts of your dataset before fitting (even with simple stuff like clearly erroneous outliers).
• Modeling the best functional form to the data which is plainly visible from plotting.
• Getting counsel from your supervisor, peers, or subject matter experts when possible (they may have their own background knowledge that is helpful).

There is likely more to it than that, but those immediately stick out to me as things you should be doing already.

#### Interactions Galore!

On to another important point, I agree with Robert that having a five-way interaction is difficult to justify, and including a combination of several main effects and interactions isn't inherently better.

First, you simply aren't going to have any statistical power with just $$n = 65$$, as the effects from interactions are necessarily weak. Generally you need 4 times the sample size to estimate an interaction that is the same size as the main effect, because the standard error effectively doubles for the coefficient (Gelman et al., 2022). Even if you get a statistically significant effect, it could be completely attributable to noise in the model. Random slope models (like yours) with insufficient variance structures to support them also tend to have less power (Matuschek et al., 2017).

Second, what is the theoretical justification for such an interaction? Is it really quite clear that these five variables are going to have a strong relationship that merits their inclusion? If there is a relationship between them, is it even strong enough to merit consideration? You should carefully consider these points before moving further with your analysis, but my hunch is that including a five-way interaction like yours here isn't going to be helpful. In fact, including all the variables as just main effects can be equally problematic if they don't make sense to add (especially given you are estimating so many coefficients with so few people). To that purpose, I recommend reading through Robert's useful post on causal diagrams like DAGs to determine what really should be included. As Cohen (1990) once noted, "less is more, except of course for sample size." Knowing which effects to model and simplifying down to the core theoretical question you are trying to answer is often better than slamming a bunch of predictors into a model and hoping for statistically significant slopes.

P.S. As Peter and Robert already noted, stepwise has been historically known as a bad idea. For a simulation study on the matter, check out Smith (2018) below. Just don't bother using it, as it too often capitalizes on chance and doesn't account enough for relationships in your model (such as predictors that are correlated with each other).