How to evaluate mixed model interactions? (Writing Custom tidier functions)

I wanted to do a power analysis of an interaction in a non-experimental dataset I have using lmer.

I can"tidy" my model using broom.mixed::tidy, but I think because of the extra parameters, it is throwing DD off. In essence, I think the problem is that I need to write a custom tidy function, and I wanted to get advice on the best way to do this.



pop <- declare_population(data = mtcars)
model <- 
  declare_estimator(mpg ~ wt*carb + (1|cyl),
                    model = lmer)

diagnose_design(pop + model, sims = 15)

#> Warning in nextpar(mat, cc, i, delta, lowcut, upcut): Last two rows have
#> identical or NA .zeta values: using minstep
#> Warning in FUN(X[[i]], ...): non-monotonic profile for .sig01
#> Warning in confint.thpr(pp, level = level, zeta = zeta): bad spline fit
#> for .sig01: falling back to linear interpolation
#> Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
#> collapsing to unique 'x' values
#> Error: Error in step 2 (estimator):
#>  Error in fit2tidy(results, coefficient_names): We were unable to tidy the output of the function provided to 'model'. 
#>            It is possible that the broom package has a tidier for that object type. 
#>            If not, you can use a custom estimator to 'estimator_function'.
#>            See examples in ?declare_estimator

edit: also see: Diagnosing multilevel models

Thank you so much for bringing this to our attention. I think it must be a bug because we’re supposed to recognize tidiers when they are available. Argh.

If you update to the latest version of DD on GitHub, this will work. We’ll also get to work making sure that it’ll auto-find that tidier.


pop <- declare_population(data = mtcars)
model <- 
  declare_estimator(mpg ~ wt*carb + (1|cyl),
                model = lmer, 
                model_summary = broom.mixed::tidy,
                term = TRUE)

diagnose_design(pop + model, sims = 15)
1 Like

This is a bug in some overly defensive parameter checking - it’s only looking for S3 methods, and not for S3 methods of parent classes of S4 objects.

That sounds confusing, because it is.

However, I think it’s a ~ 1 line change - I’ll PR it on github.

Unfortunately the above code doesn’t seem to work using DeclareDesign_0.23.0.

 Error in step 2 (estimator):
	Error in model(~(mpg ~ wt * carb + (1 | cyl)), model_summary = ~broom.mixed::tidy, : unused argument (model_summary = ~broom.mixed::tidy)

Is there another way of using the package for power analysis involving mixed models?

Hello my first post on this forum.

I want to revisit this question as I am fitting a similar model but I got a different error, seems like the “|” symbol for some reason doesn’t work. Below is my complete code and the error message. Any help will be greatly appreciated, and the analysis I am doing has a pretty tight deadline, so immediate help is needed as well. Thanks.
Z is my treatment variable, V29 is an interaction term, and C_COW_ALPHA is the group variable.


N ← nrow(work_data)

thepop ← declare_model(work_data)

they0y1 ← declare_potential_outcomes(Y~ 1.28372 + 0.06236Z + (-0.03071)V29 + 0.0608ZV29 + (1 | C_COW_ALPHA))

theassignment ← declare_assignment(Z=complete_ra(N),legacy=FALSE)

inquiry1 ← declare_inquiry(ATE = mean(Y_Z_1 - Y_Z_0))

revealY ← declare_reveal(Y,Z)

design ← thepop + they0y1 + theassignment + inquiry1 + revealY

E1c ← declare_estimator(Y ~ Z * V29 + (1 | C_COW_ALPHA), model = lmer, model_summary = broom.mixed::tidy, term = TRUE, REML = FALSE, label = “Test0”, inquiry = “ATE”)

des_plus_ests ← design + E1c

des1_sim ← simulate_design(des_plus_ests, sims = 15)

Error: Error in step 2 ():
** Error in 1 | C_COW_ALPHA: operations are possible only for numeric, logical or complex types**