Дата занятия: 19 сентября

Практика в R

LSDV-модель

Массив, с которым мы сейчас поработаем, включает данные по 27 посткоммунистическим странам. Каждый период – пятилетка начиная с 1989 г. Представлены показатели политического режима (\(fh\_polity\)) и государственной состоятельности (\(state\_capacity\)).

panel<-read_dta("RAPDC_lab1.dta")

LSDV <- lm(fh_polity~state_capacity + country, data = panel)
summary(LSDV)
## 
## Call:
## lm(formula = fh_polity ~ state_capacity + country, data = panel)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.2753 -0.4968  0.1919  0.5423  3.9501 
## 
## Coefficients:
##                          Estimate Std. Error t value Pr(>|t|)    
## (Intercept)             2.7102295  0.8322841   3.256 0.001512 ** 
## state_capacity          0.7845941  0.1343982   5.838 5.74e-08 ***
## countryArmenia          0.0003584  0.6999010   0.001 0.999592    
## countryAzerbaijan      -2.8211369  0.7265017  -3.883 0.000179 ***
## countryBelarus         -3.6416852  0.6855667  -5.312 5.95e-07 ***
## countryBulgaria         1.4102365  0.6856351   2.057 0.042135 *  
## countryCroatia         -1.2813502  0.7262585  -1.764 0.080532 .  
## countryCzech Republic   0.7278715  0.7478433   0.973 0.332602    
## countryEstonia          0.9959629  0.7095131   1.404 0.163296    
## countryGeorgia          1.3685292  0.7465533   1.833 0.069563 .  
## countryHungary          1.3291870  0.7347050   1.809 0.073237 .  
## countryKazakhstan      -3.4400550  0.6898298  -4.987 2.38e-06 ***
## countryKyrgyz Republic -1.2676060  0.7121942  -1.780 0.077938 .  
## countryLatvia           1.7385071  0.6844088   2.540 0.012516 *  
## countryLithuania        2.1021607  0.6912702   3.041 0.002965 ** 
## countryMacedonia       -0.0471028  0.6954089  -0.068 0.946124    
## countryMoldova          1.2013274  0.6972578   1.723 0.087790 .  
## countryMongolia         1.8325965  0.6835142   2.681 0.008500 ** 
## countryPoland           1.4646948  0.7143880   2.050 0.042782 *  
## countryRomania          0.6633450  0.6886991   0.963 0.337626    
## countryRussia          -0.4983025  0.6928297  -0.719 0.473568    
## countrySerbia          -0.3827508  0.6941581  -0.551 0.582517    
## countrySlovak Republic  0.8342225  0.7269053   1.148 0.253678    
## countrySlovenia        -0.0353723  0.8449212  -0.042 0.966685    
## countryTajikistan      -1.8143519  0.7798914  -2.326 0.021880 *  
## countryTurkmenistan    -5.1161202  0.7079871  -7.226 7.69e-11 ***
## countryUkraine          1.1699991  0.6940497   1.686 0.094757 .  
## countryUzbekistan      -4.3890763  0.7400814  -5.931 3.76e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.08 on 107 degrees of freedom
## Multiple R-squared:  0.8961, Adjusted R-squared:  0.8699 
## F-statistic: 34.18 on 27 and 107 DF,  p-value: < 2.2e-16

Оценим LSDV-модель. Предсказуемо, получили громоздкую малопрезентабельную выдачу.

(Intercept) – по умолчанию R в качестве базовой категории берет первую страну. Таким образом, 2.71 показывает среднее значение уровня демократичности в Албании при стартовом – нулевом значении государственной состоятельности.

Возьмем одну из значимых оценок коэффициентов для дамми-переменных. К примеру, оценка для Азербайджана демонстрирует, что при прочих равных условиях значение уровня демократичности меньше в среднем на 2.82.

В то время как FE-модель предполагает разные стартовые условия, взаимосвязь задается одинаковой во всех странах. Оценка коэффициента при \(state\_capacity\) демонстрирует, что в среднем при прочих равных при увеличении государственной состоятельности на единицу измерения уровень демократичности увеличивается на 0.78, что справедливо для всей выборки.

Есть также возможность изменить базовую категорию. К примеру, мы хотим сравнить результаты с Россией. Обратите внимание, что для изменения базовой категории группирующая переменная \(country\) должна быть факторной.

panel$countryfa1 <- relevel(factor(panel$country), ref ="Russia")
LSDV_Rus <- lm(fh_polity~state_capacity + countryfa1, data = panel)
summary(LSDV_Rus)
## 
## Call:
## lm(formula = fh_polity ~ state_capacity + countryfa1, data = panel)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.2753 -0.4968  0.1919  0.5423  3.9501 
## 
## Coefficients:
##                           Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                 2.2119     0.7416   2.983 0.003540 ** 
## state_capacity              0.7846     0.1344   5.838 5.74e-08 ***
## countryfa1Albania           0.4983     0.6928   0.719 0.473568    
## countryfa1Armenia           0.4987     0.6842   0.729 0.467704    
## countryfa1Azerbaijan       -2.3228     0.6958  -3.338 0.001161 ** 
## countryfa1Belarus          -3.1434     0.6857  -4.584 1.24e-05 ***
## countryfa1Bulgaria          1.9085     0.7047   2.708 0.007875 ** 
## countryfa1Croatia          -0.7830     0.7729  -1.013 0.313284    
## countryfa1Czech Republic    1.2262     0.8016   1.530 0.129032    
## countryfa1Estonia           1.4943     0.7488   1.996 0.048523 *  
## countryfa1Georgia           1.8668     0.7080   2.637 0.009620 ** 
## countryfa1Hungary           1.8275     0.7844   2.330 0.021687 *  
## countryfa1Kazakhstan       -2.9418     0.6835  -4.304 3.73e-05 ***
## countryfa1Kyrgyz Republic  -0.7693     0.6886  -1.117 0.266416    
## countryfa1Latvia            2.2368     0.7007   3.192 0.001854 ** 
## countryfa1Lithuania         2.6005     0.7179   3.623 0.000448 ***
## countryfa1Macedonia         0.4512     0.7257   0.622 0.535446    
## countryfa1Moldova           1.6996     0.6836   2.486 0.014459 *  
## countryfa1Mongolia          2.3309     0.6897   3.379 0.001014 ** 
## countryfa1Poland            1.9630     0.7561   2.596 0.010746 *  
## countryfa1Romania           1.1616     0.7124   1.631 0.105914    
## countryfa1Serbia            0.1156     0.6833   0.169 0.866021    
## countryfa1Slovak Republic   1.3325     0.7738   1.722 0.087945 .  
## countryfa1Slovenia          0.4629     0.9173   0.505 0.614849    
## countryfa1Tajikistan       -1.3160     0.7314  -1.799 0.074774 .  
## countryfa1Turkmenistan     -4.6178     0.6868  -6.723 8.98e-10 ***
## countryfa1Ukraine           1.6683     0.6832   2.442 0.016256 *  
## countryfa1Uzbekistan       -3.8908     0.7039  -5.527 2.31e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.08 on 107 degrees of freedom
## Multiple R-squared:  0.8961, Adjusted R-squared:  0.8699 
## F-statistic: 34.18 on 27 and 107 DF,  p-value: < 2.2e-16

Предсказуемо, коэффициент для дамми при России исчез из выдачи, а для Албании, напротив, появился. Интерпретация схожая с интерпретацией предыдущей модели, однако сравниваем теперь с Россией.

Визуализация

Проиллюстрируем, что мы “заложили” в модель разные стартовые условия и одинаковую взаимосвязь государственной состоятельности и политического режима. Для этого предварительно сохраним предсказанные значения, полученные по итогам оценивания модели. Данный график нельзя интерпретировать как подтверждение того, что взаимосвязь одинаковая, это только иллюстрация того, на каких допущениях построена спецификация модели.

panel$ypredicted <- LSDV$fitted
colors <- colors(distinct = TRUE)
mypalette <- sample(colors, length(unique(panel$country)))
ggplot(panel, aes(x = state_capacity, y = ypredicted, color = country))+geom_smooth(method = lm)  + scale_color_manual(values = mypalette) + ylab("Predicted Level of Democracy") + xlab("State Capacity") 

Внутригрупповое преобразование

Как мы уже видели, LSDV-модель довольно перегружена. Переоценим модель с использованием внутригруппового преобразования. Используем пакет plm:

  • в index укажем группирующие переменные, чтобы R правильно считал структуру данных (первая переменная \(country\) отвечает за пространственное измерение, вторая переменная \(period\) – соответственно, за временную перспективу),

  • effect = “individual”: FE включаются на пространственные единицы, в R по умолчанию стоит эта опция, так что можно ее в явном виде не прописывать, если Вас интересует модель с FE на пространственные единицы

  • model = “within”: указывает на то, что модель будет оценена с использованием внутригруппового преобразования

fe <- plm(fh_polity~state_capacity, data = panel, index=c("country", "period"), effect = "individual", model="within")
summary(fe)
## Oneway (individual) effect Within Model
## 
## Call:
## plm(formula = fh_polity ~ state_capacity, data = panel, effect = "individual", 
##     model = "within", index = c("country", "period"))
## 
## Balanced Panel: n = 27, T = 5, N = 135
## 
## Residuals:
##     Min.  1st Qu.   Median  3rd Qu.     Max. 
## -3.27528 -0.49682  0.19188  0.54226  3.95012 
## 
## Coefficients:
##                Estimate Std. Error t-value  Pr(>|t|)    
## state_capacity  0.78459    0.13440  5.8378 5.736e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Total Sum of Squares:    164.63
## Residual Sum of Squares: 124.86
## R-Squared:      0.24157
## Adj. R-Squared: 0.050186
## F-statistic: 34.0803 on 1 and 107 DF, p-value: 5.7355e-08

Можно выгрузить оценки констант и соответствующие стандартные ошибки для каждой страны следующим образом:

summary(fixef(fe))
##                 Estimate Std. Error t-value  Pr(>|t|)    
## Albania          2.71023    0.83228  3.2564 0.0015119 ** 
## Armenia          2.71059    0.71402  3.7962 0.0002440 ***
## Azerbaijan      -0.11091    0.64720 -0.1714 0.8642613    
## Belarus         -0.93146    0.78672 -1.1840 0.2390442    
## Bulgaria         4.12047    0.87985  4.6832 8.336e-06 ***
## Croatia          1.42888    1.04272  1.3703 0.1734480    
## Czech Republic   3.43810    1.09427  3.1419 0.0021713 ** 
## Estonia          3.70619    0.99437  3.7272 0.0003113 ***
## Georgia          4.07876    0.61266  6.6575 1.233e-09 ***
## Hungary          4.03942    1.06394  3.7966 0.0002436 ***
## Kazakhstan      -0.72983    0.75667 -0.9645 0.3369596    
## Kyrgyz Republic  1.44262    0.67864  2.1258 0.0358227 *  
## Latvia           4.44874    0.86559  5.1395 1.249e-06 ***
## Lithuania        4.81239    0.92003  5.2307 8.452e-07 ***
## Macedonia        2.66313    0.94091  2.8304 0.0055541 ** 
## Moldova          3.91156    0.72341  5.4071 3.928e-07 ***
## Mongolia         4.54283    0.81573  5.5690 1.922e-07 ***
## Poland           4.17492    1.00956  4.1354 7.071e-05 ***
## Romania          3.37357    0.90436  3.7303 0.0003079 ***
## Russia           2.21193    0.74161  2.9826 0.0035399 ** 
## Serbia           2.32748    0.73575  3.1634 0.0020299 ** 
## Slovak Republic  3.54445    1.04440  3.3938 0.0009677 ***
## Slovenia         2.67486    1.27028  2.1057 0.0375686 *  
## Tajikistan       0.89588    0.56953  1.5730 0.1186659    
## Turkmenistan    -2.40589    0.68957 -3.4890 0.0007052 ***
## Ukraine          3.88023    0.73621  5.2706 7.117e-07 ***
## Uzbekistan      -1.67885    0.62290 -2.6952 0.0081709 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Значения теперь уже интерпретируются не в терминах отклонения от базовой категории, они получены посредством суммирования оценки константы для базовой категории и оценки коэффициента при дамми для i-ой страны. Так для Азербайджана мы видели, что отклонение в уровне демократичности от Албании составляет примерно \(-2.82\). В результате суммирования с оценкой константы для Албании получили \(-0.11\), что показывает, чему в Азербайджане был бы равен показатель политического режима при стартовом значении государственной состоятельности. Общая оценка константы статистически незначима.

Процедура взвешивания

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

  1. Реализуем первый шаг: рассчитаем внутригрупповые дисперсии данной объясняющей переменной по подгруппам – странам.
var <- summarize(group_by(panel, country), var(state_capacity))
var
  1. На втором шаге оценим отдельные регрессионные модели \(fh\_polity\) на \(state\_capacity\) на каждой подвыборке – стране.
  • Точка в синтаксисе кажется непривычной, дело в том, что ее использование требует функция do. Точка – как обращение к данным, определенным на предыдущем шаге (до pipe %>% мы сгруппировали данные, каждый последующий раз будем оценивать модель на данных i-ой страны).

  • [2] – в квадратных скобках стоит номер параметра в модели, который нас интересует и который мы выгружаем в датафрейм. В данном случае 2 отвечает за коэффициент при предикторе (в отличие от 1 – константы).

a <- group_by(panel, country) %>%
  do(data.frame(beta = coef(lm(fh_polity ~ state_capacity, data = .))[2]))
a
  1. Третий шаг – непосредственно рассчитаем оценку коэффициента при предикторе посредством ранее определенной процедуры взвешивания.
m <- merge(a, var, by ="country")
m
m$coef <- m$beta*(m$`var(state_capacity)`/sum(m$`var(state_capacity)`))
sum(m$coef)
## [1] 0.7845941

И вуаля: получили то же значение оценки коэффициента при \(state\_capacity\), как и в FE-модели.

FE-модель VS pooled-модель

Проверим, а нужна ли нам вообще такая модель, которая учитывает панельную структуру, или можно оценить обычную регрессионную модель без поправок на неоднородность? Оценим для этого модель pooled. Конечно, можно было бы использовать хорошо знакомую нам функцию lm, однако так как мы сравниваем ее с FE-моделью из класса моделей для анализа панельных данных, для сопоставимости переоценим модель с помощью plm, уточнив тип модели: model = “pooling”.

Сформулируем нулевую гипотезу, которую мы проверяем: \[\gamma_{1} = \gamma_{2} \cdots = \gamma_{N-1} = 0\]

pooled <- plm(fh_polity~state_capacity, data=panel, model="pooling")
pFtest(fe, pooled)
## 
##  F test for individual effects
## 
## data:  fh_polity ~ state_capacity
## F = 14.702, df1 = 26, df2 = 107, p-value < 2.2e-16
## alternative hypothesis: significant effects

По результатам F-теста можно сделать вывод о том, что FE-модель более предпочтительна.