Interaction plots for mlwinfitMCMC fixed effects

Welcome to the forum for R2MLwiN users. Feel free to post your question about R2MLwiN here. The Centre for Multilevel Modelling take no responsibility for the accuracy of these posts, we are unable to monitor them closely. Do go ahead and post your question and thank you in advance if you find the time to post any answers!

Go to R2MLwiN: Running MLwiN from within R >> http://www.bris.ac.uk/cmm/software/r2mlwin/
Post Reply
wfleming
Posts: 6
Joined: Thu Aug 13, 2020 9:18 pm

Interaction plots for mlwinfitMCMC fixed effects

Post by wfleming »

Any recommendation for the best way to plot fixed interaction effects from mlwinfitMCMC objects?

Model is simple like:

Code: Select all

m <- runmlwin(y ~ 1 + lev1x*lev1x + (1 | lev2id) + (1 | lev1id), estoptions = list(EstM = 1))
and want to plot the lev1x*lev1x

Would usually use plot_model from 'sjPlot' for lm objects but naturally not written for mlwin objects. Would also do it manually extracting predicted values but predict() doesn't give the intervals for mlwin objects.

Any ideas greatly appreciated.
ChrisCharlton
Posts: 1267
Joined: Mon Oct 19, 2009 10:34 am

Re: Interaction plots for mlwinfitMCMC fixed effects

Post by ChrisCharlton »

You could manually create the data required for this plot by generating a chain of predictions and then calculating the mean and credible intervals from this. An example of this would be as follows:

Code: Select all

library(R2MLwiN)
library(coda)

data(tutorial, package = "R2MLwiN")
(model1 <- runMLwiN(normexam ~ 1 + standlrt * vrband + (1 | school) + (1 | student), estoptions = list(EstM = 1), data = tutorial))

# Data frame for variables used in fixed part of the model
modeldata <- model.matrix(normexam ~ 1 + standlrt * vrband, data = tutorial)
head(modeldata)

# Data frame containing parameter chains (specify variables in the order corresponding with the data frame above)
modelchains <- as.matrix(model1@chains)[, c("FP_Intercept", "FP_standlrt", "FP_vrbandvb2", "FP_vrbandvb3", "FP_standlrt:vrbandvb2", "FP_standlrt:vrbandvb3")]

# Prediction chain
predictchain <- as.mcmc(modelchains %*% t(modeldata))

# Mean prediction
prediction <- drop(apply(predictchain, 2, mean))

# 95% credible intervals
quantiles <- t(apply(predictchain, 2, quantile, c(0.025, 0.975)))
wfleming
Posts: 6
Joined: Thu Aug 13, 2020 9:18 pm

Re: Interaction plots for mlwinfitMCMC fixed effects

Post by wfleming »

Thanks Chris this worked well. Wouldn't have worked it out without you.

Only change I made was to use bayestestR::hdi() to get credible intervals instead of the usual quantiles.
Post Reply