Конечно, как и в случае с любой другой командой Stata, вы можете использовать интерфейс «укажи и щелкни» вместо ввода команд.

Байесовские модели панельных данных предназначены не только для непрерывных исходов. Вы можете так же легко ввести команду для бинарных исходов

. bayes: xtprobit y x1 x2

или для подсчета результатов

. bayes: xtpoisson y x1 x2

Или используйте любую из восьми моделей панельных данных, которые поддерживают префикс bayes, включая новую модель мультиномиального логита для панельных данных.

Давайте посмотрим, как это работает

• Оценка
• Оценка
• Пользовательские приоритеты
• Постериорные распределения панельных или групповых эффектов
• Байесовские прогнозы
• Последующие проверки предсказания — соответствие модели
• Уборка

Оценка

Рассмотрим подмножество данных Национального продольного обследования молодых женщин в возрасте от 14 до 24 лет в 1968 году, проживающих на Юге. Мы будем моделировать логарифм заработной платы как функцию образования индивида, grade; его опыта работы, ttl_exp, который входит в модель квадратично; и того, живет ли он в стандартном мегаполисе, not_smsa. Для каждого индивида имеется несколько наблюдений, обозначаемых id.

Мы используем байесовскую линейную модель панельных данных для учета индивидуальных эффектов. Мы также хотим вычислить байесовские прогнозы логарифмической заработной платы и сравнить индивидуальные (панельные) эффекты.

. webuse nlswork6
(Subsample of 1986 National Longitudinal Survey of Young Women)

Для экономии времени мы запускаем выборку Марковской цепи Монте-Карло (MCMC) в 1 000, вместо 10 000 по умолчанию. Для воспроизводимости мы также указываем семя случайного числа.

. bayes, mcmcsize(1000) rseed(17): xtreg ln_wage grade c.ttl_exp##c.ttl_exp i.not_smsa
note: Gibbs sampling is used for regression coefficients and variance
      components.

Burn-in 2500 aaaaaaaaa1000aaaaaaaaa2000aaaaa done
Simulation 1000 .........1000 done

Model summary

Likelihood: 
  ln_wage ~ normal(xb_ln_wage,{sigma2})

Priors: 
                {ln_wage:grade} ~ normal(0,10000)                          (1)
              {ln_wage:ttl_exp} ~ normal(0,10000)                          (1)
  {ln_wage:c.ttl_exp#c.ttl_exp} ~ normal(0,10000)                          (1)
           {ln_wage:1.not_smsa} ~ normal(0,10000)                          (1)
                {ln_wage:_cons} ~ normal(0,10000)                          (1)
                        {U[id]} ~ normal(0,{var_U})                        (1)
                       {sigma2} ~ igamma(0.01,0.01)

Hyperprior: 
  {var_U} ~ igamma(0.01,0.01)

(1) Parameters are elements of the linear form xb_ln_wage.

Bayesian RE normal regression                    MCMC iterations  =      3,500
Metropolis–Hastings and Gibbs sampling           Burn-in          =      2,500
                                                 MCMC sample size =      1,000
Group variable: id                               Number of groups =        831
                                                 Obs per group:
                                                              min =          1
                                                              avg =        1.4
                                                              max =          5
                                                 Number of obs    =      1,174
                                                 Acceptance rate  =       .796
                                                 Efficiency:  min =     .02411
                                                              avg =     .07145
Log marginal-likelihood                                       max =      .1676
 


         
 
                                                               
         
 
                                               Equal-tailed
         
 
      Mean   Std. dev.     MSCE     Median  [95% cred. interval]

ln_wage
 
                                                                     
          grade
 
  .0696248   .0052371    .00082    .069638   .0591421     .07999
        ttl_exp
 
  .0404646   .0070371   .000544   .0404901    .026454    .053855
               
 
                                                                
     c.ttl_exp#
 
                                                               
      c.ttl_exp
 
 -.0005053   .0003978   .000039  -.0005132  -.0012609   .0003009
               
 
                                                                
     1.not_smsa
 
 -.1502314   .0251494   .002755   -.150213   -.200979  -.1014476
          _cons
 
  .5646936   .0658804     .0099   .5680976   .4377888   .6957415

          var_U
 
  .0801363   .0073138   .001258   .0795163   .0674855   .0959143
         sigma2
 
  .0673308   .0046724   .000952   .0671238   .0590742   .0771401

Апостериорное среднее для коэффициента grade положительно, его величина составляет 7 процентов. Теория предполагает, что заработная плата увеличивается с ростом опыта, но со временем этот рост уменьшается. Это означает, что коэффициент ttl_exp должен быть положительным, а коэффициент c.ttl_exp#c.ttl_exp должен быть отрицательным. Мы наблюдаем это в наших данных. Наконец, проживание за пределами крупных городских центров отрицательно влияет на заработную плату.

Байесовское моделирование требует, чтобы вы указали приоритеты для всех параметров модели. Как и любая другая команда bayes, bayes: xtreg для удобства предоставляет приоритеты по умолчанию. Вы должны просмотреть спецификации приоритетов и указать свои собственные.

Сходимость MCMC

bayes: xtreg использует MCMC для получения результатов. Перед дальнейшим анализом необходимо проверить его сходимость. Вы можете сделать это графически или вычислить статистику сходимости Гельмана-Рубина, используя несколько цепочек.

Давайте оценим сходимость MCMC визуально, скажем, для коэффициента grade.

. bayesgraph diagnostics {ln_wage:grade}

На графике трассировки нет явной тенденции. Автокорреляция уменьшается со временем. У нас нет причин подозревать отсутствие сходимости для grade, но сходимость должна быть проверена для всех параметров модели.

Давайте сохраним наши результаты MCMC и оценки для дальнейшего использования.

. bayes, saving(bxtregsim)
note: file bxtregsim.dta saved.

. estimates store bxtreg

Пользовательские приоритеты

В байесовских моделях вы можете захотеть включить свои собственные прайоры. Такие прайоры часто берутся из исторических данных. Например, на основании предыдущих исследований можно предположить, что коэффициент для grade колеблется между 0 и 25. Поэтому мы можем использовать для него равномерный приоритет (0,25). Мы можем использовать опцию prior() в bayes для задания пользовательских приоритетов.

. bayes, mcmcsize(1000) rseed(17) prior({ln_wage:grade}, uniform(0,25)):
> xtreg ln_wage grade c.ttl_exp##c.ttl_exp i.not_smsa
note: Gibbs sampling is used for variance components.

  
Burn-in 2500 aaaaaaaaa1000aaaaaaaaa2000aaaaa done
Simulation 1000 .........1000 done


Model summary

Bayesian RE normal regression                    MCMC iterations  =      3,500
Metropolis–Hastings and Gibbs sampling           Burn-in          =      2,500
                                                 MCMC sample size =      1,000
Group variable: id                               Number of groups =        831
                                                 Obs per group:
                                                              min =          1
                                                              avg =        1.4
                                                              max =          5
                                                 Number of obs    =      1,174
                                                 Acceptance rate  =      .6181
                                                 Efficiency:  min =    .009861
                                                              avg =     .01259
Log marginal-likelihood                                       max =     .01751



         
 
                                                               
         
 
                                               Equal-tailed
         
 
      Mean   Std. dev.     MSCE     Median  [95% cred. interval]

ln_wage
 
                                                                     
          grade
 
  .0691558   .0063225   .002013   .0692614   .0567288   .0819581
        ttl_exp
 
  .0389882   .0080916   .002355   .0390646   .0228502   .0530886
               
 
                                                                
     c.ttl_exp#
 
                                                               
      c.ttl_exp
 
 -.0004548   .0004804   .000153  -.0004778  -.0013733   .0005973
               
 
                                                                
     1.not_smsa
 
 -.1452854   .0241483   .006649  -.1450175  -.1934294  -.0947368
          _cons
 
  .5751218   .0727991   .022234   .5795395   .4176046   .7018132

          var_U
 
  .0785815   .0082846    .00198   .0784331   .0615956   .0956798
         sigma2
 
  .0688682   .0059815   .001537   .0684258   .0583475   .0823034

Наш пользовательский приор не сильно изменил результаты.

Постериорные распределения панельных или групповых эффектов

Мы можем быть заинтересованы в том, чтобы сделать вывод об индивидуальных (панельных) эффектах. Байесовское моделирование обеспечивает естественный способ сделать это. В отличие от классических моделей панельных данных со случайными эффектами, байесовские модели панельных данных оценивают случайные (панельные) эффекты вместе со всеми другими параметрами модели. Таким образом, каждый случайный эффект представлен целой выборкой MCMC из его апостериорного распределения. Эта выборка может быть использована для выводов, таких как сравнение панельных или групповых эффектов. (Подумайте о сравнении показателей различных компаний или больниц).

Вернемся к нашей исходной модели, в которой использовались неинформативные прайоры по умолчанию.

. estimates restore bxtreg
(results bxtreg are active now)

Построим апостериорные распределения для первых девяти индивидуальных эффектов.

. bayesgraph histogram {U[1/9]}, byparm normal

Приведенные выше индивидуальные эффекты представляют собой сдвиги или смещения по сравнению со средней логарифмической заработной платой. Заработная плата отдельных лиц определенно различается. Например, мы видим, что зарплата человека 8 выше, чем у человека 5. Таким образом, между зарплатами отдельных людей все еще существуют различия, которые не учитываются предикторами модели.

Байесовские прогнозы

В рамках байесовской системы можно вычислять прогнозы и их неопределенности, не делая никаких асимптотических предположений. Это возможно потому, что полученные прогнозы являются выборками из «точного» апостериорного распределения прогноза новых данных с учетом наблюдаемых данных. При этом не предполагается, что апостериорное распределение предсказаний является асимптотически нормальным.

Вычислим апостериорные значения прогнозируемой логарифмической заработной платы вместе с их 95% доверительным интервалом

. bayespredict pmean, mean 

Computing predictions ...

. bayespredict cri_l cri_u, cri

Computing predictions ...

Перечислим результаты для первых 10 наблюдений.

. list ln_wage pmean cri_l cri_u in 1/10


     
  ln_wage      pmean       cri_l       cri_u 
  1.
 1.543923    1.55197    .8831815    2.184215
  2.
 1.815738   1.645933    1.006488    2.284888 
  3.
 1.870532   1.718343    1.112855    2.360129 
  4.
 2.340405    2.16819    1.523068    2.747879 
  5.
 1.545327   1.629953    1.092535    2.221188 
  6.
 1.594335   1.687503    1.089151    2.268714 
  7.
 1.848619   1.849651    1.266758    2.408907 
  8.
 1.757637    1.80307    1.140366    2.467546 
  9.
 2.346746   2.224904    1.529079    2.833333 
 10.
 3.539868   2.826903    2.270285    3.511156 

Предсказанные апостериорные средние близки к наблюдаемым значениям, за исключением, возможно, наблюдения 10. В следующем разделе мы покажем, как более формально проверить, хорошо ли подходит модель. Но для этого мы должны сначала сохранить все предсказания MCMC.

. bayespredict {_ysim1}, saving(bxtregpred)

Computing predictions ...

file bxtregpred.dta saved.
file bxtregpred.ster saved.

bxtregpred.dta содержит MCMC-выборку значений для каждого наблюдения. Постериорные средние (pmean), которые мы предсказали ранее, являются для каждого наблюдения средним значением по выборке MCMC.

Наборы данных MCMC-предсказаний обычно большие, поэтому их можно сохранять только при необходимости; подробности см. в [BAYES] bayespredict.

Последующие проверки предсказания — соответствие модели

Еще одним преимуществом байесовского анализа моделей панельных данных являются формализованные апостериорные прогностические проверки для проверки пригодности модели, которые недоступны при классическом частотном анализе.

Мы можем использовать апостериорные выборки для предсказанного исхода, выборку предсказания MCMC, для выполнения проверки пригодности модели. Например, давайте сравним минимальную и максимальную статистику из выборки предсказаний MCMC с теми, которые наблюдаются в данных. Эти статистики описывают хвосты распределения данных. Вы можете использовать любые другие статистики вместо (или в дополнение к) максимуму и минимуму.
bayesstats ppvalues выполняет проверку апостериорного предсказания, вычисляя так называемые апостериорные предсказательные p-значения (PPP). ППП описывают, как часто статистика из выборки предсказаний MCMC была такой же или более экстремальной, чем статистика в наблюдаемой выборке; см. [BAYES] bayesstats ppvalues.

. bayesstats ppvalues (pmin:@min({_ysim1})) (pmax:@max({_ysim1})) using bxtregpred

Posterior predictive summary   MCMC sample size =     1,000



              T 
 
      Mean   Std. dev.  E(T_obs)  P(T>=T_obs)

           pmin 
 
  .0299837   .1652143   .0176546         .565
           pmax 
 
  3.515637   .2291669    3.85025         .085

Note: P(T>=T_obs) close to 0 or 1 indicates lack of fit.

ППС для минимума, .565, не близок к 0 или 1, поэтому модель хорошо подходит для этой статистики. Но ППС для максимума, .085, близок к 0, что говорит о плохой подгонке для этой статистики.

Если мы изучим данные, то обнаружим кластер наблюдений со значениями между 3,5 и 4, которые выглядят как выбросы. Возможно, нам придется пересмотреть нашу модель и решить проблему выбросов.

Уборка

После анализа мы можем удалить файлы, созданные bayes: var и bayespredict, которые больше не нужны.

. erase bxtregsim.dta
. erase bxtregpred.dta
. erase bxtregpred.ster