Python: медленная регрессия по сравнению со Stata (манекены с фиксированным эффектом)

Я пытаюсь запустить регрессию в Python, но это занимает время и перестает работать. В Stata это работает и занимает всего несколько секунд.

Это связано с категориальным столбцом, включающим групповые фиксированные эффекты. Без переменной производительность Stata и Python примерно одинакова: на 200000 наблюдений требуется около 1 секунды:

Код Стата

reg income height Number_children

Код Python

model = smf.ols(income ~ height + Number_children, data=humans).fit()

Добавляя манекен, я меняю код Stata на areg:

areg income height Number_children, absorb(Village)

Это занимает всего 1-2 секунды дольше, чем без манекена.

В Python:

model = smf.ols(income ~ height + Number_children + Village, data=humans).fit()

куда:

Name: Village, dtype: category
Categories (3678, object):

Останавливаю регресс через 2 минуты ожидания. Есть ли идеи, как запустить код и увеличить скорость почти до такой же скорости, как у Stata? Проблема больше вызвана переменной или командой регрессии?

  • РЕДАКТИРОВАТЬ:

Основываясь на ответе Дмитрия, я попробовал это для всех переменных:

e.g.: 

humans["income_gr_m"]= humans["income"].groupby(humans['Village']).mean()
humans["income_star"] = humans["income"] - humans["income_gr_m"] + humans["income"].mean()

Однако это также заставляет Python работать как минимум 2 минуты (я снова остановился). Или трансформацию надо проводить по-другому? Спасибо


python regression dummy-variable stata
person user27074    schedule 29.03.2018    source источник

Ответы (1)


arrow_upward
3
arrow_downward

areg на самом деле не инвертирует эту матрицу с 3677 индикаторами деревень, как вы делаете с Python. Он преобразует данные таким образом, чтобы исключить необходимость в этом, поэтому это будет намного быстрее. Это также причина того, что константа из regress с деревенскими манекенами не будет соответствовать константе из areg, хотя коэффициенты наклона должны быть такими же, если вы дождетесь завершения Python.

Вот что делает areg для вычисления коэффициентов с regress. Стандартные ошибки будут слишком большими, так как я не выполняю регулировку степени свободы для 5 поглощенных эффектов, но я сделаю корректировку вручную в цикле, указанном ниже, путем умножения SE:

. sysuse auto, clear
(1978 Automobile Data)

. drop if missing(rep78)
(5 observations deleted)

. /* (1) transform the data by subtracting the group specific mean and */
. /* adding the grand/overall mean back in for outcome and regressors */
. foreach var of varlist price weight length foreign {
  2.         bys rep78: egen group_mean = mean(`var')
  3.         qui sum `var'
  4.         gen double `var'_star = `var' - group_mean + r(mean)
  5.         drop group_mean
  6. }

. /* (2) Fit the model on transformed data */
. regress price_star weight_star length_star foreign_star

      Source |       SS           df       MS      Number of obs   =        69
-------------+----------------------------------   F(3, 65)        =     26.99
       Model |   315296838         3   105098946   Prob > F        =    0.0000
    Residual |   253139578        65  3894455.05   R-squared       =    0.5547
-------------+----------------------------------   Adj R-squared   =    0.5341
       Total |   568436416        68  8359359.06   Root MSE        =    1973.4

------------------------------------------------------------------------------
  price_star |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
 weight_star |    6.15521   1.008605     6.10   0.000     4.140885    8.169534
 length_star |  -100.9268   33.82508    -2.98   0.004    -168.4801   -33.37341
foreign_star |   3394.052    782.454     4.34   0.000     1831.383     4956.72
       _cons |   5453.782   3829.487     1.42   0.159    -2194.232     13101.8
------------------------------------------------------------------------------

. /* (3) Adjust the SEs for DoF */
. foreach coef in weight_star length_star foreign_star _cons {
  2.         di "Adjusted SE for `coef': " %9.8gc _se[`coef']*sqrt(65/61)
  3. }
Adjusted SE for weight_star:  1.041149
Adjusted SE for length_star:  34.91649
Adjusted SE for foreign_star:  807.7009
Adjusted SE for _cons:   3953.05

. /* (4) Make sure areg gives the same output */
. areg price weight length foreign, absorb(rep78)

Linear regression, absorbing indicators         Number of obs     =         69
                                                F(   3,     61)   =      25.33
                                                Prob > F          =     0.0000
                                                R-squared         =     0.5611
                                                Adj R-squared     =     0.5108
                                                Root MSE          =  2037.1129

------------------------------------------------------------------------------
       price |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      weight |    6.15521   1.041149     5.91   0.000     4.073303    8.237116
      length |  -100.9268   34.91649    -2.89   0.005    -170.7466   -31.10692
     foreign |   3394.052   807.7009     4.20   0.000     1778.954    5009.149
       _cons |   5453.782    3953.05     1.38   0.173    -2450.831    13358.39
-------------+----------------------------------------------------------------
       rep78 |          F(4, 61) =      0.261   0.902           (5 categories)

Код статистики:

cls
sysuse auto, clear
drop if missing(rep78)
/* (1) transform the data by subtracting the group specific mean and */
/* adding the grand/overall mean back in for outcome and regressors */
foreach var of varlist price weight length foreign {
    bys rep78: egen group_mean = mean(`var')
    qui sum `var'
    gen double `var'_star = `var' - group_mean + r(mean)
    drop group_mean
}
/* (2) Fit the model on transformed data */
regress price_star weight_star length_star foreign_star
/* (3) Adjust the SEs for DoF */
foreach coef in weight_star length_star foreign_star _cons {
    di "Adjusted SE for `coef': " %9.8gc _se[`coef']*sqrt(65/61)
}
/* (4) Make sure areg gives the same output */
areg price weight length foreign, absorb(rep78)
person Dimitriy V. Masterov    schedule 29.03.2018
comment
Это здорово, но у вас есть идея, как быстро выполнить аналогичный расчет в Python? Преобразование всех предполагаемых переменных замедлит весь процесс. - person user27074; 30.03.2018
comment
Вычисление нескольких средних и сложение / вычитание не должно занимать так много времени, так как это можно очень легко распараллелить. К сожалению, я новичок в Python, поэтому я не могу с готовностью помочь, кроме описания алгоритма. Вы можете выполнить поиск пакетов регрессии с фиксированными эффектами большой размерности в Python. - person Dimitriy V. Masterov; 30.03.2018