Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

Setting the right random effect part in mixed effect models can be tricky in many applied situation. I will not talk here about choosing wether a grouping variable (sites, individuals …) should be included as a fixed term or as a random term, please see Gelman and Hill (2006) and Zuur et al (2009) for informations. Here I will present the use of the bootMer function in the package lme4 to compare two models with different random effect terms specification and decide wether one model do a (significantly) better job at fitting the data. The standard way to compare two model is to derive the likelihood ratio test (LRT) value and since these should follow a chi-square distribution derive a p-value corresponding to the probability to observe such an extreme LRT under the null hypothesis that both model perform equally well. This approach works relatively fine for GLM but for (G)LMM several problem arises due mainly to boundary effects (the null hypothesis being in this case that the variance of the random effects is 0) see Bolker et al (2009). One way to do model comparison in (G)LMM is to derive bootstrapped likelihood values from the two competing models and to draw confidence intervals around the observed values to decide if one model perform better than the other. Below are some code with simulated data (a cleaner version with more graphs can be found here: http://rpubs.com/hughes/22059):

 

library(lme4)
library(arm)
library(RColorBrewer)

##### work on model comparison using bootMer ##### simulate some data and fit a
##### random intercept model to them
x <- runif(100, 0, 10)
# the grouping variable
site <- gl(n = 10, k = 10)
# the random intercept effect, the simulated standard deviation around the
# intercept is 1
rnd <- rnorm(10, 0, 1)
# the simulated resposne variable, note that the fixed effect coefficient
# are 1 for the intercept and 3 for the slope. Also the simulated residuals
# will have a standard deviation of one
y <- rep(1 + rnd, each = 10) + 3 * x + rnorm(100, 0, 1)
# fit the model using Maximum Likelihood to be able to use the LRT
m1 <- lmer(y ~ x + (1 | site), REML = FALSE)

# simulate to generate credible intervals
simu <- sim(m1, n.sims = 1000)
# a new model matrix with ordered and equally spaced predictor values
new.x <- model.matrix(~x, data = data.frame(x = seq(0, 10, length.out = 100)))
new.y <- matrix(ncol = 1000, nrow = 100)
# get the predicted response values for each 1000 simulations of the fixed
# effect model parameters
new.y <- apply([email protected], 1, function(x) new.x %*% x)
# compute the lower/upper quantile
lower <- apply(new.y, 1, function(x) quantile(x, prob = 0.025))
upper <- apply(new.y, 1, function(x) quantile(x, prob = 0.975))
median <- apply(new.y, 1, function(x) quantile(x, prob = 0.5))

# nice plot
pal <- brewer.pal(10, "RdYlBu")
plot(y ~ x, col = rep(pal, each = 10), pch = 16)
lines(new.x[, 2], median, col = "blue", lwd = 2)
lines(new.x[, 2], lower, col = "red", lwd = 2, lty = 2)
lines(new.x[, 2], upper, col = "red", lwd = 2, lty = 2)

# fit a second model with a random slope effect
m2 <- lmer(y ~ x + (x | site), REML = FALSE)

# using bootMer to compute 100 bootstrapped log-likelihood
b1 <- bootMer(m1, FUN = function(x) as.numeric(logLik(x)), nsim = 100)
b2 <- bootMer(m2, FUN = function(x) as.numeric(logLik(x)), nsim = 100)

# the observed LRT value
lrt <- as.numeric(-2 * logLik(m1) + 2 * logLik(m2))
# the 100 bootstrapped LRT
lrt.b <- -2 * b1$t + 2 * b2$t
# plot
quant <- quantile(lrt.b, probs = c(0.025, 0.975))
plot(1, lrt, xlab = "", ylab = "Likelihood ratio test", xaxt = "n", ylim = c(quant[1] +
1, quant[2] + 1))
abline(h = 0, lty = 2, lwd = 2, col = "red")
segments(1, quant[1], 1, quant[2], lend = 1)


 In the above example the 95% CI of the bootstrapped LRT cross the 0 line which means that one model do not fit the data better than the other. In this case the rule use would be to use the most simple model (the one with the lower number of parameters) which is the random-intercept model. Let's make another example: # now simulate data from random intercept/ slope rnd.slope <- rnorm(10, 0, 0.5) y <- rep(1 + rnd, each = 10) + rep(3 + rnd.slope, each = 10) * x + rnorm(100, 0, 1) # the new models m3 <- lmer(y ~ x + (x | site), REML = FALSE) m4 <- lmer(y ~ x + (1 | site), REML = FALSE) # LRT the observed values lrt <- -2 * logLik(m4) + 2 * logLik(m3) # the bootstrap b3 <- bootMer(m3, FUN = function(x) as.numeric(logLik(x)), nsim = 100) b4 <- bootMer(m4, FUN = function(x) as.numeric(logLik(x)), nsim = 100) # the 100 bootstrapped LRT lrt.b <- -2 * b4$t + 2 * b3$t # the nice plot quant <- quantile(lrt.b, probs = c(0.025, 0.975)) plot(1, lrt, xlab = "", ylab = "Likelihood ratio test", xaxt = "n", ylim = c(0, quant[2] + 1)) abline(h = 0, lty = 2, lwd = 2, col = "red") segments(1, quant[1], 1, quant[2], lend = 1) In this second example the random intercept/slope model fits much better to the data than the random intercept. This random effect structure should be kept. As mentionned in Bolker et al (2009) the LRT will be relevant depending on the design and the interest that is put on the random terms. In the case were random terms are due to the particular design of the study (site, blocks …) and when there are considered as a “nuisance” they may be included in the models without testing for the increase in fitness that their inclusion provide. In the case where the random term effects is of interest (individual sampling units …) then using LRT is a sensible way to detect and interprete the effect of the random terms. Biblio: Bolker, B. M., Brooks, M. E., Clark, C. J., Geange, S. W., Poulsen, J. R., Stevens, M. H. H., & White, J. S. S. (2009). Generalized linear mixed models: a practical guide for ecology and evolution. Trends in ecology & evolution, 24(3), 127-135. Gelman, A., & Hill, J. (2006). Data analysis using regression and multilevel/hierarchical models. Cambridge University Press. Zuur, A., Ieno, E. N., Walker, N., Saveliev, A. A., & Smith, G. M. (2009). Mixed effects models and extensions in ecology with R. Springer.Filed under: R and Stat Tagged: GLMM, LRT, R var vglnk = {key: '949efb41171ac6ec1bf7f206d57e90b8'}; (function(d, t) { var s = d.createElement(t); s.type = 'text/javascript'; s.async = true; // s.defer = true; // s.src = '//cdn.viglink.com/api/vglnk.js'; s.src = 'https://www.r-bloggers.com/wp-content/uploads/2020/08/vglnk.js'; var r = d.getElementsByTagName(t)[0]; r.parentNode.insertBefore(s, r); }(document, 'script')); Related ShareTweet To leave a comment for the author, please follow the link and comment on their blog: biologyforfun » R. R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job. Want to share your content on R-bloggers? click here if you have a blog, or here if you don't. 
 
 
 function init() { var vidDefer = document.getElementsByTagName('iframe'); for (var i=0; i<vidDefer.length; i++) { if(vidDefer[i].getAttribute('data-src')) { vidDefer[i].setAttribute('src',vidDefer[i].getAttribute('data-src')); } } } window.onload = init; R bloggers Facebook page Most viewed posts (weekly) Introduction to Machine Learning with TensorFlow PCA vs Autoencoders for Dimensionality Reduction 5 Ways to Subset a Data Frame in R How to write the first for loop in R Date Formats in R How to Remove Outliers in R R – Sorting a data frame by the contents of a column Sponsors // https://support.cloudflare.com/hc/en-us/articles/200169436-How-can-I-have-Rocket-Loader-ignore-my-script-s-in-Automatic-Mode- // this must be placed higher. Otherwise it doesn't work. // data-cfasync="false" is for making sure cloudflares' rocketcache doesn't interfeare with this // in this case it only works because it was used at the original script in the text widget function createCookie(name,value,days) { var expires = ""; if (days) { var date = new Date(); date.setTime(date.getTime() + (days*24*60*60*1000)); expires = "; expires=" + date.toUTCString(); } document.cookie = name + "=" + value + expires + "; path=/"; } function readCookie(name) { var nameEQ = name + "="; var ca = document.cookie.split(';'); for(var i=0;i < ca.length;i++) { var c = ca[i]; while (c.charAt(0)==' ') c = c.substring(1,c.length); if (c.indexOf(nameEQ) == 0) return c.substring(nameEQ.length,c.length); } return null; } function eraseCookie(name) { createCookie(name,"",-1); } // no longer use async because of google // async async function readTextFile(file) { // Helps people browse between pages without the need to keep downloading the same // ads txt page everytime. This way, it allows them to use their browser's cache. var random_number = readCookie("ad_random_number_cookie"); if(random_number == null) { var random_number = Math.floor(Math.random()*100*(new Date().getTime()/10000000000)); createCookie("ad_random_number_cookie",random_number,1) } file += '?t='+random_number; var rawFile = new XMLHttpRequest(); rawFile.onreadystatechange = function () { if(rawFile.readyState === 4) { if(rawFile.status === 200 || rawFile.status == 0) { // var allText = rawFile.responseText; // document.write(allText); document.write(rawFile.responseText); } } } rawFile.open("GET", file, false); rawFile.send(null); } // readTextFile('https://raw.githubusercontent.com/Raynos/file-store/master/temp.txt'); readTextFile("https://www.r-bloggers.com/wp-content/uploads/text-widget_anti-cache.txt"); Recent Posts Using Java functions from JAR file in R using rJava September 2021 ISC Call for Proposals – Now Open! Crash Course in R Model Deployment with Docker and friends Introduction to Deep Learning A Simple Two-Stage Stochastic Linear Programming using R RObservations #13: Simulating FSAs in lieu of real postal code data. Running R code for all combinations of some parameters with lapply karate error: JAVA_HOME cannot be determined from the Registry Competition to win free training closes today Solving Einstein’s Puzzle with Constraint Programming Advancements in Symbolic Data Analysis Using bootstrapped sampling to assess variability in score predictions An RStudio Table Contest for 2021 Advances in Difference-in-Differences in Econometrics Optimal disclosure risk assessment Jobs for R-usersJunior Data Scientist / Quantitative economistSenior Quantitative AnalystR programmerData Scientist – CGIAR Excellence in Agronomy (Ref No: DDG-R4D/DS/1/CG/EA/06/20)Data Analytics Auditor, Future of Audit Lead @ London or Newcastle python-bloggers.com (python/data-science news)3 Ways To Perform Quick Exploratory Data Analysis in PythonUsing the data algebra for Statistics and Data Sciencecalmcode.io > video tutorials for open source toolsApache Kafka in Python: How to Stream Data With Producers and ConsumersClassification using linear regressionTechnical skills or business skills… why not both? (Pragmatic Institute blog post)Roll up, roll up the NHS-R Community Conference 2021 is coming to town Full list of contributing R-bloggers Archives Archives Select Month October 2021  (3) September 2021  (184) August 2021  (145) July 2021  (167) June 2021  (195) May 2021  (196) April 2021  (183) March 2021  (215) February 2021  (218) January 2021  (242) December 2020  (284) November 2020  (256) October 2020  (268) September 2020  (242) August 2020  (218) July 2020  (259) June 2020  (227) May 2020  (326) April 2020  (327) March 2020  (279) February 2020  (261) January 2020  (253) December 2019  (244) November 2019  (219) October 2019  (233) September 2019  (233) August 2019  (273) July 2019  (260) June 2019  (245) May 2019  (277) April 2019  (293) March 2019  (308) February 2019  (262) January 2019  (287) December 2018  (257) November 2018  (289) October 2018  (308) September 2018  (291) August 2018  (270) July 2018  (333) June 2018  (298) May 2018  (321) April 2018  (301) March 2018  (291) February 2018  (241) January 2018  (330) December 2017  (261) November 2017  (270) October 2017  (290) September 2017  (294) August 2017  (340) July 2017  (283) June 2017  (317) May 2017  (349) April 2017  (324) March 2017  (365) February 2017  (317) January 2017  (367) December 2016  (347) November 2016  (294) October 2016  (306) September 2016  (254) August 2016  (287) July 2016  (327) June 2016  (263) May 2016  (292) April 2016  (260) March 2016  (302) February 2016  (268) January 2016  (337) December 2015  (304) November 2015  (234) October 2015  (259) September 2015  (238) August 2015  (264) July 2015  (243) June 2015  (213) May 2015  (235) April 2015  (211) March 2015  (259) February 2015  (212) January 2015  (245) December 2014  (236) November 2014  (221) October 2014  (218) September 2014  (259) August 2014  (217) July 2014  (235) June 2014  (241) May 2014  (243) April 2014  (260) March 2014  (289) February 2014  (269) January 2014  (263) December 2013  (264) November 2013  (241) October 2013  (234) September 2013  (215) August 2013  (223) July 2013  (254) June 2013  (272) May 2013  (260) April 2013  (279) March 2013  (277) February 2013  (294) January 2013  (343) December 2012  (308) November 2012  (277) October 2012  (308) September 2012  (270) August 2012  (263) July 2012  (247) June 2012  (298) May 2012  (287) April 2012  (295) March 2012  (304) February 2012  (264) January 2012  (280) December 2011  (251) November 2011  (261) October 2011  (281) September 2011  (187) August 2011  (258) July 2011  (219) June 2011  (225) May 2011  (239) April 2011  (268) March 2011  (249) February 2011  (205) January 2011  (209) December 2010  (188) November 2010  (172) October 2010  (219) September 2010  (185) August 2010  (203) July 2010  (175) June 2010  (167) May 2010  (164) April 2010  (152) March 2010  (165) February 2010  (135) January 2010  (121) December 2009  (126) November 2009  (66) October 2009  (87) September 2009  (65) August 2009  (56) July 2009  (64) June 2009  (54) May 2009  (35) April 2009  (38) March 2009  (40) February 2009  (33) January 2009  (42) December 2008  (16) November 2008  (14) October 2008  (10) September 2008  (8) August 2008  (11) July 2008  (7) June 2008  (8) May 2008  (8) April 2008  (4) March 2008  (5) February 2008  (6) January 2008  (10) December 2007  (3) November 2007  (5) October 2007  (9) September 2007  (7) August 2007  (21) July 2007  (9) June 2007  (3) May 2007  (3) April 2007  (1) March 2007  (5) February 2007  (4) November 2006  (1) October 2006  (2) August 2006  (3) July 2006  (1) June 2006  (1) May 2006  (3) April 2006  (1) March 2006  (1) February 2006  (5) January 2006  (1) October 2005  (1) September 2005  (3) May 2005  (1) /* <![CDATA[ */ (function() { var dropdown = document.getElementById( "archives-dropdown-3" ); function onSelectChange() { if ( dropdown.options[ dropdown.selectedIndex ].value !== '' ) { document.location.href = this.options[ this.selectedIndex ].value; } } dropdown.onchange = onSelectChange; })(); /* ]]> */ Other sites Jobs for R-users SAS blogs 
 
 Copyright © 2021 | MH Corporate basic by MH Themes 
 [{"@context":"https:\/\/schema.org","@graph":[{"@type":"Organization","@id":"https:\/\/www.r-bloggers.com#Organization","name":"R-bloggers","url":"https:\/\/www.r-bloggers.com","sameAs":[],"logo":{"@type":"ImageObject","url":"http:\/\/www.r-bloggers.com\/wp-content\/uploads\/2021\/05\/R_blogger_logo1_large.png","width":"1285","height":"369"},"contactPoint":{"@type":"ContactPoint","contactType":"technical support","telephone":"","url":"https:\/\/www.r-bloggers.com\/contact-us\/"}},{"@type":"WebSite","@id":"https:\/\/www.r-bloggers.com#website","headline":"R-bloggers","name":"R-bloggers","description":"R news and tutorials contributed by hundreds of R bloggers","url":"https:\/\/www.r-bloggers.com","potentialAction":{"@type":"SearchAction","target":"https:\/\/www.r-bloggers.com\/?s={search_term_string}","query-input":"required name=search_term_string"},"publisher":{"@id":"https:\/\/www.r-bloggers.com#Organization"}},{"@context":"https:\/\/schema.org","@type":"WebPage","@id":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/#webpage","name":"Using bootMer to do model comparison in R","url":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/","lastReviewed":"2014-07-13T09:14:48-06:00","inLanguage":"en-US","description":"Setting the right random effect part in mixed effect models can be tricky in many applied situation. I will not talk here about choosing wether a grouping variable (sites, individuals &hellip;) should be included as a fixed term or as a random term, please see Gelman and Hill (2006) and Zuur et al (2009) for ","reviewedBy":{"@type":"Organization","name":"R-bloggers","url":"https:\/\/www.r-bloggers.com","logo":{"@type":"ImageObject","url":"http:\/\/www.r-bloggers.com\/wp-content\/uploads\/2021\/05\/R_blogger_logo1_large.png","width":"1285","height":"369"}},"primaryImageOfPage":{"@id":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/#primaryimage"},"mainContentOfPage":[[{"@context":"https:\/\/schema.org","@type":"SiteNavigationElement","@id":"https:\/\/www.r-bloggers.com\/#top nav","name":"Home","url":"https:\/\/www.r-bloggers.com"},{"@context":"https:\/\/schema.org","@type":"SiteNavigationElement","@id":"https:\/\/www.r-bloggers.com\/#top nav","name":"About","url":"http:\/\/www.r-bloggers.com\/about\/"},{"@context":"https:\/\/schema.org","@type":"SiteNavigationElement","@id":"https:\/\/www.r-bloggers.com\/#top nav","name":"RSS","url":"https:\/\/feeds.feedburner.com\/RBloggers"},{"@context":"https:\/\/schema.org","@type":"SiteNavigationElement","@id":"https:\/\/www.r-bloggers.com\/#top nav","name":"add your blog!","url":"http:\/\/www.r-bloggers.com\/add-your-blog\/"},{"@context":"https:\/\/schema.org","@type":"SiteNavigationElement","@id":"https:\/\/www.r-bloggers.com\/#top nav","name":"Learn R","url":"https:\/\/www.r-bloggers.com\/how-to-learn-r-2\/"},{"@context":"https:\/\/schema.org","@type":"SiteNavigationElement","@id":"https:\/\/www.r-bloggers.com\/#top nav","name":"R jobs","url":"https:\/\/www.r-users.com\/"},{"@context":"https:\/\/schema.org","@type":"SiteNavigationElement","@id":"https:\/\/www.r-bloggers.com\/#top nav","name":"Submit a new job (it's free)","url":"https:\/\/www.r-users.com\/submit-job\/"},{"@context":"https:\/\/schema.org","@type":"SiteNavigationElement","@id":"https:\/\/www.r-bloggers.com\/#top nav","name":"Browse latest jobs (also free)","url":"https:\/\/www.r-users.com\/"},{"@context":"https:\/\/schema.org","@type":"SiteNavigationElement","@id":"https:\/\/www.r-bloggers.com\/#top nav","name":"Contact us","url":"http:\/\/www.r-bloggers.com\/contact-us\/"}]],"isPartOf":{"@id":"https:\/\/www.r-bloggers.com#website"},"breadcrumb":{"@id":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/#breadcrumb"}},{"@type":"BreadcrumbList","@id":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/#breadcrumb","itemListElement":[{"@type":"ListItem","position":1,"item":{"@id":"https:\/\/www.r-bloggers.com","name":"R-bloggers"}},{"@type":"ListItem","position":2,"item":{"@id":"https:\/\/www.r-bloggers.com\/category\/r-bloggers\/","name":"R bloggers"}},{"@type":"ListItem","position":3,"item":{"@id":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/","name":"Using bootMer to do model comparison in R"}}]},{"@type":"Article","@id":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/#Article","url":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/","inLanguage":"en-US","mainEntityOfPage":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/#webpage","headline":"Using bootMer to do model comparison in R","description":"Setting the right random effect part in mixed effect models can be tricky in many applied situation. I will not talk here about choosing wether a grouping variable (sites, individuals &hellip;) should be included as a fixed term or as a random term, please see Gelman and Hill (2006) and Zuur et al (2009) for ","articleBody":"Setting the right random effect part in mixed effect models can be tricky in many applied situation. I will not talk here about choosing wether a grouping variable (sites, individuals \u2026) should be included as a fixed term or as a random term, please see Gelman and Hill (2006) and Zuur et al (2009) for informations. Here I will present the use of the bootMer function in the package lme4 to compare two models with different random effect terms specification and decide wether one model do a (significantly) better job at fitting the data. The standard way to compare two model is to derive the likelihood ratio test (LRT) value and since these should follow a chi-square distribution derive a p-value corresponding to the probability to observe such an extreme LRT under the null hypothesis that both model perform equally well. This approach works relatively fine for GLM but for (G)LMM several problem arises due mainly to boundary effects (the null hypothesis being in this case that the variance of the random effects is 0) see Bolker et al (2009). One way to do model comparison in (G)LMM is to derive bootstrapped likelihood values from the two competing models and to draw confidence intervals around the observed values to decide if one model perform better than the other. Below are some code with simulated data (a cleaner version with more graphs can be found here: http:\/\/rpubs.com\/hughes\/22059): library(lme4) library(arm) library(RColorBrewer) ##### work on model comparison using bootMer ##### simulate some data and fit a ##### random intercept model to them x &lt;- runif(100, 0, 10) # the grouping variable site &lt;- gl(n 10, k 10) # the random intercept effect, the simulated standard deviation around the # intercept is 1 rnd &lt;- rnorm(10, 0, 1) # the simulated resposne variable, note that the fixed effect coefficient # are 1 for the intercept and 3 for the slope. Also the simulated residuals # will have a standard deviation of one y &lt;- rep(1 + rnd, each 10) + 3 * x + rnorm(100, 0, 1) # fit the model using Maximum Likelihood to be able to use the LRT m1 &lt;- lmer(y ~ x + (1 | site), REML FALSE) # simulate to generate credible intervals simu &lt;- sim(m1, n.sims 1000) # a new model matrix with ordered and equally spaced predictor values new.x &lt;- model.matrix(~x, data data.frame(x seq(0, 10, length.out 100))) new.y &lt;- matrix(ncol 1000, nrow 100) # get the predicted response values for each 1000 simulations of the fixed # effect model parameters new.y &lt;- apply(simu@fixef, 1, function(x) new.x %*% x) # compute the lower\/upper quantile lower &lt;- apply(new.y, 1, function(x) quantile(x, prob 0.025)) upper &lt;- apply(new.y, 1, function(x) quantile(x, prob 0.975)) median &lt;- apply(new.y, 1, function(x) quantile(x, prob 0.5)) # nice plot pal &lt;- brewer.pal(10, &quot;RdYlBu&quot;) plot(y ~ x, col rep(pal, each 10), pch 16) lines(new.x, median, col &quot;blue&quot;, lwd 2) lines(new.x, lower, col &quot;red&quot;, lwd 2, lty 2) lines(new.x, upper, col &quot;red&quot;, lwd 2, lty 2) # fit a second model with a random slope effect m2 &lt;- lmer(y ~ x + (x | site), REML FALSE) # using bootMer to compute 100 bootstrapped log-likelihood b1 &lt;- bootMer(m1, FUN function(x) as.numeric(logLik(x)), nsim 100) b2 &lt;- bootMer(m2, FUN function(x) as.numeric(logLik(x)), nsim 100) # the observed LRT value lrt &lt;- as.numeric(-2 * logLik(m1) + 2 * logLik(m2)) # the 100 bootstrapped LRT lrt.b &lt;- -2 * b1$t + 2 * b2$t # plot quant &lt;- quantile(lrt.b, probs c(0.025, 0.975)) plot(1, lrt, xlab &quot;&quot;, ylab &quot;Likelihood ratio test&quot;, xaxt &quot;n&quot;, ylim c(quant + 1, quant + 1)) abline(h 0, lty 2, lwd 2, col &quot;red&quot;) segments(1, quant, 1, quant, lend 1) In the above example the 95% CI of the bootstrapped LRT cross the 0 line which means that one model do not fit the data better than the other. In this case the rule use would be to use the most simple model (the one with the lower number of parameters) which is the random-intercept model. Let's make another example: # now simulate data from random intercept\/ slope rnd.slope &lt;- rnorm(10, 0, 0.5) y &lt;- rep(1 + rnd, each 10) + rep(3 + rnd.slope, each 10) * x + rnorm(100, 0, 1) # the new models m3 &lt;- lmer(y ~ x + (x | site), REML FALSE) m4 &lt;- lmer(y ~ x + (1 | site), REML FALSE) # LRT the observed values lrt &lt;- -2 * logLik(m4) + 2 * logLik(m3) # the bootstrap b3 &lt;- bootMer(m3, FUN function(x) as.numeric(logLik(x)), nsim 100) b4 &lt;- bootMer(m4, FUN function(x) as.numeric(logLik(x)), nsim 100) # the 100 bootstrapped LRT lrt.b &lt;- -2 * b4$t + 2 * b3$t # the nice plot quant &lt;- quantile(lrt.b, probs c(0.025, 0.975)) plot(1, lrt, xlab &quot;&quot;, ylab &quot;Likelihood ratio test&quot;, xaxt &quot;n&quot;, ylim c(0, quant + 1)) abline(h 0, lty 2, lwd 2, col &quot;red&quot;) segments(1, quant, 1, quant, lend 1) In this second example the random intercept\/slope model fits much better to the data than the random intercept. This random effect structure should be kept. As mentionned in Bolker et al (2009) the LRT will be relevant depending on the design and the interest that is put on the random terms. In the case were random terms are due to the particular design of the study (site, blocks \u2026) and when there are considered as a \u201cnuisance\u201d they may be included in the models without testing for the increase in fitness that their inclusion provide. In the case where the random term effects is of interest (individual sampling units \u2026) then using LRT is a sensible way to detect and interprete the effect of the random terms. Biblio: Bolker, B. M., Brooks, M. E., Clark, C. J., Geange, S. W., Poulsen, J. R., Stevens, M. H. H., &amp; White, J. S. S. (2009). Generalized linear mixed models: a practical guide for ecology and evolution. Trends in ecology &amp; evolution, 24(3), 127-135. Gelman, A., &amp; Hill, J. (2006). Data analysis using regression and multilevel\/hierarchical models. Cambridge University Press. Zuur, A., Ieno, E. N., Walker, N., Saveliev, A. A., &amp; Smith, G. M. (2009). Mixed effects models and extensions in ecology with R. Springer.Filed under: R and Stat Tagged: GLMM, LRT, R","keywords":"","datePublished":"2014-07-13T09:14:48-06:00","dateModified":"2014-07-13T09:14:48-06:00","author":{"@type":"Person","name":"grumble10","description":"","url":"https:\/\/www.r-bloggers.com\/author\/grumble10\/","sameAs":["http:\/\/biologyforfun.wordpress.com"],"image":{"@type":"ImageObject","url":"https:\/\/secure.gravatar.com\/avatar\/ed235250c11f89cf0e6aa8e167251e53?s=96&d=mm&r=g","height":96,"width":96}},"editor":{"@type":"Person","name":"grumble10","description":"","url":"https:\/\/www.r-bloggers.com\/author\/grumble10\/","sameAs":["http:\/\/biologyforfun.wordpress.com"],"image":{"@type":"ImageObject","url":"https:\/\/secure.gravatar.com\/avatar\/ed235250c11f89cf0e6aa8e167251e53?s=96&d=mm&r=g","height":96,"width":96}},"publisher":{"@id":"https:\/\/www.r-bloggers.com#Organization"},"image":[{"@type":"ImageObject","url":"http:\/\/biologyforfun.files.wordpress.com\/2014\/07\/fig11.png?w=497&#038;h=326","width":497,"height":326,"@id":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/#primaryimage"},{"@type":"ImageObject","url":"https:\/\/feeds.wordpress.com\/1.0\/comments\/biologyforfun.wordpress.com\/295\/","width":0,"height":0},{"@type":"ImageObject","url":"https:\/\/feeds.wordpress.com\/1.0\/delicious\/biologyforfun.wordpress.com\/295\/","width":0,"height":0},{"@type":"ImageObject","url":"https:\/\/feeds.wordpress.com\/1.0\/facebook\/biologyforfun.wordpress.com\/295\/","width":0,"height":0},{"@type":"ImageObject","url":"https:\/\/feeds.wordpress.com\/1.0\/twitter\/biologyforfun.wordpress.com\/295\/","width":0,"height":0},{"@type":"ImageObject","url":"https:\/\/feeds.wordpress.com\/1.0\/stumble\/biologyforfun.wordpress.com\/295\/","width":0,"height":0},{"@type":"ImageObject","url":"https:\/\/feeds.wordpress.com\/1.0\/digg\/biologyforfun.wordpress.com\/295\/","width":0,"height":0},{"@type":"ImageObject","url":"https:\/\/feeds.wordpress.com\/1.0\/reddit\/biologyforfun.wordpress.com\/295\/","width":0,"height":0},{"@type":"ImageObject","url":"http:\/\/pixel.wp.com\/b.gif?host=biologyforfun.wordpress.com&#038;blog=42461605&#038;%23038;post=295&#038;%23038;subd=biologyforfun&#038;%23038;ref=&#038;%23038;feed=1","width":1,"height":1}],"isPartOf":{"@id":"https:\/\/www.r-bloggers.com\/2014\/07\/using-bootmer-to-do-model-comparison-in-r\/#webpage"}}]}] var snp_f = []; var snp_hostname = new RegExp(location.host); var snp_http = new RegExp("^(http|https)://", "i"); var snp_cookie_prefix = ''; var snp_separate_cookies = false; var snp_ajax_url = 'https://www.r-bloggers.com/wp-admin/admin-ajax.php'; var snp_ajax_nonce = '590e19d30d'; var snp_ignore_cookies = false; var snp_enable_analytics_events = true; var snp_enable_mobile = false; var snp_use_in_all = false; var snp_excluded_urls = []; Never miss an update! Subscribe to R-bloggers to receive e-mails with the latest R posts. (You will not see this message again.) Click here to close (This popup will not appear again) .snp-pop-109583 .snp-theme6 { max-width: 700px;} .snp-pop-109583 .snp-theme6 h1 {font-size: 17px;} .snp-pop-109583 .snp-theme6 { color: #a0a4a9;} .snp-pop-109583 .snp-theme6 .snp-field ::-webkit-input-placeholder { color: #a0a4a9;} .snp-pop-109583 .snp-theme6 .snp-field :-moz-placeholder { color: #a0a4a9;} .snp-pop-109583 .snp-theme6 .snp-field :-ms-input-placeholder { color: #a0a4a9;} .snp-pop-109583 .snp-theme6 .snp-field input { border: 1px solid #a0a4a9;} .snp-pop-109583 .snp-theme6 .snp-field { color: #000000;} .snp-pop-109583 .snp-theme6 { background: #f2f2f2;} jQuery(document).ready(function() { }); var CaptchaCallback = function() { jQuery('.g-recaptcha').each(function(index, el) { grecaptcha.render(el, { 'sitekey' : '' }); }); }; /* <![CDATA[ */!function(e,n){var r={"selectors":{"block":"pre","inline":"code"},"options":{"indent":4,"ampersandCleanup":true,"linehover":true,"rawcodeDbclick":false,"textOverflow":"scroll","linenumbers":false,"theme":"enlighter","language":"r","retainCssClasses":false,"collapse":false,"toolbarOuter":"","toolbarTop":"{BTN_RAW}{BTN_COPY}{BTN_WINDOW}{BTN_WEBSITE}","toolbarBottom":""},"resources":["https:\/\/www.r-bloggers.com\/wp-content\/plugins\/enlighter\/cache\/enlighterjs.min.css?4WJrVky+dDEQ83W","https:\/\/www.r-bloggers.com\/wp-content\/plugins\/enlighter\/resources\/enlighterjs\/enlighterjs.min.js"]},o=document.getElementsByTagName("head")[0],t=n&&(n.error||n.log)||function(){};e.EnlighterJSINIT=function(){!function(e,n){var r=0,l=null;function c(o){l=o,++r==e.length&&(!0,n(l))}e.forEach(function(e){switch(e.match(/\.([a-z]+)(?:[#?].*)?$/)[1]){case"js":var n=document.createElement("script");n.onload=function(){c(null)},n.onerror=c,n.src=e,n.async=!0,o.appendChild(n);break;case"css":var r=document.createElement("link");r.onload=function(){c(null)},r.onerror=c,r.rel="stylesheet",r.type="text/css",r.href=e,r.media="all",o.appendChild(r);break;default:t("Error: invalid file extension",e)}})}(r.resources,function(e){e?t("Error: failed to dynamically load EnlighterJS resources!",e):"undefined"!=typeof EnlighterJS?EnlighterJS.init(r.selectors.block,r.selectors.inline,r.options):t("Error: EnlighterJS resources not loaded yet!")})},(document.querySelector(r.selectors.block)||document.querySelector(r.selectors.inline))&&e.EnlighterJSINIT()}(window,console); /* ]]> */ window.FPConfig= { delay: 0, ignoreKeywords: ["\/wp-admin","\/wp-login.php","\/cart","add-to-cart","logout","#","?",".png",".jpeg",".jpg",".gif",".svg"], maxRPS: 3, hoverDelay: 50 }; _stq = window._stq || []; _stq.push([ 'view', {v:'ext',j:'1:7.3.3',blog:'11524731',post:'77502',tz:'-6',srv:'www.r-bloggers.com'} ]); _stq.push([ 'clickTrackerInit', '11524731', '77502' ]); jQuery(document).ready(function ($) { for (let i = 0; i < document.forms.length; ++i) { let form = document.forms[i]; if ($(form).attr("method") != "get") {$(form).append('<input type="hidden" name="LHfQMdbo_Al" value="8I4awxr7kvVDt" />'); } if ($(form).attr("method") != "get") {$(form).append('<input type="hidden" name="zxVu-DITS" value="Nce7adstg1[H]" />'); } if ($(form).attr("method") != "get") {$(form).append('<input type="hidden" name="fFcgQ_OCSUuEzl" value="DRK3bIW0w8L" />'); } } $(document).on('submit', 'form', function () { if ($(this).attr("method") != "get") { $(this).append('<input type="hidden" name="LHfQMdbo_Al" value="8I4awxr7kvVDt" />'); } if ($(this).attr("method") != "get") { $(this).append('<input type="hidden" name="zxVu-DITS" value="Nce7adstg1[H]" />'); } if ($(this).attr("method") != "get") { \$(this).append('<input type="hidden" name="fFcgQ_OCSUuEzl" value="DRK3bIW0w8L" />'); } return true; }); jQuery.ajaxSetup({ beforeSend: function (e, data) { if (data.type !== 'POST') return; if (typeof data.data === 'object' && data.data !== null) { data.data.append("LHfQMdbo_Al", "8I4awxr7kvVDt"); data.data.append("zxVu-DITS", "Nce7adstg1[H]"); data.data.append("fFcgQ_OCSUuEzl", "DRK3bIW0w8L"); } else { data.data = data.data + '&LHfQMdbo_Al=8I4awxr7kvVDt&zxVu-DITS=Nce7adstg1[H]&fFcgQ_OCSUuEzl=DRK3bIW0w8L'; } } }); });