Predicción con factores pseudo-medios en los nuevos datos

Question

R y Stata tienen diferentes comportamientos por defecto cuando se hacen predicciones de un modelo que utiliza covariables categóricas/factoriales. Por ejemplo, si quiero predecir los resultados para ambos niveles de un factor covariable de dos niveles (en este caso, coches extranjeros y nacionales), manteniendo todos los demás valores en sus medias, el programa de Stata margins [varname], atmeans hace cosas raras con los factores, calculando el valor medio 0/1 para cada nivel:

. sysuse auto2

. reg price mpg i.foreign i.rep78

...

------------------------------------------------------------------------------
       price |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
         mpg |  -299.6068   63.34525    -4.73   0.000    -426.2322   -172.9815
             |
     foreign |
    Foreign  |   1102.334   901.7772     1.22   0.226    -700.2928    2904.961
             |
       rep78 |
       Fair  |   841.3622   2055.452     0.41   0.684    -3267.428    4950.153
    Average  |   1285.116   1901.486     0.68   0.502    -2515.901    5086.132
       Good  |   1155.571   1984.561     0.58   0.562     -2811.51    5122.652
  Excellent  |   2353.179   2130.577     1.10   0.274    -1905.784    6612.142
             |
       _cons |   10856.24   2266.757     4.79   0.000      6325.06    15387.43
------------------------------------------------------------------------------

. margins foreign, atmeans

Adjusted predictions                            Number of obs     =         69
Model VCE    : OLS

Expression   : Linear prediction, predict()
at           : mpg             =    21.28986 (mean)
               0.foreign       =    .6956522 (mean)
               1.foreign       =    .3043478 (mean)
               1.rep78         =    .0289855 (mean)
               2.rep78         =     .115942 (mean)
               3.rep78         =    .4347826 (mean)
               4.rep78         =    .2608696 (mean)
               5.rep78         =    .1594203 (mean)

------------------------------------------------------------------------------
             |            Delta-method
             |     Margin   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
     foreign |
   Domestic  |    5810.55    415.892    13.97   0.000     4979.194    6641.907
    Foreign  |   6912.884   700.8393     9.86   0.000     5511.927    8313.842
------------------------------------------------------------------------------

R, por otro lado, no puede calcular la media de un factor (ya que eso es matemáticamente imposible de todos modos), y no divide los factores en proporciones numéricas como Stata. En su lugar, al crear un nuevo marco de datos de covariables para pasar al modelo, tengo que elegir uno de los niveles del factor:

library(haven)
auto <- read_stata("http://www.stata-press.com/data/r13/auto2.dta")

model <- lm(price ~ mpg + as.factor(foreign) + as.factor(rep78), data=auto)
summary(model)
#> Coefficients:
#>                     Estimate Std. Error t value Pr(>|t|)    
#> (Intercept)         10856.24    2266.76   4.789 1.08e-05 ***
#> mpg                  -299.61      63.35  -4.730 1.34e-05 ***
#> as.factor(foreign)1  1102.33     901.78   1.222    0.226    
#> as.factor(rep78)2     841.36    2055.45   0.409    0.684    
#> as.factor(rep78)3    1285.12    1901.49   0.676    0.502    
#> as.factor(rep78)4    1155.57    1984.56   0.582    0.562    
#> as.factor(rep78)5    2353.18    2130.58   1.104    0.274    
#> ---

# Create new data with average values of all covariates for both foreign and
# domestic cars
newdata <- expand.grid(mpg = mean(auto$mpg, na.rm=TRUE),
                       foreign = c(0, 1),
                       rep78 = 3)  # One of the factor levels

# Not the same as Stata, obviously
predict(model, newdata=newdata)
#>        1        2 
#> 5760.544 6862.878

Estoy usando R para replicar un estudio que fue hecho originalmente en Stata que usó margins [varname], atmeans para generar resultados predichos a partir de un modelo con varias covariables categóricas. ¿Existe una forma de replicar el valor pseudo-medio del factor como hace Stata (descomponiendo el factor en sus niveles individuales, codificados como valores ficticios 0/1), o hay una forma más precisa de utilizar predict() con categorías "promedio" en R (que no sea elegir arbitrariamente uno de los niveles)? ¿Qué enfoque (la media de cada nivel de Stata frente a la elección de un nivel de R) es más preciso/apropiado?

Answer 1

La forma correcta de utilizar la función de Stata margins en este contexto sería añadir el comando over(foreign) opción. Stata tiene una larga tradición de tomar al usuario literalmente en lugar de tratar de interpretar lo que realmente se quiere. Así que si se añade la opción atmeans opción de pedir una evaluación en los medios, y se obtiene precisamente eso. Si quieres que se separe por un factor, añades la opción over() opción.

Estoy de acuerdo en que si se olvida el over() obtendrá resultados que no me gustaría presentar, pero es por eso que Stata le advierte dando esa salida adicional bastante "antiestática" que indica a qué valores se evalúan los márgenes.

Answer 2

Mi preferencia personal a la hora de calcular los efectos marginales en la media (EEM) sería utilizar el valor base para las variables factoriales y la media para las continuas. A veces también se denominan efectos marginales en valores representativos (MER).

Por ejemplo,

sysuse auto2, clear
reg price mpg i.foreign i.rep78, coefl
margins foreign, at((means) mpg (base) rep78)

Esto equivale a:

sum mpg if e(sample)
di _b[_cons] + _b[mpg]*r(mean)  + _b[1.foreign]
di _b[_cons] + _b[mpg]*r(mean) 

Esencialmente, usted está prediciendo una vez como si todos los coches eran extranjeros y luego otra vez como si todos los coches eran nacionales, con el mpg ajustado al en general media y rep78 ajustado a "Pobre" para todos los vehículos. Esto es diferente de la predicción para los coches extranjeros y nacionales por separado, con las covariables establecidas a las medias específicas de origen, que es lo que over() lo hace.

En esta época de big data, si no quieres escribir todos los nombres de las variables, puedes utilizar simplemente esta abreviatura:

márgenes extranjeros, at((medios) _continuos (base) _factor)

Puedes replicar lo que hace R con:

margins foreign, at((means) _continuous rep78=3)

o especificando previamente la base de esta manera:

fvset base 3 rep78
reg price mpg i.foreign i.rep78, coefl
margins foreign, at((means) _continuous (base) _factor)

En cuanto a las mejores prácticas, no creo que haya un consenso, ni siquiera dentro de un campo concreto. Personalmente, prefiero los Efectos Marginales Medios (EMA) a los EEM y los MER. Eso es lo que Stata calcula por defecto. Hay un documento más antiguo de SJ de Tamas Bartus que repasa las distintas opciones de la sección 2.