Multilevel Descriptive Statistics

Reading the Matrix

Upper triangle — within-person correlations: Each coefficient answers: for a given person, on days when they score higher than their own average on variable X, do they also tend to score higher than their own average on variable Y? For example, the cell where stress meets wellbeing in the upper triangle tells us that days with higher stress than usual also tend to be days with lower wellbeing for the same person. For time-invariant variables (e.g., trait self-control), within-group correlations cannot be estimated and appear as NA in the upper triangle.

Lower triangle — between-person correlations: Each coefficient answers: do persons who score higher on variable X on average also tend to score higher on variable Y on average? For example, the cell where stress meets wellbeing in the lower triangle tells us that people who are more stressed on average also tend to have worse wellbeing overall.

When Within and Between Diverge

Importantly, within-group and between-group correlations can sometimes differ in direction:

within_between_correlations(
  data  = media_diary,
  group = "person",
  vars  = c("wellbeing", "screen_time")
)
#> # Within- and Between-Person Correlations
#>   =========== ====== =====
#>   variable       `1`   `2`
#>   ----------- ------ -----
#> 1 wellbeing        – 0.42*
#> 2 screen_time -0.34*     –
#>   =========== ====== =====
#> # ℹ Within-person correlations above, between-person correlations below the
#> #   diagonal.
#> # ℹ All correlations marked with a star are significant at p < .05.
#> # ℹ Correlations estimated via variance decomposition.

The two correlations for screen_time and wellbeing tell different stories:

• Within-person (upper triangle): Positive. On days when a person watches more entertainment media than they typically do, their wellbeing tends to be slightly higher than their own average — consistent with short-term escapism or mood repair.
• Between-person (lower triangle): Negative. People who watch more entertainment media on average tend to have lower average wellbeing — perhaps because chronic heavy media use reflects lower trait self-control.

Now look at what a naive pooled correlation gives:

cor(media_diary$screen_time, media_diary$wellbeing)
#> [1] 0.04597101

The pooled correlation is near zero. A researcher relying only on this number would conclude that entertainment media use has no relationship with wellbeing — missing a real positive within-person effect and a real negative between-person association that point in opposite directions.

Estimation Method and Weighting

Three estimation methods for within- and between-group correlations are available via the method argument:

• method = "decomposition" (default): A fast, closed-form variance-decomposition approach (Pedhazur, 1997). Within-group correlations are computed from group-mean-centred scores; between-group correlations are computed from the group means.
• method = "sem": Fits a two-level structural equation model using robust maximum likelihood. This handles very unequal group sizes more rigorously for significance testing.
• method = "bayes": Fits Bayesian multivariate models via brms, reporting credible intervals instead of p-values. See vignette("correlation-methods") for the statistical rationale and “Bayesian Estimation” below for the additional arguments it requires.

within_between_correlations(
  data   = media_diary,
  group  = "person",
  vars   = vars,
  method = "sem" # or "bayes"
)

The weight argument (decomposition method only) controls how much influence each group has on the between-group estimate:

• weight = TRUE (default): Between-group correlations and descriptive statistics are weighted by group size. Recommended for unbalanced data to recover the correct population-level estimate (Snijders & Bosker, 2012).
• weight = FALSE: Every group counts equally regardless of size.

Note that significance tests for between-group correlations always use unweighted group means regardless of the weight setting — see vignette("correlation-methods") for the statistical rationale.

within_between_correlations(
  data   = media_diary,
  group  = "person",
  vars   = vars,
  weight = FALSE
)

Customising the Output

Several options control the appearance:

mldesc(
  data                = media_diary,
  group               = "person",
  vars                = vars,
  significance        = "detailed",  # *, **, *** for p < .05, .01, .001
  flip                = TRUE,        # between above diagonal, within below
  remove_leading_zero = FALSE        # keep "0.45" instead of ".45"
)
#> # Multilevel Descriptive Statistics
#>   ============ ===== ====== ===== ===== ===== ======== ======== ======== =====
#>   variable     n_obs      m    sd range   `1`      `2`      `3`      `4`   icc
#>   ------------ ----- ------ ----- ----- ----- -------- -------- -------- -----
#> 1 Self control 1,400   4.03  0.83   2–6     –  0.61*** -0.67*** -0.53***  1.00
#> 2 Wellbeing    1,400   4.45  0.87   2–7    NA        – -0.34*** -0.38***  0.46
#> 3 Screen time  1,400 128.66 42.29 0–272    NA  0.42***        –  0.38***  0.45
#> 4 Stress       1,400   3.81  0.91   1–7    NA -0.43***  0.29***        –  0.33
#>   ============ ===== ====== ===== ===== ===== ======== ======== ======== =====
#> # ℹ Between-person correlations above, within-person correlations below the
#> #   diagonal.
#> # ℹ Correlations marked with * are significant at p < .05, ** at p < .01, and
#> #   *** at p < .001.
#> # ℹ Correlations estimated via variance decomposition.
#> # ℹ Group-weighted multilevel descriptive statistics computed with mlstats.

• significance = "detailed": Adds one star for p < .05, two for p < .01, and three for p < .001. The default ("basic") marks only p < .05.
• flip = TRUE: Swaps the positions of within- and between-person correlations (between above the diagonal, within below).
• remove_leading_zero = FALSE: Keeps the leading zero in decimal numbers. The default removes it for APA formatting (.45 instead of 0.45).

The method and weight arguments can also be used and work exactly as for within_between_correlations() (see above).

Pipe-Friendly Output

Although the main purpose is to produce publication-ready tables, the output of both mldesc() and within_between_correlations() are ordinary tibbles. The correlation and ICC columns use the class mlstats_stat (abbreviated mls) to store formatted strings. To use these values in further calculations, cast them to numeric, which strips formatting and significance stars:

result_num <- result[c(1, 6:10)]
result_num[-1] <- lapply(result[6:10], as.numeric)
as_tibble(result_num)
#> # A tibble: 4 × 6
#>   variable       `1`   `2`   `3`   `4`   icc
#>   <chr>        <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 Self control NA    NA    NA    NA     1   
#> 2 Wellbeing     0.61 NA     0.42 -0.43  0.46
#> 3 Screen time  -0.67 -0.34 NA     0.29  0.45
#> 4 Stress       -0.53 -0.38  0.38 NA     0.33