Comparison: Mandala vs lme4

Field trial analysis vs general-purpose mixed models

Overview

This tutorial compares the mandala package with lme4, the most widely-used R package for mixed models. While both can fit linear mixed models, they serve different purposes.

Background

Package	Primary Focus	Best For
lme4	General-purpose	Any discipline
Mandala	Agricultural trials	Field experiments

Loading Packages

# library(mandala)  # Install Mandala first (request download: https://listoagriculture.com/products)
library(lme4)
library(lmerTest)
library(dplyr)
library(ggplot2)

---

Example 1: Basic RCBD

Data Setup

set.seed(321)
n_genotypes <- 12
n_blocks <- 5

rcbd_data <- expand.grid(
  genotype = factor(paste0("G", sprintf("%02d", 1:n_genotypes))),
  block = factor(paste0("B", 1:n_blocks))
)

rcbd_data$yield <- 4500 +
  rnorm(n_genotypes, 0, 20)[as.numeric(rcbd_data$genotype)] +
  rnorm(n_blocks, 0, 10)[as.numeric(rcbd_data$block)] +
  rnorm(nrow(rcbd_data), 0, 15)

head(rcbd_data)

  genotype block    yield
1      G01    B1 4542.286
2      G02    B1 4496.306
3      G03    B1 4503.019
4      G04    B1 4513.249
5      G05    B1 4497.798
6      G06    B1 4522.072

Syntax Comparison

Mandala:

mandala(
  fixed  = yield ~ genotype,
  random = ~ block,
  data   = rcbd_data
)

lme4:

lmer(
  yield ~ genotype + (1|block),
  data = rcbd_data
)

lme4 Analysis

lme4_model <- lmer(
  yield ~ genotype + (1|block),
  data = rcbd_data,
  REML = TRUE
)

summary(lme4_model)

Linear mixed model fit by REML. t-tests use Satterthwaite's method [
lmerModLmerTest]
Formula: yield ~ genotype + (1 | block)
   Data: rcbd_data

REML criterion at convergence: 423

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.9341 -0.7533  0.1652  0.6741  1.7233 

Random effects:
 Groups   Name        Variance Std.Dev.
 block    (Intercept) 125.5    11.20   
 Residual             221.6    14.89   
Number of obs: 60, groups:  block, 5

Fixed effects:
            Estimate Std. Error       df t value Pr(>|t|)    
(Intercept) 4529.901      8.333   19.684 543.637  < 2e-16 ***
genotypeG02  -47.166      9.415   44.000  -5.010 9.34e-06 ***
genotypeG03  -33.401      9.415   44.000  -3.548 0.000938 ***
genotypeG04  -19.832      9.415   44.000  -2.106 0.040913 *  
genotypeG05  -28.048      9.415   44.000  -2.979 0.004693 ** 
genotypeG06  -16.668      9.415   44.000  -1.770 0.083595 .  
genotypeG07  -13.595      9.415   44.000  -1.444 0.155835    
genotypeG08  -24.406      9.415   44.000  -2.592 0.012896 *  
genotypeG09  -17.345      9.415   44.000  -1.842 0.072182 .  
genotypeG10  -43.653      9.415   44.000  -4.636 3.17e-05 ***
genotypeG11  -25.541      9.415   44.000  -2.713 0.009488 ** 
genotypeG12   -1.386      9.415   44.000  -0.147 0.883608    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) gntG02 gntG03 gntG04 gntG05 gntG06 gntG07 gntG08 gntG09
genotypeG02 -0.565                                                        
genotypeG03 -0.565  0.500                                                 
genotypeG04 -0.565  0.500  0.500                                          
genotypeG05 -0.565  0.500  0.500  0.500                                   
genotypeG06 -0.565  0.500  0.500  0.500  0.500                            
genotypeG07 -0.565  0.500  0.500  0.500  0.500  0.500                     
genotypeG08 -0.565  0.500  0.500  0.500  0.500  0.500  0.500              
genotypeG09 -0.565  0.500  0.500  0.500  0.500  0.500  0.500  0.500       
genotypeG10 -0.565  0.500  0.500  0.500  0.500  0.500  0.500  0.500  0.500
genotypeG11 -0.565  0.500  0.500  0.500  0.500  0.500  0.500  0.500  0.500
genotypeG12 -0.565  0.500  0.500  0.500  0.500  0.500  0.500  0.500  0.500
            gntG10 gntG11
genotypeG02              
genotypeG03              
genotypeG04              
genotypeG05              
genotypeG06              
genotypeG07              
genotypeG08              
genotypeG09              
genotypeG10              
genotypeG11  0.500       
genotypeG12  0.500  0.500

# Variance components
as.data.frame(VarCorr(lme4_model))

       grp        var1 var2     vcov    sdcor
1    block (Intercept) <NA> 125.5453 11.20470
2 Residual        <NA> <NA> 221.6143 14.88671

---

Example 2: Nested Design

set.seed(654)
n_genotypes <- 30
n_complete <- 3

nested_data <- expand.grid(
  genotype = factor(paste0("G", sprintf("%02d", 1:n_genotypes))),
  complete_block = factor(1:n_complete)
)

nested_data$incomplete_block <- factor(
  paste0(nested_data$complete_block, ".",
         rep(1:5, length.out = nrow(nested_data)))
)

nested_data$yield <- 5500 +
  rnorm(n_genotypes, 0, 30)[as.numeric(nested_data$genotype)] +
  rnorm(n_complete, 0, 15)[as.numeric(nested_data$complete_block)] +
  rnorm(nrow(nested_data), 0, 20)

lme4_nested <- lmer(
  yield ~ genotype + (1|complete_block) + (1|complete_block:incomplete_block),
  data = nested_data,
  REML = TRUE
)

as.data.frame(VarCorr(lme4_nested))

                              grp        var1 var2        vcov      sdcor
1 complete_block:incomplete_block (Intercept) <NA>  24.4848159  4.9482134
2                  complete_block (Intercept) <NA>   0.7170808  0.8468062
3                        Residual        <NA> <NA> 301.0144193 17.3497671

---

Diagnostics

diag_data <- data.frame(
  Fitted = fitted(lme4_model),
  Residuals = residuals(lme4_model)
)

p1 <- ggplot(diag_data, aes(x = Fitted, y = Residuals)) +
  geom_point(alpha = 0.6, color = "#5D4E6D") +
  geom_hline(yintercept = 0, linetype = "dashed", color = "#9B5B5B") +
  geom_smooth(se = FALSE, color = "#8B9A6D") +
  theme_minimal() +
  labs(title = "Residuals vs Fitted")

p2 <- ggplot(diag_data, aes(sample = Residuals)) +
  stat_qq(color = "#5D4E6D") +
  stat_qq_line(color = "#9B5B5B") +
  theme_minimal() +
  labs(title = "Normal Q-Q Plot")

gridExtra::grid.arrange(p1, p2, ncol = 2)

---

Key Differences

Task	Mandala	lme4
Random intercept	random = ~ block	(1\|block)
Variance components	$varcomp	VarCorr()
BLUPs	mandala_predict()	ranef()
Heritability	Built-in	Manual
Spatial models	Built-in	Extensions needed

When to Use Each

Use Mandala for:

• Agricultural field trials
• Spatial correlation modeling
• Quick heritability estimates

Use lme4 for:

• General mixed models
• Extensive diagnostics
• Large documentation/community

Session Information

sessionInfo()

R version 4.5.2 (2025-10-31)
Platform: aarch64-apple-darwin20
Running under: macOS Tahoe 26.1

Matrix products: default
BLAS:   /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib 
LAPACK: /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.1

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

time zone: America/New_York
tzcode source: internal

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] ggplot2_4.0.0  dplyr_1.1.4    lmerTest_3.1-3 lme4_1.1-37    Matrix_1.7-4  

loaded via a namespace (and not attached):
 [1] gtable_0.3.6        jsonlite_2.0.0      compiler_4.5.2     
 [4] tidyselect_1.2.1    Rcpp_1.1.0          gridExtra_2.3      
 [7] splines_4.5.2       scales_1.4.0        boot_1.3-32        
[10] yaml_2.3.10         fastmap_1.2.0       lattice_0.22-7     
[13] R6_2.6.1            labeling_0.4.3      generics_0.1.4     
[16] knitr_1.50          rbibutils_2.3       htmlwidgets_1.6.4  
[19] MASS_7.3-65         tibble_3.3.0        nloptr_2.2.1       
[22] minqa_1.2.8         pillar_1.11.1       RColorBrewer_1.1-3 
[25] rlang_1.1.6         xfun_0.54           S7_0.2.0           
[28] cli_3.6.5           mgcv_1.9-3          withr_3.0.2        
[31] magrittr_2.0.4      Rdpack_2.6.4        digest_0.6.37      
[34] grid_4.5.2          lifecycle_1.0.4     nlme_3.1-168       
[37] reformulas_0.4.2    vctrs_0.6.5         evaluate_1.0.5     
[40] glue_1.8.0          numDeriv_2016.8-1.1 farver_2.1.2       
[43] rmarkdown_2.30      pkgconfig_2.0.3     tools_4.5.2        
[46] htmltools_0.5.8.1