Links to sources:

https://cran.r-project.org/web/packages/pdynmc/vignettes/pdynmc-introLong.pdf

https://cran.r-project.org/web/packages/pdynmc/pdynmc.pdf

The following research is conducted by Maria R. Koldasheva and Nikolai A. Popov.

Libraries

library(dplyr) 
library(psych)
library(readxl)
library(writexl)
library(kableExtra)

# Econometrics
library(tidyverse)
library(plm) # panel models
library(sandwich) #covariance matrixes
library(lmtest) # tests
library(xtable)# latex tables
library(stargazer) # latex regression tables
library(ggpubr) # correlation test


# Dynamic models
library(pdynmc)

# To fix RPubs
# install.packages('rsconnect')
library(rsconnect)

Downloading the data

# Downloading
long_data <- read_csv('C:/Users/Kolya/Documents/New_studies/Diploma/Data/Data_use/Long_united_data.csv',
                      show_col_types = FALSE)

long_data

Dataset for 2015-2019

We exclude 2014 and 2020 as crisis years, to check whether our results are robust.

long_data <- subset(long_data, year %in% c(2015:2019))

Model construction

Specifications, Panel models

Theil index specification

Theil <- log(GRP) ~ lag(log(GRP)) + lag(Div) + Edu + Pop + lag(Inv) + lag(FDI) + log(Imp)

EM and IM specification

Margins <- log(GRP) ~ lag(log(GRP)) + lag(EM) + lag(IM) + Edu + Pop + lag(Inv) +lag(FDI)

Manually making logarithms for “pdynmc”

# new dataset
long_data_log = data.frame(long_data)

# Check if memory addresses are the same
tracemem(long_data_log) == tracemem(long_data)
# FALSE expected

# logarithm of appropriate columns
long_data_log$GRP=log(long_data_log$GRP)
long_data_log$Imp=log(long_data_log$Imp)
long_data_log$GRPL=log(long_data_log$GRPL)
colnames(long_data_log)

Model setting

Estimation: system GMM

Theil

Theil <- log(GRP) ~ lag(log(GRP)) + lag(Div) + Edu + Pop + lag(Inv) + lag(FDI) + log(Imp)

# Theil
Theil_d <- pdynmc(dat = long_data_log, varname.i = "region", varname.t = "year",
use.mc.diff = TRUE, use.mc.lev = FALSE, use.mc.nonlin = FALSE, inst.stata = TRUE,
include.y = TRUE, varname.y = "GRP", lagTerms.y = 1,
fur.con = TRUE, fur.con.diff = TRUE, fur.con.lev = FALSE,
varname.reg.fur = c("DivL", "Edu", "Pop", "InvL", "FDIL", "Imp"),
lagTerms.reg.fur = c(0,0,0,0,0,0),
include.dum = TRUE, dum.diff = TRUE, dum.lev = FALSE, varname.dum = "year",
w.mat.stata = TRUE, std.err = "dbl.corrected", # Hwang et al. (2021)
estimation = "onestep", opt.meth = "none")
## tracemem[0x00000262d854f510 -> 0x00000262d8554790]: pdynmc eval eval eval_with_user_handlers withVisible withCallingHandlers handle timing_fn evaluate_call <Anonymous> evaluate in_dir in_input_dir eng_r block_exec call_block process_group.block process_group withCallingHandlers process_file <Anonymous> <Anonymous> 
## tracemem[0x00000262d8554790 -> 0x00000262d85537d0]: $<-.data.frame $<- pdynmc eval eval eval_with_user_handlers withVisible withCallingHandlers handle timing_fn evaluate_call <Anonymous> evaluate in_dir in_input_dir eng_r block_exec call_block process_group.block process_group withCallingHandlers process_file <Anonymous> <Anonymous>
# w.mat = "zero.cov", 

Theil_d_sum <- summary(Theil_d)
Theil_d_sum
## 
## Dynamic linear panel estimation (onestep)
## Estimation steps: 1
## 
## Coefficients:
##           Estimate    Std.Err z-value Pr(>|z|)   
## L1.GRP   0.0459660  0.2908427   0.158  0.87446   
## L0.DivL  0.0095977  0.0131267   0.731  0.46478   
## L0.Edu  -0.0008814  0.0006573  -1.341  0.17992   
## L0.Pop  -0.0009117  0.0029735  -0.307  0.75884   
## L0.InvL  0.1776548  0.1858710   0.956  0.33907   
## L0.FDIL -0.1018030  0.0894338  -1.138  0.25512   
## L0.Imp  -0.0090019  0.0171567  -0.525  0.59958   
## 2017     0.0299152  0.0175772   1.702  0.08876 . 
## 2018     0.0849823  0.0317961   2.673  0.00752 **
## 2019     0.1054123  0.0434944   2.424  0.01535 * 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##  15 total instruments are employed to estimate 10 parameters
##  6 linear (DIF) 
##  6 further controls (DIF) 
##  3 time dummies (DIF) 
##  
## J-Test (overid restrictions):  7 with 5 DF, pvalue: 0.2209
## F-Statistic (slope coeff):  0.72 with 7 DF, pvalue: 0.9981
## F-Statistic (time dummies):  1.6 with 3 DF, pvalue: 0.6587

J-Test: null hypothesis is not rejected - > our set of IVs is jointly zero F-Statistic (slope coeff): null hypothesis is not rejected - > our betas are jointly zero (bad model) F-Statistic (time dummies): null hypothesis is rejected - > our betas are not jointly zero (good dummies model)

Margins

Margins <- log(GRP) ~ lag(log(GRP)) + lag(EM) + lag(IM) + Edu + Pop + lag(Inv) +lag(FDI)

# Margins
Margins_d <- pdynmc(dat = long_data_log, varname.i = "region", varname.t = "year",
use.mc.diff = TRUE, use.mc.lev = FALSE, use.mc.nonlin = FALSE, inst.stata = TRUE,
include.y = TRUE, varname.y = "GRP", lagTerms.y = 1,
fur.con = TRUE, fur.con.diff = TRUE, fur.con.lev = FALSE,
varname.reg.fur = c("EML", "IML", "Edu", "Pop", "InvL", "FDIL", "Imp"),
lagTerms.reg.fur = c(0,0,0,0,0,0,0),
include.dum = TRUE, dum.diff = TRUE, dum.lev = FALSE, varname.dum = "year",
w.mat.stata = TRUE, std.err = "dbl.corrected", # Hwang et al. (2021)
estimation = "onestep", opt.meth = "none")
## tracemem[0x00000262d854f510 -> 0x00000262d85480d0]: pdynmc eval eval eval_with_user_handlers withVisible withCallingHandlers handle timing_fn evaluate_call <Anonymous> evaluate in_dir in_input_dir eng_r block_exec call_block process_group.block process_group withCallingHandlers process_file <Anonymous> <Anonymous> 
## tracemem[0x00000262d85480d0 -> 0x00000262d8547d10]: $<-.data.frame $<- pdynmc eval eval eval_with_user_handlers withVisible withCallingHandlers handle timing_fn evaluate_call <Anonymous> evaluate in_dir in_input_dir eng_r block_exec call_block process_group.block process_group withCallingHandlers process_file <Anonymous> <Anonymous>
# w.mat = "zero.cov",

Margins_d_sum <- summary(Margins_d)
Margins_d_sum
## 
## Dynamic linear panel estimation (onestep)
## Estimation steps: 1
## 
## Coefficients:
##           Estimate    Std.Err z-value Pr(>|z|)   
## L1.GRP   0.0981078  0.2517465   0.390  0.69654   
## L0.EML   0.0086024  0.0231197   0.372  0.70989   
## L0.IML  -0.0005594  0.0007777  -0.719  0.47214   
## L0.Edu  -0.0008666  0.0006543  -1.324  0.18550   
## L0.Pop  -0.0025527  0.0043447  -0.588  0.55653   
## L0.InvL  0.2044539  0.1603735   1.275  0.20231   
## L0.FDIL -0.1063605  0.0926739  -1.148  0.25097   
## L0.Imp  -0.0126473  0.0181670  -0.696  0.48643   
## 2017     0.0282322  0.0160022   1.764  0.07773 . 
## 2018     0.0827926  0.0299156   2.768  0.00564 **
## 2019     0.1032482  0.0398312   2.592  0.00954 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##  16 total instruments are employed to estimate 11 parameters
##  6 linear (DIF) 
##  7 further controls (DIF) 
##  3 time dummies (DIF) 
##  
## J-Test (overid restrictions):  6.89 with 5 DF, pvalue: 0.2288
## F-Statistic (slope coeff):  0.91 with 8 DF, pvalue: 0.9987
## F-Statistic (time dummies):  1.65 with 3 DF, pvalue: 0.6485

J-Test: null hypothesis is not rejected - > our set of IVs is jointly zero F-Statistic (slope coeff): null hypothesis is not rejected - > our betas are jointly zero (bad model) F-Statistic (time dummies): null hypothesis is rejected - > our betas are not jointly zero (good dummies model) # United table creating

# Table for Theil
united_table_dynamic_panels_Theil <-
cbind(Theil_d_sum$coefficients[c(1:(nrow(Theil_d_sum$coefficients)-3)), c(1, 2, 4)])

# Table for Margins
united_table_dynamic_panels_Margins <- 
cbind(Margins_d_sum$coefficients[c(1:(nrow(Margins_d_sum$coefficients)-3)), c(1, 2, 4)])

# matrixes into dataframe
united_table_dynamic_panels_Theil <- as.data.frame(united_table_dynamic_panels_Theil)

united_table_dynamic_panels_Margins <- as.data.frame(united_table_dynamic_panels_Margins)

# renaming
names(united_table_dynamic_panels_Theil) <- c('Est_Theil', 'SE_Theil', 'Pvalue_Theil')

names(united_table_dynamic_panels_Margins) <- c('Est_Margins', 'SE_Margins', 'Pvalue_Margins')

united_table_dynamic_panels_Theil
united_table_dynamic_panels_Margins

Table

library(modelsummary)
## 
## Присоединяю пакет: 'modelsummary'
## Следующий объект скрыт от 'package:psych':
## 
##     SD
ti_Theil <- data.frame(
  term = c('L1.GRP', 'L1.Div', 'L0.Edu', 'L0.Pop', 'L1.Inv', 'L1.FDI', 'L0.Imp'),
  estimate = united_table_dynamic_panels_Theil$Est_Theil,
  std.error = united_table_dynamic_panels_Theil$SE_Theil,
  p.value = united_table_dynamic_panels_Theil$Pvalue_Theil)

ti_Margins <- data.frame(
  term = c('L1.GRP', 'L1.IM', 'L1.EM', 'L0.Edu', 'L0.Pop', 'L1.Inv', 'L1.FDI', 'L0.Imp'),
  estimate = united_table_dynamic_panels_Margins$Est_Margins,
  std.error = united_table_dynamic_panels_Margins$SE_Margins,
  p.value = united_table_dynamic_panels_Margins$Pvalue_Margins)

gl <- data.frame(
  Observations = "546",
  Model = "System GMM")

mod_Theil <- list(tidy = ti_Theil, glance = gl)
mod_Margins <- list(tidy = ti_Margins, glance = gl)

# Model class
class(mod_Theil) <- "modelsummary_list"
class(mod_Margins) <- "modelsummary_list"


# created named list
models <- list()
models[['log(GRP) (1)']] <- mod_Theil
models[['log(GRP) (2)']] <- mod_Margins

# Model
modelsummary(models, stars = T, title = 'Dynamic model', fmt = 4, coef_map = c('L1.GRP', 'L1.Div',
                                                                               'L1.EM', 'L1.IM',
                                                                               'L0.Edu', 'L0.Pop',
                                                                               'L1.Inv','L1.FDI',
                                                                               'L0.Imp'))
Dynamic model
log(GRP) (1)  log(GRP) (2)
L1.GRP 0.0460 0.0981
(0.2908) (0.2517)
L1.Div 0.0096
(0.0131)
L1.EM -0.0006
(0.0008)
L1.IM 0.0086
(0.0231)
L0.Edu -0.0009 -0.0009
(0.0007) (0.0007)
L0.Pop -0.0009 -0.0026
(0.0030) (0.0043)
L1.Inv 0.1777 0.2045
(0.1859) (0.1604)
L1.FDI -0.1018 -0.1064
(0.0894) (0.0927)
L0.Imp -0.0090 -0.0126
(0.0172) (0.0182)
Observations 546 546
Model System GMM System GMM
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001
#output = 'latex'
# modelsummary(mod)

SAME CODE, BUT FOR ROBUST DATASET (NO MSK, SPB)

Downloading the data, robust

# Downloading
long_data_robust <- read_csv(
  'C:/Users/Kolya/Documents/New_studies/Diploma/Data/Data_use/Long_united_data_robust.csv',
  show_col_types = FALSE)

long_data_robust

Dataset for 2015-2019, robust

We exclude 2014 and 2020 as crisis years, to check whether our results are robust.

long_data_robust <- subset(long_data_robust, year %in% c(2015:2019))

Model construction, robust

Manually making logarithms for “pdynmc”, robust

# new dataset
long_data_robust_log = data.frame(long_data_robust)

# Check if memory addresses are the same
tracemem(long_data_robust_log) == tracemem(long_data_robust)
# FALSE expected

# logarithm of appropriate columns
long_data_robust_log$GRP=log(long_data_robust_log$GRP)
long_data_robust_log$Imp=log(long_data_robust_log$Imp)
long_data_robust_log$GRPL=log(long_data_robust_log$GRPL)
colnames(long_data_robust_log)

Model setting, robust

Theil_robust, robust

Theil_robust <- log(GRP) ~ lag(log(GRP)) + lag(Div) + Edu + Pop + lag(Inv) + lag(FDI) + log(Imp)

# Theil_robust
Theil_robust_d <- pdynmc(dat = long_data_robust_log, varname.i = "region", varname.t = "year",
use.mc.diff = TRUE, use.mc.lev = FALSE, use.mc.nonlin = FALSE, inst.stata = TRUE,
include.y = TRUE, varname.y = "GRP", lagTerms.y = 1,
fur.con = TRUE, fur.con.diff = TRUE, fur.con.lev = FALSE,
varname.reg.fur = c("DivL", "Edu", "Pop", "InvL", "FDIL", "Imp"),
lagTerms.reg.fur = c(0,0,0,0,0,0),
include.dum = TRUE, dum.diff = TRUE, dum.lev = FALSE, varname.dum = "year",
w.mat.stata = TRUE, std.err = "dbl.corrected", # Hwang et al. (2021)
estimation = "onestep", opt.meth = "none")
## tracemem[0x00000262d8556350 -> 0x00000262d8556a10]: pdynmc eval eval eval_with_user_handlers withVisible withCallingHandlers handle timing_fn evaluate_call <Anonymous> evaluate in_dir in_input_dir eng_r block_exec call_block process_group.block process_group withCallingHandlers process_file <Anonymous> <Anonymous> 
## tracemem[0x00000262d8556a10 -> 0x00000262d85564d0]: $<-.data.frame $<- pdynmc eval eval eval_with_user_handlers withVisible withCallingHandlers handle timing_fn evaluate_call <Anonymous> evaluate in_dir in_input_dir eng_r block_exec call_block process_group.block process_group withCallingHandlers process_file <Anonymous> <Anonymous>
# w.mat = "zero.cov", 

Theil_robust_d_sum <- summary(Theil_robust_d)
Theil_robust_d_sum
## 
## Dynamic linear panel estimation (onestep)
## Estimation steps: 1
## 
## Coefficients:
##           Estimate    Std.Err z-value Pr(>|z|)  
## L1.GRP   0.0913607  0.3167996   0.288   0.7733  
## L0.DivL  0.0058302  0.0158164   0.369   0.7121  
## L0.Edu  -0.0009219  0.0006888  -1.338   0.1809  
## L0.Pop  -0.0006918  0.0014506  -0.477   0.6334  
## L0.InvL  0.1978262  0.2047430   0.966   0.3340  
## L0.FDIL -0.1033001  0.0943469  -1.095   0.2735  
## L0.Imp  -0.0099062  0.0178129  -0.556   0.5782  
## 2017     0.0276186  0.0182355   1.515   0.1298  
## 2018     0.0811092  0.0341460   2.375   0.0175 *
## 2019     0.0995659  0.0466591   2.134   0.0328 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##  15 total instruments are employed to estimate 10 parameters
##  6 linear (DIF) 
##  6 further controls (DIF) 
##  3 time dummies (DIF) 
##  
## J-Test (overid restrictions):  7.5 with 5 DF, pvalue: 0.1861
## F-Statistic (slope coeff):  0.68 with 7 DF, pvalue: 0.9985
## F-Statistic (time dummies):  1.33 with 3 DF, pvalue: 0.7212

J-Test: null hypothesis is not rejected - > our set of IVs is jointly zero F-Statistic (slope coeff): null hypothesis is not rejected - > our betas are jointly zero (bad model) F-Statistic (time dummies): null hypothesis is rejected - > our betas are not jointly zero (good dummies model)

Margins_robust

Margins_robust <- log(GRP) ~ lag(log(GRP)) + lag(EM) + lag(IM) + Edu + Pop + lag(Inv) +lag(FDI)

# Margins_robust
Margins_robust_d <- pdynmc(dat = long_data_robust_log, varname.i = "region", varname.t = "year",
use.mc.diff = TRUE, use.mc.lev = FALSE, use.mc.nonlin = FALSE, inst.stata = TRUE,
include.y = TRUE, varname.y = "GRP", lagTerms.y = 1,
fur.con = TRUE, fur.con.diff = TRUE, fur.con.lev = FALSE,
varname.reg.fur = c("EML", "IML", "Edu", "Pop", "InvL", "FDIL", "Imp"),
lagTerms.reg.fur = c(0,0,0,0,0,0,0),
include.dum = TRUE, dum.diff = TRUE, dum.lev = FALSE, varname.dum = "year",
w.mat.stata = TRUE, std.err = "dbl.corrected", # Hwang et al. (2021)
estimation = "onestep", opt.meth = "none")
## tracemem[0x00000262d8556350 -> 0x00000262d8548cd0]: pdynmc eval eval eval_with_user_handlers withVisible withCallingHandlers handle timing_fn evaluate_call <Anonymous> evaluate in_dir in_input_dir eng_r block_exec call_block process_group.block process_group withCallingHandlers process_file <Anonymous> <Anonymous> 
## tracemem[0x00000262d8548cd0 -> 0x00000262d8548e50]: $<-.data.frame $<- pdynmc eval eval eval_with_user_handlers withVisible withCallingHandlers handle timing_fn evaluate_call <Anonymous> evaluate in_dir in_input_dir eng_r block_exec call_block process_group.block process_group withCallingHandlers process_file <Anonymous> <Anonymous>
# w.mat = "zero.cov",

Margins_robust_d_sum <- summary(Margins_robust_d)
Margins_robust_d_sum
## 
## Dynamic linear panel estimation (onestep)
## Estimation steps: 1
## 
## Coefficients:
##           Estimate    Std.Err z-value Pr(>|z|)  
## L1.GRP   0.1251321  0.2666220   0.469   0.6391  
## L0.EML   0.0080904  0.0242749   0.333   0.7391  
## L0.IML  -0.0006511  0.0007936  -0.820   0.4122  
## L0.Edu  -0.0009081  0.0006827  -1.330   0.1835  
## L0.Pop  -0.0012888  0.0017751  -0.726   0.4678  
## L0.InvL  0.2141202  0.1747507   1.225   0.2206  
## L0.FDIL -0.1082829  0.0967774  -1.119   0.2631  
## L0.Imp  -0.0134586  0.0189388  -0.711   0.4771  
## 2017     0.0268686  0.0164012   1.638   0.1014  
## 2018     0.0800734  0.0313938   2.551   0.0107 *
## 2019     0.0986844  0.0416427   2.370   0.0178 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##  16 total instruments are employed to estimate 11 parameters
##  6 linear (DIF) 
##  7 further controls (DIF) 
##  3 time dummies (DIF) 
##  
## J-Test (overid restrictions):  7.29 with 5 DF, pvalue: 0.1998
## F-Statistic (slope coeff):  0.94 with 8 DF, pvalue: 0.9986
## F-Statistic (time dummies):  1.4 with 3 DF, pvalue: 0.7058

J-Test: null hypothesis is not rejected - > our set of IVs is jointly zero F-Statistic (slope coeff): null hypothesis is not rejected - > our betas are jointly zero (bad model) F-Statistic (time dummies): null hypothesis is rejected - > our betas are not jointly zero (good dummies model) # United table creating, robust

# Table for Theil_robust
united_table_dynamic_panels_Theil_robust <-
cbind(Theil_robust_d_sum$coefficients[c(1:(nrow(Theil_robust_d_sum$coefficients)-3)), c(1, 2, 4)])

# Table for Margins_robust
united_table_dynamic_panels_Margins_robust <- 
cbind(Margins_robust_d_sum$coefficients[c(1:(nrow(Margins_robust_d_sum$coefficients)-3)), c(1, 2, 4)])

# matrixes into dataframe
united_table_dynamic_panels_Theil_robust <- as.data.frame(united_table_dynamic_panels_Theil_robust)

united_table_dynamic_panels_Margins_robust <- as.data.frame(united_table_dynamic_panels_Margins_robust)

# renaming
names(united_table_dynamic_panels_Theil_robust) <- c('Est_Theil_robust', 'SE_Theil_robust', 'Pvalue_Theil_robust')

names(united_table_dynamic_panels_Margins_robust) <- c('Est_Margins_robust', 'SE_Margins_robust', 'Pvalue_Margins_robust')

united_table_dynamic_panels_Theil_robust
united_table_dynamic_panels_Margins_robust

Table, robust

library(modelsummary)
ti_Theil_robust <- data.frame(
  term = c('L1.GRP', 'L1.Div', 'L0.Edu', 'L0.Pop', 'L1.Inv', 'L1.FDI', 'L0.Imp'),
  estimate = united_table_dynamic_panels_Theil_robust$Est_Theil_robust,
  std.error = united_table_dynamic_panels_Theil_robust$SE_Theil_robust,
  p.value = united_table_dynamic_panels_Theil_robust$Pvalue_Theil_robust)

ti_Margins_robust <- data.frame(
  term = c('L1.GRP', 'L1.IM', 'L1.EM', 'L0.Edu', 'L0.Pop', 'L1.Inv', 'L1.FDI', 'L0.Imp'),
  estimate = united_table_dynamic_panels_Margins_robust$Est_Margins_robust,
  std.error = united_table_dynamic_panels_Margins_robust$SE_Margins_robust,
  p.value = united_table_dynamic_panels_Margins_robust$Pvalue_Margins_robust)

gl <- data.frame(
  Observations = "532",
  Model = "System GMM")

mod_Theil_robust <- list(tidy = ti_Theil_robust, glance = gl)
mod_Margins_robust <- list(tidy = ti_Margins_robust, glance = gl)

# Model class
class(mod_Theil_robust) <- "modelsummary_list"
class(mod_Margins_robust) <- "modelsummary_list"


# created named list
models <- list()
models[['log(GRP) (1)']] <- mod_Theil_robust
models[['log(GRP) (2)']] <- mod_Margins_robust

# Model
modelsummary(models, stars = T, title = 'Dynamic model', fmt = 4, coef_map = c('L1.GRP', 'L1.Div',
                                                                               'L1.EM', 'L1.IM',
                                                                               'L0.Edu', 'L0.Pop',
                                                                               'L1.Inv','L1.FDI',
                                                                               'L0.Imp'))
Dynamic model
log(GRP) (1)  log(GRP) (2)
L1.GRP 0.0914 0.1251
(0.3168) (0.2666)
L1.Div 0.0058
(0.0158)
L1.EM -0.0007
(0.0008)
L1.IM 0.0081
(0.0243)
L0.Edu -0.0009 -0.0009
(0.0007) (0.0007)
L0.Pop -0.0007 -0.0013
(0.0015) (0.0018)
L1.Inv 0.1978 0.2141
(0.2047) (0.1748)
L1.FDI -0.1033 -0.1083
(0.0943) (0.0968)
L0.Imp -0.0099 -0.0135
(0.0178) (0.0189)
Observations 532 532
Model System GMM System GMM
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001
#output = 'latex'
# modelsummary(mod)

Final table

# created named list
models <- list()
models[['log(GRP) (1)']] <- mod_Theil
models[['log(GRP) (2)']] <- mod_Theil_robust
models[['log(GRP) (3)']] <- mod_Margins
models[['log(GRP) (4)']] <- mod_Margins_robust
# Model
modelsummary(models, stars = T, title = 'Dynamic models, robust (no crisis years)',
             fmt = 4, coef_map = c('L1.GRP', 'L1.Div','L1.EM', 'L1.IM','L0.Edu', 
                                   'L0.Pop', 'L1.Inv','L1.FDI','L0.Imp'))
Dynamic models, robust (no crisis years)
log(GRP) (1)  log(GRP) (2)  log(GRP) (3)  log(GRP) (4)
L1.GRP 0.0460 0.0914 0.0981 0.1251
(0.2908) (0.3168) (0.2517) (0.2666)
L1.Div 0.0096 0.0058
(0.0131) (0.0158)
L1.EM -0.0006 -0.0007
(0.0008) (0.0008)
L1.IM 0.0086 0.0081
(0.0231) (0.0243)
L0.Edu -0.0009 -0.0009 -0.0009 -0.0009
(0.0007) (0.0007) (0.0007) (0.0007)
L0.Pop -0.0009 -0.0007 -0.0026 -0.0013
(0.0030) (0.0015) (0.0043) (0.0018)
L1.Inv 0.1777 0.1978 0.2045 0.2141
(0.1859) (0.2047) (0.1604) (0.1748)
L1.FDI -0.1018 -0.1033 -0.1064 -0.1083
(0.0894) (0.0943) (0.0927) (0.0968)
L0.Imp -0.0090 -0.0099 -0.0126 -0.0135
(0.0172) (0.0178) (0.0182) (0.0189)
Observations 546 532 546 532
Model System GMM System GMM System GMM System GMM
+ p < 0.1, * p < 0.05, ** p < 0.01, *** p < 0.001

Coefficients are different, significance was only lost.