I want to compute predictions for a fixed effects ols regression by group with the marginaleffects package. When I generate new data for the prediction and do not specify the specific level of the fixed effects variables that I want the prediction to use, the prediction function calculates the values at am==0, gear==0 and carb==4 (see below).
Where do these values come from? What's the most sensible way to calculate average predictions in the presence of fixed effects in this example?
library(tidyverse)
library(fixest)
library(marginaleffects)
dat <- mtcars
m1 <- dat %>%
mutate(cyl = as_factor(cyl)) %>%
feols(mpg ~ hp*cyl | am + gear + carb)
etable(m1)
#> m1
#> Dependent Var.: mpg
#>
#> hp -0.1316 (0.0504)
#> cyl6 -47.87 (27.42)
#> cyl8 -21.31 (14.80)
#> hp x cyl6 0.4434 (0.2619)
#> hp x cyl8 0.1471 (0.0665)
#> Fixed-Effects: ----------------
#> am Yes
#> gear Yes
#> carb Yes
#> _______________ ________________
#> S.E.: Clustered by: am
#> Observations 32
#> R2 0.88892
#> Within R2 0.38565
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
pred1 <- predictions(m1,
newdata = datagrid(
hp = c(min(dat$hp):max(dat$hp)),
cyl = levels(dat$cyl)))
pred1
#>
#> Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 % cyl am gear carb hp
#> 9.76 8.72 1.119 0.263 1.9 -7.33 26.9 8 0 3 4 52
#> 9.78 8.60 1.136 0.256 2.0 -7.09 26.6 8 0 3 4 53
#> 9.79 8.49 1.154 0.249 2.0 -6.84 26.4 8 0 3 4 54
#> 9.81 8.37 1.172 0.241 2.1 -6.60 26.2 8 0 3 4 55
#> 9.82 8.25 1.190 0.234 2.1 -6.35 26.0 8 0 3 4 56
#> --- 274 rows omitted. See ?avg_predictions and ?print.marginaleffects ---
#> 14.09 23.89 0.590 0.555 0.8 -32.73 60.9 8 0 3 4 331
#> 14.11 24.01 0.588 0.557 0.8 -32.94 61.2 8 0 3 4 332
#> 14.12 24.12 0.586 0.558 0.8 -33.15 61.4 8 0 3 4 333
#> 14.14 24.24 0.583 0.560 0.8 -33.37 61.6 8 0 3 4 334
#> 14.16 24.36 0.581 0.561 0.8 -33.58 61.9 8 0 3 4 335
#> Columns: rowid, estimate, std.error, statistic, p.value, s.value, conf.low, conf.high, mpg, cyl, am, gear, carb, hp
Created on 2023-08-22 with reprex v2.0.2
You ask:
When I generate new data for the prediction and do not specify the specific level of the fixed effects variables that I want the prediction to use, the prediction function calculates the values at am==0, gear==0 and carb==4
(I think you meant gear=6. Must be a typo.)
When a variable is not explicitly specified by name in datagrid, the value is picked by a function which depends on the variable class. Please read ?datagrid and check out the FUN_character, FUN_numeric, and other arguments to that function.
In fixest models, the variables after the | are treated as factor variables, so the default is to take the mode. When there is more than one mode, the choice is arbitrary (I think it’s the first one to appear in the dataset, but not 100% sure). The user can of course customize this by specifying the relevant FUN_* argument.
You ask:
What’s the most sensible way to calculate average predictions in the presence of fixed effects in this example?
In my opinion, it is typically most sensible to not specify the newdata argument at all when computing average predictions by group. When you leave, newdata=NULL, the function will compute predictions for each row of the original data, and then average them by group. This default behavior gives you predicted values averaged across the empirical distribution of the covariates.
That said, this is a statistical question whose answer may vary depending on your specific research context and application. If I were you, I would consult a local statistician or open a new question on Cross Validated (not Stack Overflow) with a lot more detail about the specific research question that you hope to answer.
library(fixest)
library(marginaleffects)
dat <- mtcars
mod <- feols(mpg ~ hp*cyl | am + gear + carb, data = dat)
avg_predictions(mod, by = c("am", "cyl"))
#
# am cyl Estimate Std. Error z Pr(>|z|) S 2.5 % 97.5 %
# 0 4 23.8 12.1 1.973 0.0484 4.4 0.163 47.4
# 0 6 18.8 14.2 1.321 0.1863 2.4 -9.088 46.7
# 0 8 14.9 17.2 0.867 0.3858 1.4 -18.807 48.7
# 1 4 27.8 11.8 2.360 0.0183 5.8 4.720 51.0
# 1 6 20.4 15.4 1.327 0.1845 2.4 -9.746 50.6
# 1 8 16.5 22.1 0.747 0.4550 1.1 -26.859 59.9
#
# Columns: am, cyl, estimate, std.error, statistic, p.value, s.value, conf.low, conf.high