# Why individual fixed effect method is not estimating average marginal effects (panel data)?

I am using an individual fixed effect method in a panel data. I look whether the working hours changed differently between men and women following the 2008 financial crisis.

Here is a simple model in Stata and R

xtreg workhours c.age i.sex i.education i.sex##i.crisis, fe # Stata

plm(workhours ~ age + sex + education + sex*crisis, index=c("id", "year"), data=df) # R


In the model above sex and education are time fixed and hence constant over time. The crisis variable is a dummy variable (0= before 2008, 1= post 2008).

The interaction of interest sex*crisis tells me whether the gap in working hours has increased or decreased following the financial crisis. The findings show that the gap increases following the crisis, however, I am interested to see whether this increase is driven by a decrease in working hours for men or by an increase in working hours for women or whether it is driven by both cases. Thus, I would like to run the average marginal effects to see graphically the association between both variables. As there is no easy way to do it in R, I run it in Stata. After running margins crisis##sex, I get the following text: . (not estimable).

I tried to see whether it is possible to run it using a random effect model instead of FE. Both in R and Stata I can run the average marginal effects for the random effect models with no problem. This leaves me with the conclusion that the error I get for average marginal effects using FE is related to the FE itself. Maybe my question is stupid, but could someone explain to me if I am missing something here? Why the FE does not give me average marginal effects? I don't want to get into the discussion of why it is better to use RE over FE or vice versa. I am interested in understanding why this issue could occur if I use FE.

EDIT 20/10/2021:

Below you find the Stata output using the code suggested by Dimitri. Education and sex are dropped out in FE model because they are time-invariant. Then I include the average marginal effects using Random Effects and Dimitriy's code.

 xtset id year

Panel variable: id (unbalanced)
Time variable: year, 2005 to 2015, but with gaps
Delta: 1 unit

. xtreg workhours c.age i.sex i.education i.sex##i.crisis, fe
note: 2.sex omitted because of collinearity.
note: 2.education omitted because of collinearity.
note: 3.education omitted because of collinearity.

Fixed-effects (within) regression               Number of obs     =    472,757
Group variable: id                              Number of groups  =    196,088

R-squared:                                      Obs per group:
Within  = 0.0014                                         min =          1
Between = 0.0020                                         avg =        2.4
Overall = 0.0011                                         max =          8

F(3,276666)       =     128.29
corr(u_i, Xb) = -0.2447                         Prob > F          =     0.0000

----------------------------------------------------------------------------------
workhours | Coefficient  Std. err.      t    P>|t|     [95% conf. interval]
-----------------+----------------------------------------------------------------
age |  -.1713355   .0103908   -16.49   0.000    -.1917011   -.1509698
|
sex |
Female  |          0  (omitted)
|
education |
Upper-secondary  |          0  (omitted)
Tertiary  |          0  (omitted)
|
1.crisis |  -.1401835    .050655    -2.77   0.006     -.239466    -.040901
|
sex#crisis |
Female#1  |   .2906291   .0690111     4.21   0.000     .1553693     .425889
|
_cons |   46.62693   .4306217   108.28   0.000     45.78292    47.47093
-----------------+----------------------------------------------------------------
sigma_u |  9.7706582
sigma_e |  5.5553808
rho |  .75569746   (fraction of variance due to u_i)
----------------------------------------------------------------------------------
F test that all u_i=0: F(196087, 276666) = 5.87              Prob > F = 0.0000

.
end of do-file

. do "/var/folders/wb/2v3hpch94wd3_r_6q8ssng300000gp/T//SD86690.000000"

. margins sex, dydx(crisis)

Average marginal effects                               Number of obs = 472,757
Model VCE: Conventional

Expression: Linear prediction, predict()
dy/dx wrt:  1.crisis

------------------------------------------------------------------------------
|            Delta-method
|      dy/dx   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
0.crisis     |  (base outcome)
-------------+----------------------------------------------------------------
1.crisis     |
sex |
Male  |          .  (not estimable)
Female  |          .  (not estimable)
------------------------------------------------------------------------------
Note: dy/dx for factor levels is the discrete change from the base level.


Random effects models using Dimitriy's code:

     xtset id year

Panel variable: id (unbalanced)
Time variable: year, 2005 to 2015, but with gaps
Delta: 1 unit

. xtreg workhours c.age i.sex i.education i.sex##i.crisis, re

Random-effects GLS regression                   Number of obs     =    472,757
Group variable: id                              Number of groups  =    196,088

R-squared:                                      Obs per group:
Within  = 0.0003                                         min =          1
Between = 0.0885                                         avg =        2.4
Overall = 0.0764                                         max =          8

Wald chi2(6)      =   19144.98
corr(u_i, X) = 0 (assumed)                      Prob > chi2       =     0.0000

----------------------------------------------------------------------------------
workhours | Coefficient  Std. err.      z    P>|z|     [95% conf. interval]
-----------------+----------------------------------------------------------------
age |   .0209332   .0017377    12.05   0.000     .0175275    .0243389
|
sex |
Female  |  -5.770971   .0491049  -117.52   0.000    -5.867215   -5.674727
|
education |
Upper-secondary  |   .0708841      .0486     1.46   0.145    -.0243701    .1661382
Tertiary  |  -.2056841   .0543011    -3.79   0.000    -.3121124   -.0992559
|
1.crisis |  -.7906896   .0369781   -21.38   0.000    -.8631654   -.7182138
|
sex#crisis |
Female#1  |   .5714118   .0544573    10.49   0.000     .4646775    .6781461
|
_cons |   41.22829   .0851835   483.99   0.000     41.06133    41.39525
-----------------+----------------------------------------------------------------
sigma_u |  8.0439189
sigma_e |  5.5553808
rho |   .6770612   (fraction of variance due to u_i)
----------------------------------------------------------------------------------

. margins sex, dydx(crisis)

Average marginal effects                               Number of obs = 472,757
Model VCE: Conventional

Expression: Linear prediction, predict()
dy/dx wrt:  1.crisis

------------------------------------------------------------------------------
|            Delta-method
|      dy/dx   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
0.crisis     |  (base outcome)
-------------+----------------------------------------------------------------
1.crisis     |
sex |
Male  |  -.7906896   .0369781   -21.38   0.000    -.8631654   -.7182138
Female  |  -.2192778   .0402583    -5.45   0.000    -.2981825   -.1403731
------------------------------------------------------------------------------
Note: dy/dx for factor levels is the discrete change from the base level.

.


It’s not clear if the sex variable corresponds to being male or female. Let’s say it’s female. In that case, the effect for men is the crisis coefficient and the effect for women is the sum of crisis and the sex-crisis interaction coefficients. If sex corresponds to male, it's reversed.

You are asking Stata to calculate predictions for the four combinations of sex and crisis, which it cannot do since the effect of sex is not separately identified from the FEs.

You should use margins sex, dydx(crisis), which will give you the AME for each group.

Here's a reproducible example with a continuous-binary interaction:

. webuse pig, clear
(Longitudinal analysis of pig weights)

. xtset id week

Panel variable: id (strongly balanced)
Time variable: week, 1 to 9
Delta: 1 unit

. gen group = mod(id,2) + 1

. replace weight = weight*(1 + .029*week/9) if group == 2

. xtreg weight c.week c.week#i.group, fe vce(cluster id)

Fixed-effects (within) regression               Number of obs     =        432
Group variable: id                              Number of groups  =         48

R-squared:                                      Obs per group:
Within  = 0.9852                                         min =          9
Between = 0.0358                                         avg =        9.0
Overall = 0.9321                                         max =          9

F(2,47)           =    2294.34
corr(u_i, Xb) = -0.0006                         Prob > F          =     0.0000

(Std. err. adjusted for 48 clusters in id)
------------------------------------------------------------------------------
|               Robust
weight | Coefficient  std. err.      t    P>|t|     [95% conf. interval]
-------------+----------------------------------------------------------------
week |   6.179861   .1278409    48.34   0.000     5.922678    6.437044
|
group#c.week |
2  |   .3240612   .1874244     1.73   0.090     -.052988    .7011104
|
_cons |   19.17099   .4685609    40.91   0.000     18.22837    20.11361
-------------+----------------------------------------------------------------
sigma_u |  3.9891874
sigma_e |  2.1314439
rho |  .77791824   (fraction of variance due to u_i)
------------------------------------------------------------------------------

. margins group, dydx(week)

Average marginal effects                                   Number of obs = 432
Model VCE: Robust

Expression: Linear prediction, predict()
dy/dx wrt:  week

------------------------------------------------------------------------------
|            Delta-method
|      dy/dx   std. err.      z    P>|z|     [95% conf. interval]
-------------+----------------------------------------------------------------
week         |
group |
1  |   6.179861   .1278409    48.34   0.000     5.929297    6.430425
2  |   6.503922   .1370569    47.45   0.000     6.235296    6.772549
------------------------------------------------------------------------------


You can follow that up with marginsplot if you want to see these effects graphed.

New Code:

webuse pig, clear
xtset id week
gen group = mod(id,2)
replace weight = weight*(1 + .029*week/9) if group == 2
gen dummy = week >5
xtreg weight i.dummy i01.dummy#i.group, fe vce(cluster id)
margins group, dydx(dummy)
xtreg weight c.dummy c.dummy#i.group, fe vce(cluster id)
margins group, dydx(dummy)

• Thanks a lot for the response, i use the method suggested but I still get the same response. In order for things to be more clear, I include now both the average marginal effects of FE and RE using your suggested code in my original post. You can check them there
– Jack
Oct 20, 2021 at 7:04
• I am not sure exactly what is happening with your data. Sometimes you need to tweak/trick the factor variable notation as I did above. If that does not work, then try simplifying your specification by removing the extraneous variables. I might also add how your variables of interest are coded. Oct 20, 2021 at 7:44