Visualization: Collapse by [random] group

[mattansb]

Are we allowed to steal from ourselves?

https://strengejacke.github.io/ggeffects/reference/collapse_by_group.html

[DominiqueMakowski]

library(ggeffects)
library(lme4)
library(dplyr)


data(efc)
efc$e15relat <- as.factor(efc$e15relat)
efc$c161sex <- as.factor(efc$c161sex)
levels(efc$c161sex) <- c("male", "female")
model <- lmer(neg_c_7 ~ c161sex + (1 | e15relat), data = efc)

# ggpredict
me <- ggpredict(model, terms = "c161sex")

collapse_by_group(me, model, "e15relat")
#>    x group_col facet random  response
#> 1  1         1     1      1 12.297872
#> 2  2         1     1      1 13.347107
#> 3  1         1     1      2 11.585586
#> 4  2         1     1      2 12.118310
#> 5  1         1     1      3 12.166667
#> 6  2         1     1      3 10.545455
#> 7  1         1     1      4 10.750000
#> 8  2         1     1      4 11.726027
#> 9  1         1     1      5 11.333333
#> 10 2         1     1      5 10.235294
#> 11 1         1     1      6  8.200000
#> 12 2         1     1      6  9.235294
#> 13 1         1     1      7 13.000000
#> 14 2         1     1      7 10.400000
#> 15 1         1     1      8  9.666667
#> 16 2         1     1      8 10.955882



# Modelbased
pr <- modelbased::estimate_expectation(model)

dat <- insight::get_data(attributes(pr)$model)
summarise_if(group_by(dat, e15relat), is.numeric, mean)
#> # A tibble: 8 x 2
#>   e15relat neg_c_7
#>   <fct>      <dbl>
#> 1 1           13.1
#> 2 2           12.0
#> 3 3           10.9
#> 4 4           11.6
#> 5 5           10.5
#> 6 6            9  
#> 7 7           10.8
#> 8 8           10.6

^{Created on 2021-07-06 by the reprex package (v2.0.0)}

lol what is x in ggeffects' output?

[mattansb]

library(ggeffects)
library(lme4)
#> Loading required package: Matrix
library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union


data(efc)
efc$e15relat <- as.factor(efc$e15relat)
efc$c161sex <- as.factor(efc$c161sex)
levels(efc$c161sex) <- c("male", "female")
model <- lmer(neg_c_7 ~ c161sex + (1 | e15relat), data = efc)

# ggpredict
me <- ggpredict(model, terms = "c161sex")
collapse_by_group(me, model, "e15relat")
#>    x group_col facet random  response
#> 1  1         1     1      1 12.297872
#> 2  2         1     1      1 13.347107
#> 3  1         1     1      2 11.585586
#> 4  2         1     1      2 12.118310
#> 5  1         1     1      3 12.166667
#> 6  2         1     1      3 10.545455
#> 7  1         1     1      4 10.750000
#> 8  2         1     1      4 11.726027
#> 9  1         1     1      5 11.333333
#> 10 2         1     1      5 10.235294
#> 11 1         1     1      6  8.200000
#> 12 2         1     1      6  9.235294
#> 13 1         1     1      7 13.000000
#> 14 2         1     1      7 10.400000
#> 15 1         1     1      8  9.666667
#> 16 2         1     1      8 10.955882

dat <- insight::get_data(model)
summarise_if(group_by(dat, e15relat, c161sex), is.numeric, mean) # by random group ALSO
#> # A tibble: 16 x 3
#> # Groups:   e15relat [8]
#>    e15relat c161sex neg_c_7
#>    <fct>    <fct>     <dbl>
#>  1 1        male      12.3 
#>  2 1        female    13.3 
#>  3 2        male      11.6 
#>  4 2        female    12.1 
#>  5 3        male      12.2 
#>  6 3        female    10.5 
#>  7 4        male      10.8 
#>  8 4        female    11.7 
#>  9 5        male      11.3 
#> 10 5        female    10.2 
#> 11 6        male       8.2 
#> 12 6        female     9.24
#> 13 7        male      13   
#> 14 7        female    10.4 
#> 15 8        male       9.67
#> 16 8        female    11.0

This makes more sense for plotting:

plot(me, add.data = TRUE)
#> Loading required namespace: ggplot2

plot(me, add.data = TRUE, collapse.group = "e15relat")

can also be used with residualized plots (here the random effect is residualized)

plot(me, residuals = TRUE)

plot(me, residuals = TRUE, collapse.group = "e15relat")
#> Error in FUN(X[[i]], ...): object 'group_col' not found

Well… it should…
@strengejacke why is this failing? :(

^{Created on 2021-07-06 by the reprex package (v2.0.0)}

[DominiqueMakowski]

I wonder (also for #129 ) how would a modelbased/easystats API look for that kind of stuff

[mattansb]

(This is also only for data overlay)

[bwiernik]

@DominiqueMakowski What does the API look like for estimate_prediction()? I think collapse could be an argument in those functions.

For #129, I think a different function might be best? How do we do effects plots here generally?

[DominiqueMakowski]

What does the API look like for estimate_prediction()?

For plotting? Currently data overlay is added fairly straightforwardly here (allowing for multiple data overlays like a 2D density + points):

modelbased/R/visualisation_recipe.estimate_predicted.R

Lines 177 to 186 in 5d81124

	for (i in show_data) {
	if (i %in% c("point", "points", "jitter")) {
	layers[[paste0("l", l)]] <- .visualisation_predicted_points(rawdata, x1, y, color, shape = shape, stroke = stroke, type = i, point = point)
	} else if (i %in% c("density_2d", "density_2d_filled", "density_2d_polygon", "density_2d_raster")) {
	layers[[paste0("l", l)]] <- .visualisation_predicted_density2d(rawdata, x1, y, type = i, density_2d = density_2d)
	} else {
	stop("'show_data' can only be some of 'points', 'density_2d', 'density_2d_filled', density_2d_polygon', 'density_2d_raster'. Check spelling.")
	}
	l <- l + 1
	}

Any adjustments to the raw-data should be made before that directly to the rawdata obtained currently via insight::get_data(), but how? using the model to make controlled predictions (adjusted for some variables)? by using datawizard::adjust()?

[mattansb]

Why would you need adjust() here?

A. For collapsing you just add another step to insight::get_data().
B. For #129 you add a different step to to insight::get_data() which also accounts for the data grid being used.
C. For #129 + collapsing, you do (B) and then (A).

[mattansb]

Sorry if I'm talking like I have no idea how this package works - I don't.... Waiting for that vignette...

[DominiqueMakowski]

We need to translate summarise_if(group_by_(data, "variable"), is.numeric, mean) into base R / easystats is that it ^^?

[bwiernik]

Yeah, basically.

Something like:

lapply(df, df$group, function(x) if (is.numeric(x)) mean(x, na.rm = TRUE) else NA_real_)