Introduction

This project focuses on observation of economic convergence, explanation of its determinants, and assessment of deductive mathematical models of Solow (1956) and Romer (1986). It combines multiple methods and brings focus onto merging variable both occurent in existing research and those accessible through IMF, the World Bank, and the UNDP in order to confirm or reject their relevance.

Disclaimer

This document accompanies a final year Economics Project paper, which is to be submitted by the author. Any wording match between this page and the final submission, which is not found elsewhere is hereby asked to be ignored. Everything contained herein is solely product of my work, unless cited.

Workspace preparation

Inclusion of packages to be used and definition of functions

setwd("C:/Users/matt/MEGA/School/ECN3020/CW02")
#Data handling packages
require(readr)
require(reshape2)
#Mathematical packages
require(plm)
#Result presentation packages
require(stargazer)
require(knitr)
#DATA API packages
require(WDI)
require(maddison)

#import of variables crucial for most operations
ccdat <- function(all){
    #Specify countries of interest
    countries <<- c("GBR", "AUT", "BEL", "DNK", "FRA", "DEU", "ITA", "LUX",
                   "NLD", "NOR", "SWE", "CHE", "FIN", "GRC", "ISL", "IRL",
                   "PRT", "ESP", "BLR", "ALB", "BGR", "MDA", "UKR", "CZE",
                   "SVK", "EST", "LVA", "SRB", "MNE", "HUN", "LTU", "HRV",
                   "SVN", "MKD", "BIH", "POL", "ROU")
    #Load dataset of unique country identifiers 
    #created for purpose of this exercise
    CCX <<- read_csv("./DAT/CCX.csv")
    if(!all){CCX <<- CCX[CCX$ccx %in% countries,]}
                  
}
#World Bank data cleaning function to make the format match
wbimport <- function(icode, nic=icode, cntrs=CCX$ci){
    WD0 <- WDI(cntrs, icode, 1985, 2015)
    iname0 <- WDIsearch(icode, "indicator")
    iname <- ifelse(is.vector(iname0), iname0[2], ifelse(length(iname0)==0, nic, iname0[1,2]))
    WD1 <- data.frame(ci=WD0$iso2c, indicator=iname, 
                       ic=nic, year=WD0$year, value=as.numeric(WD0[,3]))
    WD2 <- merge(CCX[,c(1, 3, 5)], WD1, by="ci")[,-1]
    return(WD2)
}
#Might be useful, eventually
do.matrix <- function(x){
    reshape(data=pPAD[pPAD$year %in% 2001:2010 ,c("year", "ccx", x)],
        idvar="ccx",
        v.names = x,
        timevar = "year",
        direction="wide",
        new.row.names = unique(pPAD$ccx)
        )[,-1]
}
#Get coefficient column from a model, add significance stars
getcof <- function(model){
    paste0(
        ifelse(model$coefficients<0, "", "+"),
        round(model$coefficients, 2),
        ifelse(summary(model)$coefficients[4]>0.05, ".", "*")
        )
}

Variable preparation

Selection of appropriate variables from appropriate sources has proven to be rather complex task, for which reason the chunks of code shown within this section had been repeated multiple times before an appropriate dataset in an appropriate form was created. While research suggested what the appropriate variables should be, choice of their measures was especially hard, for they seem to be hardly ever collected.

Macroeconomic Indicators

Because the data to be used was taken from multiple sources in various formats, conversion into a single dataset featuring countries of concern only shall be extracted. The following process runs through all the datasets, applying appropriate structural conversion, finally composing one single long-format dataset.

IMF

The database access API provided by IMF is fairly unreliable and not the easiest to deal with, which is why CSV import is used instead.

#Load countries of concern and a list of countries by country code
ccdat(TRUE)

#CDIS <- read_csv("./DAT/CDIS.csv")
#    colnames(CDIS) <- c("cntry", "cc", "indicator", "ic", "counterpart", "ccc", "year", "value", "status", "blank")
#    CDIS <- CDIS[CDIS$ccc==1,c(1:4, 7:8)] #Only select counterpart World
DOT <- read_csv("./DAT/DOT.csv")
    colnames(DOT) <- c("cntry", "cc", "indicator", "ic", "counterpart", "ccc", "year", "value", "status", "blank")
    DOT <- DOT[DOT$ccc==1,c(1:4, 7:8)] #Only select counterpart World
ED <- read_csv("./DAT/ED.csv")
    colnames(ED) <- c("cntry", "cc", "indicator", "ic", "year", "value", "blank")
    ED <- ED[,1:6]
FM <- read_csv("./DAT/FM.csv")
    colnames(FM) <- c("cntry", "cc", "indicator", "ic", "year", "value", "status", "blank")
    FM <- FM[,1:6] 
HPDD <- read_csv("./DAT/HPDD.csv")
    colnames(HPDD) <- c("cntry", "cc", "indicator", "ic", "year", "value", "blank")
    HPDD <- HPDD[,1:6]
IRFCL <- read_csv("./DAT/IRFCL.csv")
    colnames(IRFCL) <- c("cntry", "cc", "indicator", "ic", "sec_name", "sec_code", "year", "value", "status", "blank")
    IRFCL<-IRFCL[,c(1:4, 7:8)]
    IRFCL$indicator <- "Nominal Value of Foreign Reserves [USD]"

#Slightly more complex
EQ <- read_csv("./DAT/EQ.csv")
    colnames(EQ) <- c("cntry", "cc", "indicator", "ic", "prod_name", "prod_code", "year", "value", "blank")
    EQ$prod_name <- ifelse(EQ$prod_name!="All products", substring(EQ$prod_name, 7), EQ$prod_name)
    EQ$prod_code <- ifelse(EQ$prod_code!="All products", EQ$prod_code, "ALL")
    EQ <- data.frame(cntry=EQ$cntry, cc=EQ$cc,
                     indicator=paste0(EQ$indicator, ": ", EQ$prod_name), #
                     ic=paste0(EQ$ic, ":", EQ$prod_code), #8*4=32 values; BEC1-7 and all
                     year=EQ$year, value=EQ$value)
    EQ <- EQ[!(substring(as.character(EQ$ic),1,1) %in% c("X", "q", "u")) ,]

#Consolidating IMF data
IMFBULK <- rbind(DOT, ED, FM, HPDD, IRFCL, EQ)
    #Replace IMF cc with more comprehensible ISO3166 ccx
    IMFBULK <- merge(CCX[,1:3], IMFBULK[,-1], by="cc")[,-1]
    #Pick the rows
    IMFBULK <- IMFBULK[IMFBULK$ccx %in% countries & IMFBULK$year %in% 1980:2015,]
#Save CSV
write_csv(IMFBULK, "./DAT/IMFBULK.csv")
rm(list = setdiff(ls(), lsf.str()))

WorldBank

This is an excellent API. Easy to use.

#Load countries of concern and a list of countries by country code
ccdat(FALSE)

WBBULK <- rbind(
    wbimport("NY.GDP.PCAP.KD", "WBGDP"),
    wbimport("NY.GDP.PCAP.KD.ZG", "WBGDPYOY"),
    wbimport("SE.TER.ENRR", "EDU3er"),
    wbimport("SE.SEC.ENRR", "EDU2er"),
    wbimport("IP.PAT.RESD", "PAT"),
    wbimport("SP.POP.DPND", "DEPRATIO"),  
    wbimport("NY.GNP.PCAP.PP.KD", "GNIpcPPP"),
    wbimport("SP.POP.TOTL", "POPULATION"),
    wbimport("SL.UEM.TOTL.ZS", "UNEMPLOYMENT"),
    wbimport("NE.GDI.TOTL.ZS", "GKF"),
    wbimport("BX.KLT.DINV.CD.WD", "FDI"),
    wbimport("FI.RES.TOTL.CD", "WBFRUSD"),
    wbimport("GC.XPN.TOTL.GD.ZS", "GOVT"),
    wbimport("NE.IMP.GNFS.ZS", "IMPORTS"),
    wbimport("NE.EXP.GNFS.ZS", "EXPORTS"),
    wbimport("NV.SRV.TETC.ZS", "SERVICEsva"),
    wbimport("FP.CPI.TOTL.ZG", "INFLATION")
)

#Save CSV
write_csv(WBBULK, "./DAT/WBBULK.csv")
#Clear environment
rm(list = setdiff(ls(), lsf.str()))

Miscellaneous

#Load countries of concern and a list of countries by country code
ccdat(TRUE)
#Human Development Index - UNDP
#Single indicator, file import
HDI <- read_csv("./DAT/HDI_W.csv")
    HDI$Year <- factor(HDI$Year)
    HDI <- gather(HDI, "cntry", "value" , 2:38, factor_key=TRUE)
    HDI <- data.frame(cntry=HDI$cntry, indicator="HDI", ic="HDI", year=as.integer(as.character(HDI$Year)), value=HDI$value)
    HDI <- merge(CCX[,c(1,3)], HDI, by="cntry")
    
#GDP Per Capita; Historical - Maddison project
MDS <- maddison[(maddison$iso3c %in% countries) & (as.numeric(as.character.Date(maddison$year, "%Y")) %in% 1985:2010),c(4,6,1,3)]
MDS <- data.frame(cntry=MDS$country, ccx=MDS$iso3c, indicator="Maddison GDP per Capita (2010$)", ic="MDS", 
                  year=as.integer(as.character.Date(MDS$year, "%Y")), value=MDS$gdp_pc)
    
#OECD data on Currency Values
FX<-read_csv("./DAT/FX.csv")
    colnames(FX) <- c("ccx", "indicator", "subj", "ic", "freq", "year", "value", "fc")
    FX <- merge(CCX[,c(1,3)], FX, by="ccx")
    FX$indicator <- "Domestic currency [DOM]"
    FX <- FX[,c(2, 1, 3, 5, 7:8)]

#Stick it together
MSCBULK <- rbind(HDI, MDS, FX)
#Pick the rows
MSCBULK <- MSCBULK[MSCBULK$ccx %in% countries & MSCBULK$year %in% 1980:2015,]
#Save CSV
write_csv(MSCBULK, "./DAT/MSCBULK.csv")
#Clear the environment
rm(list = setdiff(ls(), lsf.str()))

Dummy variables

Now that all scale variables have been collected, dummy variables can be created to indicate whether given country was a member of given international organization in given year. In addition, dummies for period of great recession is to be added.

This piece of code is awful and I am not proud of it.

#Load countries and CCX
ccdat(FALSE)
#Import gathered data for dummies 
DUMW <- read_csv("./DAT/DUM.csv")
#Prep a data frame
DUMBULK <- expand.grid(ccx=CCX$ccx, ic=colnames(DUMW)[5:14], year=1985:2015)
DUMBULK$indicator <- paste(DUMBULK$ic, "membership [1/0]")
#Run a for loop to create long data for the orgs
for(j in 1:length(DUMW$ccx)){#Row -> Country
    for(i in 5:length(DUMW)){#Columns -> Organization
        for(yr in 1985:2015){#return membership binary for each year, country, and organization
            DUMBULK$value[DUMBULK$ic==colnames(DUMW)[i] & DUMBULK$year==yr & DUMBULK$ccx==DUMW$ccx[j]] <-
                as.numeric(ifelse(is.na(DUMW[j,i]), 0, ifelse(DUMW[j,i]<=yr,1,0)))
        }}}
#Create a dummy for The Great Recession
RNG <- 1985:2015
VAL <- ifelse(RNG %in% 2008:2009, 1,0)
for(k in 1:length(countries)){#Row -> Country
    DUMADD2 <- data.frame(ccx=countries[k], indicator=paste("Great Recession [1/0]"), ic="GR08", year=RNG, value=VAL)
    DUMBULK<- rbind(DUMBULK, DUMADD2)
}
DUMBULK <- merge(CCX[,c(1,3)], DUMBULK, by="ccx")[,c(2,1,5,3,4,6)]
DUMBULK$ic <- paste0("x_",DUMBULK$ic)
#Write to CSV
write_csv(DUMBULK, "./DAT/DUMBULK.csv")
#Clear environment
rm(list = setdiff(ls(), lsf.str()))

Data consolidation

With pieces of data from various sources successfuly saved in separate uniformely structured data files, we can put them together into a complete source data file.

#Stick the bulks together and save them
CAD <- rbind(read_csv("./DAT/IMFBULK.csv"),
             read_csv("./DAT/WBBULK.csv"),
             read_csv("./DAT/MSCBULK.csv"),
             read_csv("./DAT/DUMBULK.csv"))
write_csv(CAD, "./DAT/CAD.CSV")

#Clear environment
rm(list = setdiff(ls(), lsf.str()))

Transformation

In order to make the data work with plm package, it needs to be reshaped into wide panel data format.

#Specify countries of interest and CCX
ccdat(FALSE)
#import complete dataset
CAD <- read_csv("./DAT/CAD.csv")
#extract list of indicators and save to file
    ICX <- data.frame(ic=unique(CAD$ic), indicator=unique(CAD$indicator))
    write_csv(ICX, "./DAT/ICX.csv")
    
#Pick relevant columns
CAD <- CAD[,c("ccx", "year", "ic", "value")]
    CAD$ccx <- factor(CAD$ccx)
#Turn into a panel
PAD <- dcast(CAD, ccx + year ~ ic, value.var="value")
#Add country group identifier
PAD <- merge(CCX[,c("ccx", "cg")], PAD, by="ccx")
PAD <- PAD[PAD$year %in% 1990:2015,]

#Save to file
write_csv(PAD, "./DAT/PAD.csv")
#Environmentalists must love me. 
rm(list = setdiff(ls(), lsf.str()))

Calculation of secondary variables

#Specify countries of interest
ccdat(FALSE)
#Take the panel, interpret it as a panel
PAD <- pdata.frame(read_csv("./DAT/PAD.csv"), 
                   index=c("ccx", "year"))
#Take list of indicators
ICX <- read_csv("./DAT/ICX.csv")

#CALCULATE whatever additional variables you need here
PAD$POPg <- 100*diff(PAD$POPULATION)/lag(PAD$POPULATION)
ICX <- rbind(ICX, c("POPg", "Annual population growth rate (%)"))

PAD$GNIpcPPPg <- 100*diff(PAD$GNIpcPPP)/lag(PAD$GNIpcPPP)
ICX <- rbind(ICX, c("GNIpcPPPg", "GNI PPP 2000USD growth rate (%)"))

PAD$PATpcc <- 10000*PAD$PAT/PAD$POPULATION
ICX <- rbind(ICX, c("PATpcc", "Patent applications per 1000 citizens"))

PAD$FDIpcGDP <- 100*PAD$FDI/PAD$WBGDP
ICX <- rbind(ICX, c("FDIpcGDP", "FDI net inflows, %GDP"))

PAD$lWBGDP <- log(PAD$WBGDP)
ICX <- rbind(ICX, c("lWBGDP", "Log of WBGDP"))

#Save
write_csv(PAD, "./DAT/PAD.csv")
write_csv(ICX, "./DAT/ICX.csv")
#Clear environment
rm(list = setdiff(ls(), lsf.str()))

Integrity check

#Specify countries of interest
ccdat(FALSE)
#Take the panel, interpret it as panel
PAD <- pdata.frame(read_csv("./DAT/PAD.csv"), index=c("ccx", "year"))
#Create variable list
ICX <- read_csv("./DAT/ICX.csv")
#Determine number of missing observations per indicator
ICY <- ICX; ICY$missing <- NA
for(a in 1:length(ICY$ic)){
    ICY$missing[a] <- round(100*length(which(is.na(PAD[,ICY$ic[a]])))/length(PAD$ccx))
}

#Express proportion of missing data on the whole dataset
cat(paste0(round(mean(ICY$missing), 2),"% NA initially"))
#Remove variable whenever more than 25% of observations lack its value
PAD <- PAD[,!colnames(PAD) %in% ICY$ic[ICY$missing > 25]]
ICY <- ICY[ICY$missing < 25,]
cat(paste0(round(mean(ICY$missing), 2),"% NA after automatic removal"))

#Remove ugly variables manually
varsToRemove <- c("MDS", "PAT", "FDI", "value:5BEC", "value:1BEC", 
                  "value:4BEC", "value:2BEC", "value:6BEC", "value:3BEC",
                  "value:7BEC", "G_X_G01_GDP_PT")
PAD <- PAD[,!colnames(PAD) %in% varsToRemove]
ICY <- ICY[!ICY$ic %in% varsToRemove,]
cat(paste0(round(mean(ICY$missing), 2),"% NA after manual removal"))

#display a table
kable(ICY[order(ICY$missing, decreasing = TRUE),])

#Write new indicator list that to a file
write_csv(ICY, "./DAT/ICY.csv")
write_csv(PAD, "./DAT/PAD.csv")

#Clear the environment
rm(list = setdiff(ls(), lsf.str()))

Data inspection

Here we go, finally! First, we load the file we are going to work with

# Load country metadata
ccdat(FALSE)
# Load panel data
PAD <- pdata.frame(read_csv("./DAT/PAD.csv"), index=c("ccx", "year"))
# Load indicator metadata
ICY <- read_csv("./DAT/ICY.csv")
kable(ICY[order(ICY$missing, decreasing = TRUE),])

ic	indicator	missing
GGR_G01_GDP_PT	Revenue (% of GDP)	23
NATUSD	Domestic currency [DOM]	23
GOVT	Expense (% of GDP)	20
GNIpcPPPg	GNI PPP 2000USD growth rate (%)	17
value:ALL	Trade value of exports: All products	15
GNIpcPPP	GNIpcPPP	13
PATpcc	Patent applications per 1000 citizens	12
EDU3er	School enrollment, tertiary (% gross)	11
SERVICEsva	Services, etc., value added (% of GDP)	11
EDU2er	School enrollment, secondary (% gross)	10
FDIpcGDP	FDI net inflows, %GDP	10
GGXWDG_GDP	Debt to GDP Ratio	9
between_theil	Extensive Margin	9
within_theil	Intensive Margin	9
total_theil	Export Diversification Index	9
TXG_FOB_USD	Goods, Value of Exports, Free on board (FOB), US Dollars	8
TBG_USD	Goods, Value of Trade Balance, US Dollars	8
INFLATION	Inflation, consumer prices (annual %)	8
WBGDPYOY	GDP per capita growth (annual %)	7
WBFRUSD	Total reserves (includes gold, current US$)	7
WBGDP	GDP per capita (constant 2000 US$)	5
GKF	Gross capital formation (% of GDP)	5
IMPORTS	Imports of goods and services (% of GDP)	5
EXPORTS	Exports of goods and services (% of GDP)	5
HDI	HDI	5
lWBGDP	Log of WBGDP	5
UNEMPLOYMENT	Unemployment, total (% of total labor force)	4
POPg	Annual population growth rate (%)	4
DEPRATIO	Age dependency ratio (% of working-age population)	0
POPULATION	Population, total	0
x_GR08	Great Recession [1/0]	0
x_EU	EU membership [1/0]	0
x_CISFTA	CISFTA membership [1/0]	0
x_CEFTA	CEFTA membership [1/0]	0
x_WB	WB membership [1/0]	0
x_IMF	IMF membership [1/0]	0
x_VG	VG membership [1/0]	0
x_UN	UN membership [1/0]	0
x_NATO	NATO membership [1/0]	0
x_EURO	EURO membership [1/0]	0
x_OECD	OECD membership [1/0]	0

# Generate region metadata
REG <- list(codes=sort(unique(PAD$cg)),
            names=c("Balkan", "Eastern Bloc", "Post-Soviet", "W Europe"),
            pairs=list("EB", "BK", "PS", "WE",
                       c("EB", "WE"), c("BK", "WE"), c("PS", "WE"),
                       c("EB", "WE", "BK", "PS"))
            )

Plots

Dependent variables

GDP Growth and HDI

The two variables of interest, as mentioned in the Method section are annual GDP per capita growth, as a rather crude proxy for the rate of improvement of well-being in given country, and the Human Development Index UNDP(2017), which includes level of GDP per capita in addition to further development-related components.

Firstly, visualizing the cross-region differences in distribution of these in all concerned economies between years 1985 and 2015 is the figure below.

In case of GDP Growth, some suggestion of long-run cross-regional convergence exists, as the median growth rate in each region is higher than that of the Westewrn Europe. It is however also important to note that the variable seems far more stable in the west, while there are numerous, mainly negative outliers in all the other regions. It might be either due to the initial period of shock caused by political revolutions, or due to especially underperforming economy in the region.

For this reason, both options are considered further

par(mar=c(0,0,1.1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(2,0,0,0),
    cex.axis=0.7, #reduce size of axis numbers 
    cex.main=0.7, #Increase title size
    mfrow=c(1,2))

boxplot(WBGDPYOY ~ cg, data=PAD,
        horizontal=TRUE,  outline=TRUE,  axes=FALSE,
        notch=TRUE,
        whisklty=1, 
        staplelty=0,
        col="grey70",
        main="Real GDP per capita Growth [%]",
        boxwex=0.7,
        ylim=c(-15, 25),
        names=REG$names,
        pch=20,
        outwex=0.1
        )
abline(v=0, col="grey70")
axis(1, las=0, at=10*-2:3)
axis(4, lty=0, pos=1, las=2, at=1.5:4.5, labels=REG$names)
box()
boxplot(GNIpcPPPg ~ cg, data=PAD,
        horizontal=TRUE,  outline=TRUE,  axes=FALSE,
        notch=TRUE,
        whisklty=1, 
        staplelty=0,
        col="grey70",
        main="PPP GNI per capita growth [%]",
        boxwex=0.7,
        #ylim=c(0.5, 1),
        names=REG$names,
        pch=20,
        outwex=0.1
        )
abline(v=0, col="grey70")
axis(1, las=0, at=10*-3:2)
axis(2, lty=0, pos=-1, las=2, at=1.5:4.5, labels=REG$names)
box()

As it can be seen, the distributions are fairly similar within regions. Although Ukraine drags the PostSoviet Region median down, the high volatility of growth is common among these countries. Eastern Bloc experiences growth of lower volatility, but Romania and Bulgaria both seem to fit the Post-Soviet Country pattern, which makes sense due to the geographic proximity to this region. Otherwise, growth values in countries of EB, especially Poland, Slovakia, Czechia, and Hungary were closest to that of WE over this time. In case of the Balkan region, the differences in GDP growth distribution is larger.

As for HDI, Portugal may be seen as unfitting to the WE region, which has significantly larger median over the time. As is generally true, the maximum values may also be values from the last year, as the inspected variable here is a level as opposed to the percentage growth in case of the GDP. One would then assume that the larger the range, the more progress the country has experienced. The Eastern Bloc and the Post-Soviet region both have ranges generally greater than Western Europe, while some countries of the Balkan seem to stagnate.

It took me an hour and half to order and label the boxes correctly. Such stupid solution is one I haven’t seen before :D

par(mar=c(0,0,1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(1.5,0,0,0),
    cex.axis=0.6, #reduce size of axis numbers 
    cex.main=0.7, #Increase title size
    mfrow=c(1,2))

boxhue <- c(rep("grey70",length(unique(PAD$ccx[PAD$cg=="BK"]))),
            rep("grey100",length(unique(PAD$ccx[PAD$cg=="EB"]))),
            rep("grey70",length(unique(PAD$ccx[PAD$cg=="PS"]))),
            rep("grey100",length(unique(PAD$ccx[PAD$cg=="WE"]))))

boxplot(WBGDPYOY ~ paste0(cg,ccx), data=PAD, horizontal=TRUE, 
        axes=FALSE, ylim=c(-21, 30), boxwex=0, outline=FALSE)
        axis(1, las=0, pos=0, at=10*-2.5:2.5)
        axis(4, lty=3, pos=25, las=2, at=1:37, labels=unique(PAD$ccx[order(PAD$cg)]), tck=1, col="grey20")
        axis(4, lty=1, pos=25, las=2, at=1:37, labels=NA, tck=-0.05)

boxplot(WBGDPYOY ~ paste0(cg,ccx), data=PAD,
        horizontal=TRUE,  outline=TRUE,  axes=FALSE,
        whisklty=1, staplewex=0.9,
        main="Real GDP Growth [%]",
        boxwex=0.7, pch=20,
        outwex=0.1, col=boxhue,
        add=TRUE)
abline(v=0, col="grey70")

boxplot(GNIpcPPPg ~ paste0(cg,ccx), data=PAD, horizontal=TRUE, 
        axes=FALSE, ylim=c(-30, 21), boxwex=0, outline=FALSE)
    axis(1, las=0, pos=0, at=10*-2.5:2.5)
    axis(2, lty=3, pos=-25, las=2, at=1:37, labels=NA, tck=1, col="grey20")
    axis(2, lty=1, pos=-25, las=2, at=1:37, labels=NA, tck=-0.05)

boxplot(GNIpcPPPg ~ paste0(cg,ccx), data=PAD,
        horizontal=TRUE,  outline=TRUE,  axes=FALSE,
        whisklty=1, staplewex=0.9,
        main="PPP GNI per capita growth [%]",
        boxwex=0.7, pch=20,
        outwex=0.1, col=boxhue,
        add=TRUE)
abline(v=0 , col="grey70")

#xyplot(WBGDPYOY+GKF~year,groups=ccx,data=PAD[PAD$year>1990,], scales = list(y = list(relation = "free")), type='l')
rm(list=c("boxhue"))

Now, to inspect the over-time changes in distribution of values for both indicators, candlestick charts are shown below.

par(mar=c(0,0,0,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(3,2.5,1.1,0),
    cex.axis=0.9, #reduce size of axis numbers 
    cex.main=1, #Increase title size
    mfrow=c(4,2))
#Draw a boxplot series for each region
for(g in 1:length(REG$codes)){
    boxplot(WBGDPYOY ~ year, data=PAD[PAD$cg == REG$codes[g],],
        horizontal=FALSE,  outline=TRUE,  axes=FALSE,
        whisklty=1, 
        staplelty=0, 
        boxwex=0.75,
        ylim=c(-20, 25),
        pch=16)
    axis(2, las=2)
    text(0,22, REG$names[g], adj=0)
    abline(0,0, col="orangered")
    box()
    if(g==length(REG$codes)){axis(1, las=2, at=(2*1:13)-1, labels=seq(1990,2015,by=2))}

    
    boxplot(GNIpcPPPg ~ year, data=PAD[PAD$cg == REG$codes[g],],
        horizontal=FALSE,  outline=TRUE,  axes=FALSE,
        whisklty=1, 
        staplelty=0, 
        boxwex=0.75,
        ylim=c(-20, 25),
        pch=16)
    text(0,22, REG$names[g], adj=0)
    abline(0,0, col="orangered")
    box()
    if(g==length(REG$codes)){axis(1, las=2, at=(2*1:13)-1, labels=seq(1990,2015,by=2))}
    
    if(g==1){title(main="PPP GNI per capita growth [%]", outer=TRUE, adj=0.9)
            title(main="Real GDP Growth [%]", outer=TRUE, adj=0.15)}
}

rm(list=c("g"))

Testing for stationarity of the annual GDP growth is then the table below, showing p-value for null hypothesis of non stationarity. This implies that in most cases, the variable presents with unit root, for which reason the general linear models must utilize its second difference.

Latvia, Greece, Iceland, and Portugal… Screwit! :D

Independent Variables

ICY$unitRootp <- NA
ICY$unitRootdp <- NA

dPAD <- PAD[,1:3]
dPAD[,4:length(PAD)] <-NA
colnames(dPAD) <- colnames(PAD)

par(mar=c(0,0,1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(0,0,0,0),
    cex.axis=0.5, #reduce size of axis numbers 
    cex.main=0.8, #Increase title size
    mfrow=c(4,8))

plot(c(1,2,3), c(0,0,0),
         type="l", 
         ylim=c(0,1),
         main="Variable levels plot",
         axes=FALSE,
         frame.plot=TRUE)
text(1, 0.95, "x: time (1985 - 2015)", 0, cex=0.8)
text(1, 0.80, "y: variable (relative)", 0, cex=0.8)
text(2, 0.5, "Hadri test p-value", cex=0.8)
text(2, 0.35, "shown in corner", cex=0.8)
text(2, 0.2, "NULL: No unit root", cex=0.8)

#4:length(colnames(PAD)) &

for(i in colnames(PAD)[!grepl("x_", colnames(PAD))][-(1:3)]){
    ICY$unitRootp[ICY$ic == i] <- round(
        purtest(PAD[i], data=PAD,
        index = c("ccx", "year"),
        test= "hadri", 
        lags = "AIC",
        exo = "intercept")$statistic$p.value,
        2)
    
    dPAD[,i] <- diff(PAD[,i])
    
    for(j in 1:length(CCX$ccx)){
        if(j==1){
            plot(PAD[PAD$ccx==CCX$ccx[1],i], 
                     type="l", 
                     pch=20, 
                     ylim=c(min(na.omit(PAD[,i])), 
                            max(na.omit(PAD[,i]))), 
                     main=i,
                     axes=FALSE,
                     frame.plot=TRUE)
        } else {
            lines(PAD[PAD$ccx==CCX$ccx[j],i], type="l", pch=20, lty=1, col=paste0("grey",2*j))
        }
        
    }
    legend("topleft", legend=ICY$unitRootp[ICY$ic == i],
           x.intersp=0, y.intersp=0, 
           bg=ifelse(ICY$unitRootp[ICY$ic == i]>0.05, "palegreen", "snow"))
}

par(mar=c(0,0,1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(0,0,0,0),
    cex.axis=0.5, #reduce size of axis numbers 
    cex.main=0.8, #Increase title size
    mfrow=c(4,8))

plot(c(1,2,3), c(0,0,0),
         type="l", 
         ylim=c(0,1),
         main="First difference plot",
         axes=FALSE,
         frame.plot=TRUE)

text(1, 0.95, "x: time (1985 - 2015)", 0, cex=0.8)
text(1, 0.80, "y: variable (relative)", 0, cex=0.8)
text(2, 0.5, "Hadri test p-value", cex=0.8)
text(2, 0.35, "shown in corner", cex=0.8)
text(2, 0.2, "NULL: No unit root", cex=0.8)

for(i in colnames(dPAD)[!grepl("x_", colnames(dPAD))][-(1:3)]){
    
    ICY$unitRootdp[ICY$ic == i] <- round(
        purtest(dPAD[i], data=dPAD,
        index = c("ccx", "year"),
        test= "hadri", 
        lags = "AIC",
        exo = "intercept")$statistic$p.value,
        2)
    
    for(j in 1:length(CCX$ccx)){
        if(j==1){
            plot(dPAD[dPAD$ccx==CCX$ccx[1],i], 
                     type="l", 
                     pch=20, 
                     ylim=c(min(na.omit(dPAD[,i])), 
                            max(na.omit(dPAD[,i]))), 
                     main=paste0("d_",i),
                     axes=FALSE,
                     frame.plot=TRUE)
        } else {
            lines(dPAD[dPAD$ccx==CCX$ccx[j],i], type="l", pch=20, lty=1, col=paste0("grey",2*j))
        }
        legend("topleft", legend=ICY$unitRootdp[ICY$ic == i], 
               x.intersp=0, y.intersp=0, 
               bg=ifelse(ICY$unitRootdp[ICY$ic == i]>0.05, "palegreen", "snow"))
        abline(0,0, col="orangered")
    }
}
rm(list=c("dPAD", "j", "i"))

Economic Convergence

The very first question to be considered is whether there is convergence across the region considered.

Beta

GDP per capita

par(mar=c(0.1,0,0.1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(3,0,1,0),
    cex.axis=1, #reduce size of axis numbers 
    cex.main=1, #Increase title size
    mfrow=c(4,2))

for(w in 1:length(REG$pairs)){
    r <- plm(WBGDPYOY ~ log(lag(WBGDP))*year - log(lag(WBGDP)),
                data=PAD[PAD$ccx %in% CCX$ccx[CCX$cg %in% do.call(paste, c(as.list(REG$pairs[w]), sep = ""))],],
                model="within", effect="time", 
                index=c("ccx", "year"), na.action = na.exclude)
    Braw <- summary(r)$coefficients[,1]
    names(Braw) <- NULL
    Bmis <- 26-length(Braw)
    Bplt <- c(rep(NA, Bmis), Braw)
    
    eBraw <- summary(r)$coefficients[,2]
    names(eBraw) <- NULL
    eBmis <- 26-length(eBraw)
    eBplt <- c(rep(NA, eBmis), eBraw)
    
    plot(Bplt, type="l", ylim=c(-10, 10), xlim=c(3,25), lwd=2, axes=FALSE)
    lines(Bplt+2*eBplt, type="l", lty=3, col="grey30")
    lines(Bplt-2*eBplt, type="l", lty=3, col="grey30")
    abline(0,0, col="orangered")
    text(26, 8, "D", pos=2)
    text(26, -8, "C", pos=2)
    text(3, -8, gsub("[^A-Z]", "", paste0(as.character(REG$pairs[w]))), cex=1.25, pos=4)
    rug(which(!(Bplt-2*eBplt)<0), ticksize=0.1, lwd=2, col="grey60", side=3)
    rug(which(!(Bplt+2*eBplt)>0), ticksize=0.1, lwd=2, col="grey60")
    box()
    if(w %in% c(length(REG$pairs)-1, length(REG$pairs))){axis(1, las=2, labels=seq(1992, 2014, by=2), at=seq(3, 25, by=2))}
}
title(main="Beta convergence - GDP per capita", outer=TRUE)

rm(list=c("w", "r", "Braw", "Bmis", "Bplt", "eBraw", "eBmis", "eBplt"))

Betasim GDP

#Run average convergence model
betaGDPmodel <- plm(WBGDPYOY ~ log(lag(WBGDP))
                                                #Labor Force
            + lag(diff(EDU3er))
            + lag(diff(DEPRATIO))
            # + lag(diff(UNEMPLOYMENT))
            + lag(POPg)
            #+ lag(diff(EDU2er))
        #Physical Capital
            + lag(diff(GKF))
            # + lag(PATpcc) #patents per 10k citizens
        #International Trade
            + lag(EXPORTS-IMPORTS)
        #Fiscal and monetary climate
            + lag(diff(INFLATION))
            + lag(GOVT)
            # + lag(GGXWDG_GDP)
        #International Organizations
            + x_EU
            # + x_EURO
        , data=PAD,
            model="pooling", 
            index=c("ccx", "year"), na.action = na.exclude)
summary(betaGDPmodel)

## Pooling Model
## 
## Call:
## plm(formula = WBGDPYOY ~ log(lag(WBGDP)) + lag(diff(EDU3er)) + 
##     lag(diff(DEPRATIO)) + lag(POPg) + lag(diff(GKF)) + lag(EXPORTS - 
##     IMPORTS) + lag(diff(INFLATION)) + lag(GOVT) + x_EU, data = PAD, 
##     na.action = na.exclude, model = "pooling", index = c("ccx", 
##         "year"))
## 
## Unbalanced Panel: n = 35, T = 6-24, N = 644
## 
## Residuals:
##      Min.   1st Qu.    Median   3rd Qu.      Max. 
## -19.52998  -1.22094   0.28629   1.59452  21.12929 
## 
## Coefficients:
##                          Estimate Std. Error t-value  Pr(>|t|)    
## (Intercept)            16.5483838  2.5335409  6.5317 1.335e-10 ***
## log(lag(WBGDP))        -1.3546298  0.2599798 -5.2105 2.550e-07 ***
## lag(diff(EDU3er))       0.2234885  0.0536580  4.1651 3.544e-05 ***
## lag(diff(DEPRATIO))    -1.3728329  0.2654007 -5.1727 3.099e-07 ***
## lag(POPg)              -0.5598190  0.2474709 -2.2622  0.024025 *  
## lag(diff(GKF))          0.4391534  0.0509959  8.6115 < 2.2e-16 ***
## lag(EXPORTS - IMPORTS)  0.1108287  0.0197637  5.6077 3.063e-08 ***
## lag(diff(INFLATION))   -0.0052247  0.0016505 -3.1655  0.001622 ** 
## lag(GOVT)              -0.0376382  0.0129889 -2.8977  0.003889 ** 
## x_EU                    0.7383157  0.3250358  2.2715  0.023452 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Total Sum of Squares:    9876.2
## Residual Sum of Squares: 6965.7
## R-Squared:      0.2947
## Adj. R-Squared: 0.28468
## F-statistic: 29.4337 on 9 and 634 DF, p-value: < 2.22e-16

alphaGDP <- as.numeric(betaGDPmodel$coefficients)[1]
betaGDP <- as.numeric(betaGDPmodel$coefficients)[2]

#Prepare dataframe
GDPfc <- data.frame(year=1990:2600)
#Fill with historic values
for(r in REG$codes){
    GDPfc[,r] <- 
        c(
        by(PAD$WBGDP[PAD$cg == r],
           PAD$year[PAD$cg == r],
           mean, na.rm=FALSE),
        rep(NA, 585)
        )
}
#Fill with predicted values
for(t in 2016:2600){
    for(r in REG$codes){
    GDPfc[GDPfc$year==t,r] <-
        ((log(GDPfc[GDPfc$year==(t-1),r])*betaGDP+alphaGDP)/100+1)*GDPfc[GDPfc$year==(t-1),r]
}
}

par(mar=c(0,0,1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(2,2,0,0),
    cex.axis=0.8, #reduce size of axis numbers 
    cex.main=1, #Increase title size
    mfrow=c(1,1))
plot(GDPfc$year, GDPfc$WE, type="l", xlim=c(2005, 2150), ylim=c(5000,150000), 
     lwd=2, lty=2, axes=TRUE, main="Dynamic forecast: GDP nu-convergence")
    legend("bottomright", REG$names, 
       col = c("grey80","grey40","black","black"),
       lty = c(4, 3, 1, 2), lwd=c(2, 2, 2, 2),
       ncol=2, cex=0.8, title="GDP per capita 2010$", bty="n")
    lines(GDPfc$year, GDPfc$EB, type="l", lwd=2, lty=3, col="grey40")
    lines(GDPfc$year, GDPfc$BK, type="l", lwd=2, lty=4, col="grey80")
    lines(GDPfc$year, GDPfc$PS, type="l", lwd=2, lty=1, col="black")

rm(list=c("t", "r", "alphaGDP", "betaGDP", "betaGDPmodel", "GDPfc"))

GNI PPP per capita

par(mar=c(0.1,0,0.1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(3,0,1,0),
    cex.axis=1, #reduce size of axis numbers 
    cex.main=1, #Increase title size
    mfrow=c(4,2))

for(w in 1:length(REG$pairs)){
    r <- plm(GNIpcPPPg ~ log(lag(GNIpcPPPg))*year - log(lag(GNIpcPPPg)),
                data=PAD[PAD$ccx %in% CCX$ccx[CCX$cg %in% do.call(paste, c(as.list(REG$pairs[w]), sep = ""))],],
                model="within", effect="time", 
                index=c("ccx", "year"), na.action = na.exclude)
    Braw <- summary(r)$coefficients[,1]
    names(Braw) <- NULL
    Bmis <- 26-length(Braw)
    Bplt <- c(rep(NA, Bmis), Braw)
    
    eBraw <- summary(r)$coefficients[,2]
    names(eBraw) <- NULL
    eBmis <- 26-length(eBraw)
    eBplt <- c(rep(NA, eBmis), eBraw)
    
    plot(Bplt, type="l", ylim=c(-10, 10), xlim=c(3,25), col="grey30", lwd=2, axes=FALSE)
    lines(Bplt+2*eBplt, type="l", lty=3, col="grey40")
    lines(Bplt-2*eBplt, type="l", lty=3, col="grey40")
    abline(0,0, col="orangered")
    text(26, 8, "D", pos=2)
    text(26, -8, "C", pos=2)
    text(3, -8, gsub("[^A-Z]", "", paste0(as.character(REG$pairs[w]))), cex=1.25, pos=4)
    rug(which(!(Bplt-2*eBplt)<0), ticksize=0.1, lwd=2, col="grey60", side=3)
    rug(which(!(Bplt+2*eBplt)>0), ticksize=0.1, lwd=2, col="grey60")
    box()
    if(w %in% c(length(REG$pairs)-1, length(REG$pairs))){axis(1, las=2, labels=seq(1992, 2014, by=2), at=seq(3, 25, by=2))}
}
title(main="Beta convergence - GNI PPP per capita", outer=TRUE)

rm(list=c("w", "r", "Braw", "Bmis", "Bplt", "eBraw", "eBmis", "eBplt"))

Betasim GNI

#Run average convergence model
betaGNImodel <- plm(GNIpcPPPg ~ log(lag(GNIpcPPP))
                                                #Labor Force
            + lag(diff(EDU3er))
            + lag(diff(DEPRATIO))
            # + lag(diff(UNEMPLOYMENT))
            + lag(POPg)
            #+ lag(diff(EDU2er))
        #Physical Capital
            + lag(diff(GKF))
            # + lag(PATpcc) #patents per 10k citizens
        #International Trade
            + lag(EXPORTS-IMPORTS)
        #Fiscal and monetary climate
            + lag(diff(INFLATION))
            + lag(GOVT)
            # + lag(GGXWDG_GDP)
        #International Organizations
            + x_EU
            # + x_EURO
        , data=PAD,
            model="pooling", 
            index=c("ccx", "year"), na.action = na.exclude)
summary(betaGNImodel)

## Pooling Model
## 
## Call:
## plm(formula = GNIpcPPPg ~ log(lag(GNIpcPPP)) + lag(diff(EDU3er)) + 
##     lag(diff(DEPRATIO)) + lag(POPg) + lag(diff(GKF)) + lag(EXPORTS - 
##     IMPORTS) + lag(diff(INFLATION)) + lag(GOVT) + x_EU, data = PAD, 
##     na.action = na.exclude, model = "pooling", index = c("ccx", 
##         "year"))
## 
## Unbalanced Panel: n = 35, T = 6-24, N = 615
## 
## Residuals:
##      Min.   1st Qu.    Median   3rd Qu.      Max. 
## -20.44418  -1.32419   0.21119   1.68156  15.55223 
## 
## Coefficients:
##                          Estimate Std. Error t-value  Pr(>|t|)    
## (Intercept)            27.8064394  4.2480796  6.5456 1.265e-10 ***
## log(lag(GNIpcPPP))     -2.4652565  0.4226821 -5.8324 8.900e-09 ***
## lag(diff(EDU3er))       0.2181769  0.0556859  3.9180 9.948e-05 ***
## lag(diff(DEPRATIO))    -1.2102863  0.2939589 -4.1172 4.368e-05 ***
## lag(POPg)              -0.6682192  0.2412528 -2.7698  0.005781 ** 
## lag(diff(GKF))          0.4609750  0.0532146  8.6626 < 2.2e-16 ***
## lag(EXPORTS - IMPORTS)  0.1328824  0.0204069  6.5116 1.563e-10 ***
## lag(diff(INFLATION))   -0.0044207  0.0016507 -2.6781  0.007606 ** 
## lag(GOVT)              -0.0326411  0.0140601 -2.3215  0.020588 *  
## x_EU                    0.7992050  0.3358602  2.3796  0.017642 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Total Sum of Squares:    9416.4
## Residual Sum of Squares: 6637.6
## R-Squared:      0.29511
## Adj. R-Squared: 0.28462
## F-statistic: 28.143 on 9 and 605 DF, p-value: < 2.22e-16

alphaGNI <- as.numeric(betaGNImodel$coefficients)[1]
betaGNI <- as.numeric(betaGNImodel$coefficients)[2]

#Prepare dataframe
GNIfc <- data.frame(year=1990:2600)
#Fill with historic values
for(r in REG$codes){
    GNIfc[,r] <- 
        c(
        by(PAD$GNIpcPPP[PAD$cg == r],
           PAD$year[PAD$cg == r],
           mean, na.rm=FALSE),
        rep(NA, 585)
        )
}
#Fill with predicted values
for(t in 2016:2600){
    for(r in REG$codes){
    GNIfc[GNIfc$year==t,r] <-
        ((log(GNIfc[GNIfc$year==(t-1),r])*betaGNI+alphaGNI)/100+1)*GNIfc[GNIfc$year==(t-1),r]
}
}

par(mar=c(0,0,1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(2,2,0,0),
    cex.axis=0.8, #reduce size of axis numbers 
    cex.main=1, #Increase title size
    mfrow=c(1,1))
plot(GNIfc$year, GNIfc$WE, type="l", xlim=c(2005, 2150), ylim=c(10000,80000), 
     lwd=2, lty=2, axes=TRUE, main="Dynamic forecast: GNI PPP nu-convergence")
    legend("bottomright", REG$names, 
       col = c("grey80","grey40","black","black"),
       lty = c(4, 3, 1, 2), lwd=c(2, 2, 2, 2),
       ncol=2, cex=0.8, title="GDP per capita 2010$", bty="n")
    lines(GNIfc$year, GNIfc$EB, type="l", lwd=2, lty=3, col="grey40")
    lines(GNIfc$year, GNIfc$BK, type="l", lwd=2, lty=4, col="grey80")
    lines(GNIfc$year, GNIfc$PS, type="l", lwd=2, lty=1, col="black")

rm(list=c("t", "r", "alphaGNI", "betaGNI", "betaGNImodel", "GNIfc"))

Interpreting the output of the very first model, it is obvious that due to the fact that coefficient of the logged lag of GDP is negative, convergence does exist. In other words, the higher GDP country had last year, the lower its real percentage growth tends to be.

Sigma

convS <- data.frame(year=1990:2015)

for(t in 1:length(REG$pairs)){
    convS[,paste0(gsub("[^A-Z]", "", as.character(REG$pairs[t])),"gdp")] <- 
        by(PAD$WBGDP[PAD$ccx %in% CCX$ccx[CCX$cg %in% do.call(paste, c(as.list(REG$pairs[t]), sep = ""))]],
           PAD$year[PAD$ccx %in% CCX$ccx[CCX$cg %in% do.call(paste, c(as.list(REG$pairs[t]), sep = ""))]],
           sd, na.rm=TRUE)
        convS[,paste0(gsub("[^A-Z]", "", as.character(REG$pairs[t])),"gni")] <- 
        by(PAD$GNIpcPPP[PAD$ccx %in% CCX$ccx[CCX$cg %in% do.call(paste, c(as.list(REG$pairs[t]), sep = ""))]],
           PAD$year[PAD$ccx %in% CCX$ccx[CCX$cg %in% do.call(paste, c(as.list(REG$pairs[t]), sep = ""))]],
           sd, na.rm=TRUE)
}

GDP per capita

par(mar=c(0.1,0,0.1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(3,0,1,0),
    cex.axis=1, #reduce size of axis numbers 
    cex.main=1, #Increase title size
    mfrow=c(4,2))

#GDP
for(s in seq(2, (length(colnames(convS))), by=2)){
    plot(convS$year[-1], diff(convS[,(s)]),
       type="l", axes=FALSE, 
       ylim=c(-1500, 1500), lwd=2, lty=1)
    
    rug(1990+which(!(as.vector(diff(convS[,(s)])))>0), ticksize=0.1, lwd=2, col="grey60")
    abline(0,0, col="orangered")
    text(2016, 1200, "D", pos=2)
    text(2016, -1200, "C", pos=2)
    text(1990, -1200, gsub("[^A-Z]", "", colnames(convS)[s]), cex=1.25, pos=4)
    box()
    if(s %in% c(length(colnames(convS))-1, length(colnames(convS))-3)){axis(1, las=2, at=seq(1992, 2015, by=2))}
}

title(main="Sigma convergence - GDP per capita", outer=TRUE)

GNI PPP per capita

par(mar=c(0.1,0,0.1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(3,0,1,0),
    cex.axis=1, #reduce size of axis numbers 
    cex.main=1, #Increase title size
    mfrow=c(4,2))

#GNI
for(s in seq(2, (length(colnames(convS))), by=2)){
    plot(convS$year[-1], diff(convS[,(s+1)]),
       type="l", col="grey40", axes=FALSE, 
       ylim=c(-1500, 1500), lwd=2, lty=1)
    
    rug(1990+which(!(as.vector(diff(convS[,(s+1)])))>0), ticksize=0.1, lwd=2, col="grey60")
    abline(0,0, col="orangered")
    text(2016, 1200, "D", pos=2)
    text(2016, -1200, "C", pos=2)
    text(1990, -1200, gsub("[^A-Z]", "", colnames(convS)[s]), cex=1.25, pos=4)
    box()
    if(s %in% c(length(colnames(convS))-1, length(colnames(convS))-3)){axis(1, las=2, at=seq(1992, 2015, by=2))}
}
title(main="Sigma convergence - GNI PPP per capita", outer=TRUE)

Sigma convergence does not seem to be present. The standard deviation of real GDP per capita continues to increase through time with exception of five cases. This also implies heteroskedasticity in the variable of WBGDP

Endogenous Beta

GDP per capita

par(mar=c(0.1,0,0.1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(3,0,1,0),
    cex.axis=1, #reduce size of axis numbers 
    cex.main=1, #Increase title size
    mfrow=c(4,2))

for(t in 1:length(REG$pairs)){
    r  <- plm(WBGDPYOY ~ lag(log(WBGDP))*year - lag(log(WBGDP))
        #Labor Force
            + lag(diff(EDU3er))
            + lag(diff(DEPRATIO))
            + lag(diff(UNEMPLOYMENT))
            + lag(POPg)
            #+ lag(diff(EDU2er))
        #Physical Capital
            + lag(diff(GKF))
            + lag(PATpcc) #patents per 10k citizens
        #International Trade
            + lag(EXPORTS-IMPORTS)
        #Fiscal and monetary climate
            + lag(diff(INFLATION))
            + lag(GOVT)
            + lag(GGXWDG_GDP)
        #International Organizations
            + x_EU
            + x_EURO
        , data=PAD[PAD$ccx %in% CCX$ccx[CCX$cg %in% do.call(paste, c(as.list(REG$pairs[t]), sep = ""))],],
        model="within",
        effect="time",
        index=c("ccx", "year"), 
        na.action = na.exclude,
        pooling=TRUE)
    #g<-step(g, direction = "backward", trace=0)
    
    Braw <- summary(r)$coefficients[,1][-(1:11)]
    names(Braw) <- NULL
    Bmis <- 26-length(Braw)
    Bplt <- c(rep(NA, Bmis), Braw)
    
    eBraw <- summary(r)$coefficients[,2][-(1:11)]
    names(eBraw) <- NULL
    eBmis <- 26-length(eBraw)
    eBplt <- c(rep(NA, eBmis), eBraw)
    
    plot(Bplt, type="l", ylim=c(-10, 10), xlim=c(3,25), lwd=2, axes=FALSE)
    lines(Bplt+2*eBplt, type="l", lty=3, col="grey30")
    lines(Bplt-2*eBplt, type="l", lty=3, col="grey30")
    abline(0,0, col="orangered")
    text(26, 8, "D", pos=2)
    text(26, -8, "C", pos=2)
    text(3, -8, gsub("[^A-Z]", "", paste0(as.character(REG$pairs[t]))), cex=1.25, pos=4)
    rug(which(!(Bplt-2*eBplt)<0), ticksize=0.1, lwd=2, col="grey60", side=3)
    rug(which(!(Bplt+2*eBplt)>0), ticksize=0.1, lwd=2, col="grey60")
    box()
    if(t %in% c(length(REG$pairs)-1, length(REG$pairs))){axis(1, las=2, labels=seq(1992, 2014, by=2), at=seq(3, 25, by=2))}
    
    
    assign(paste0("GDPc", gsub("[^A-Z]", "", as.character(REG$pairs[t]))), r)
}

title(main="Nu Convergence: GDP per capita", outer=TRUE)

stargazer(GDPcEB, GDPcBK, GDPcPS, GDPcWE, GDPcEBWE, GDPcBKWE, GDPcPSWE, GDPcEBWEBKPS,
          type="html", keep.stat=c("n", "adj.rsq", "rsq"), model.numbers=FALSE,
          omit=(length(coef(GDPcEBWEBKPS))-23):length(coef(GDPcEBWEBKPS)),
          column.labels = c("EB", "BK", "PS", "WE", "EBWE", "BKWE", "PSWE", "All"))


	Dependent variable:

	WBGDPYOY
	EB	BK	PS	WE	EBWE	BKWE	PSWE	All

lag(diff(EDU3er))	0.250^*	1.035	-0.188	0.184^***	0.176^***	0.186^***	0.145^***	0.156^***
	(0.135)	(Inf.000)	(0.162)	(0.050)	(0.043)	(0.046)	(0.051)	(0.042)

lag(diff(DEPRATIO))	-1.719^*	-4.759	-3.243^***	-0.307	-0.860^***	-0.548^**	-1.270^***	-1.121^***
	(0.863)	(Inf.000)	(1.105)	(0.300)	(0.240)	(0.275)	(0.314)	(0.248)

lag(diff(UNEMPLOYMENT))	-0.526^**	-0.307	-0.432	-0.371^***	-0.304^***	-0.328^***	-0.352^***	-0.342^***
	(0.260)	(Inf.000)	(0.290)	(0.118)	(0.090)	(0.100)	(0.097)	(0.080)

lag(POPg)	1.523	0.843	-1.631	0.034	0.020	0.225	-0.845^***	-0.590^***
	(1.191)	(Inf.000)	(1.246)	(0.275)	(0.239)	(0.218)	(0.222)	(0.180)

lag(diff(GKF))	0.308^**	0.194	0.139	0.443^***	0.419^***	0.382^***	0.429^***	0.398^***
	(0.149)	(Inf.000)	(0.149)	(0.084)	(0.061)	(0.072)	(0.068)	(0.053)

lag(PATpcc)	-0.083	-1.988	0.980	-0.061	-0.120	-0.044	-0.283^**	-0.222^**
	(1.711)	(Inf.000)	(1.479)	(0.098)	(0.094)	(0.095)	(0.110)	(0.101)

lag(EXPORTS - IMPORTS)	0.083	-0.077	0.043	0.179^***	0.115^***	0.142^***	0.068^***	0.069^***
	(0.104)	(Inf.000)	(0.056)	(0.027)	(0.018)	(0.021)	(0.017)	(0.015)

lag(diff(INFLATION))	-0.024^**	-0.204	0.001	-0.014	-0.017^**	-0.072	-0.006^*	-0.012^***
	(0.012)	(Inf.000)	(0.010)	(0.092)	(0.007)	(0.079)	(0.003)	(0.002)

lag(GOVT)	0.053	0.550	-0.004	0.011	-0.001	0.007	-0.006	-0.019^*
	(0.077)	(Inf.000)	(0.032)	(0.015)	(0.014)	(0.014)	(0.012)	(0.010)

lag(GGXWDG_GDP)	-0.020	-0.177	-0.067^*	0.001	-0.002	0.002	-0.012^**	-0.010^**
	(0.016)	(Inf.000)	(0.034)	(0.006)	(0.004)	(0.005)	(0.005)	(0.004)

x_EU	-3.264	-6.658	-0.277	-0.272	0.391	-0.034	0.313	0.564^**
	(3.173)	(Inf.000)	(3.737)	(0.362)	(0.269)	(0.284)	(0.334)	(0.258)


Observations	108	44	108	343	451	387	451	603
R²	0.502	1.000	0.580	0.356	0.473	0.390	0.564	0.509
Adjusted R²	-0.065		0.211	0.228	0.396	0.284	0.500	0.458

Note:	p<0.1; p<0.05; p<0.01

##Export to Tex
trash <- capture.output(
            stargazer(GDPcEB, GDPcBK, GDPcPS, GDPcWE, GDPcEBWE, GDPcBKWE, GDPcPSWE, GDPcEBWEBKPS,
                type="latex", keep.stat=c("n", "adj.rsq", "rsq"), model.numbers=FALSE,
                omit=(length(coef(GDPcEBWEBKPS))-23):length(coef(GDPcEBWEBKPS)),
                column.labels = c("EB", "BK", "PS", "WE", "EBWE", "BKWE", "PSWE", "All"),
                out="./TAB/WGDPc.tex"))

GNI PPP per capita

par(mar=c(0.1,0,0.1,0), #margins: L&B 2 for axes, T 0.5 for title
    oma=c(3,0,1,0),
    cex.axis=1, #reduce size of axis numbers 
    cex.main=1, #Increase title size
    mfrow=c(4,2))

for(t in 1:length(REG$pairs)){
    r  <- plm(GNIpcPPPg ~ lag(log(GNIpcPPP))*year - lag(log(GNIpcPPP))
        #Labor Force
            + lag(diff(EDU3er))
            + lag(diff(DEPRATIO))
            + lag(diff(UNEMPLOYMENT))
            + lag(POPg)
            #+ lag(diff(EDU2er))
        #Physical Capital
            + lag(diff(GKF))
            + lag(PATpcc) #patents per 10k citizens
        #International Trade
            + lag(EXPORTS-IMPORTS)
        #Fiscal and monetary climate
            + lag(diff(INFLATION))
            + lag(GOVT)
            + lag(GGXWDG_GDP)
        #International Organizations
            + x_EU
            + x_EURO
        , data=PAD[PAD$ccx %in% CCX$ccx[CCX$cg %in% do.call(paste, c(as.list(REG$pairs[t]), sep = ""))],],
        model="within",
        effect="time",
        index=c("ccx", "year"), 
        na.action = na.exclude,
        pooling=TRUE)
    #g<-step(g, direction = "backward", trace=0)
    
    Braw <- summary(r)$coefficients[,1][-(1:11)]
    names(Braw) <- NULL
    Bmis <- 26-length(Braw)
    Bplt <- c(rep(NA, Bmis), Braw)
    
    eBraw <- summary(r)$coefficients[,2][-(1:11)]
    names(eBraw) <- NULL
    eBmis <- 26-length(eBraw)
    eBplt <- c(rep(NA, eBmis), eBraw)
    
    plot(Bplt, type="l", ylim=c(-10, 10), xlim=c(3,25), lwd=2, axes=FALSE, col="grey30")
    lines(Bplt+2*eBplt, type="l", lty=3, col="grey40")
    lines(Bplt-2*eBplt, type="l", lty=3, col="grey40")
    abline(0,0, col="orangered")
    text(26, 8, "D", pos=2)
    text(26, -8, "C", pos=2)
    text(3, -8, gsub("[^A-Z]", "", paste0(as.character(REG$pairs[t]))), cex=1.25, pos=4)
    rug(which(!(Bplt-2*eBplt)<0), ticksize=0.1, lwd=2, col="grey60", side=3)
    rug(which(!(Bplt+2*eBplt)>0), ticksize=0.1, lwd=2, col="grey60")
    box()
    if(t %in% c(length(REG$pairs)-1, length(REG$pairs))){axis(1, las=2, labels=seq(1992, 2014, by=2), at=seq(3, 25, by=2))}
    assign(paste0("GNIc", gsub("[^A-Z]", "", as.character(REG$pairs[t]))), r)
}

title(main="NuConvergence: GNI PPP per capita", outer=TRUE)

stargazer(GNIcEB, GNIcBK, GNIcPS, GNIcWE, GNIcEBWE, GNIcBKWE, GNIcPSWE, GNIcEBWEBKPS,
          type="html", keep.stat=c("n", "adj.rsq", "rsq"), model.numbers=FALSE, 
          omit=(length(coef(GDPcEBWEBKPS))-23):length(coef(GDPcEBWEBKPS)),
          column.labels = c("EB", "BK", "PS", "WE", "EBWE", "BKWE", "PSWE", "All"))


	Dependent variable:

	GNIpcPPPg
	EB	BK	PS	WE	EBWE	BKWE	PSWE	All

lag(diff(EDU3er))	0.389^**	-0.213	-0.233	0.134^**	0.181^***	0.171^***	0.126^**	0.153^***
	(0.150)	(Inf.000)	(0.207)	(0.053)	(0.048)	(0.052)	(0.056)	(0.046)

lag(diff(DEPRATIO))	-2.606^**	-55.837	-3.494^**	0.175	-0.570^*	-0.268	-0.817^**	-0.941^***
	(1.144)	(Inf.000)	(1.543)	(0.341)	(0.307)	(0.323)	(0.368)	(0.294)

lag(diff(UNEMPLOYMENT))	-0.375	0.154	-0.314	-0.157	-0.206^**	-0.158	-0.398^***	-0.371^***
	(0.312)	(Inf.000)	(0.330)	(0.125)	(0.104)	(0.112)	(0.104)	(0.088)

lag(POPg)	2.256	0.837	-1.140	-0.547^*	-0.066	-0.079	-0.876^***	-0.579^***
	(1.447)	(Inf.000)	(1.361)	(0.289)	(0.253)	(0.230)	(0.227)	(0.188)

lag(diff(GKF))	0.104	1.864	0.289^*	0.795^***	0.594^***	0.633^***	0.559^***	0.492^***
	(0.200)	(Inf.000)	(0.166)	(0.088)	(0.074)	(0.081)	(0.073)	(0.060)

lag(PATpcc)	-1.939	-4.011	0.022	-0.145	-0.026	-0.031	-0.180	-0.135
	(2.770)	(Inf.000)	(1.543)	(0.099)	(0.097)	(0.099)	(0.111)	(0.104)

lag(EXPORTS - IMPORTS)	-0.098	3.751	0.031	0.231^***	0.153^***	0.180^***	0.094^***	0.094^***
	(0.159)	(Inf.000)	(0.068)	(0.024)	(0.021)	(0.022)	(0.019)	(0.017)

lag(diff(INFLATION))	-0.045^*	-0.820	0.010	0.074	-0.021^***	0.039	-0.006^**	-0.012^***
	(0.026)	(Inf.000)	(0.012)	(0.096)	(0.008)	(0.089)	(0.003)	(0.002)

lag(GOVT)	0.114	3.324	0.032	0.003	-0.009	-0.006	-0.011	-0.016
	(0.206)	(Inf.000)	(0.041)	(0.016)	(0.016)	(0.016)	(0.013)	(0.012)

lag(GGXWDG_GDP)	-0.007	-1.372	-0.077	-0.001	0.002	0.006	-0.005	-0.006
	(0.056)	(Inf.000)	(0.049)	(0.006)	(0.004)	(0.005)	(0.005)	(0.004)

x_EU	-1.755	-36.308	-0.866	-0.589^*	0.306	-0.030	0.279	0.717^**
	(3.112)	(Inf.000)	(3.188)	(0.354)	(0.305)	(0.315)	(0.353)	(0.278)


Observations	98	42	96	340	438	382	436	576
R²	0.533	1.000	0.659	0.484	0.486	0.459	0.570	0.511
Adjusted R²	-0.132		0.280	0.380	0.409	0.363	0.505	0.457

Note:	p<0.1; p<0.05; p<0.01

##Export to Tex

trash <- capture.output(
            stargazer(GNIcEB, GNIcBK, GNIcPS, GNIcWE, GNIcEBWE, GNIcBKWE, GNIcPSWE, GNIcEBWEBKPS,
                type="latex", keep.stat=c("n", "adj.rsq", "rsq"), model.numbers=FALSE,
                omit=(length(coef(GDPcEBWEBKPS))-23):length(coef(GDPcEBWEBKPS)),
                column.labels = c("EB", "BK", "PS", "WE", "EBWE", "BKWE", "PSWE", "All"),
                out="./TAB/WGNIc.tex"))

rm(list=c("GNIcEB", "GNIcBK", "GNIcPS", "GNIcWE", "GNIcEBWE", "GNIcBKWE", "GNIcPSWE", "GNIcEBWEBKPS",
          "GDPcEB", "GDPcBK", "GDPcPS", "GDPcWE", "GDPcEBWE", "GDPcBKWE", "GDPcPSWE", "GDPcEBWEBKPS",
          "r", "t", "trash"))
warnings()

## NULL

Economic Convergence Between Post-Communist Nations and Western Europe

Matěj Adámek

20 April 2018

Introduction

Disclaimer

Workspace preparation

Variable preparation

Macroeconomic Indicators

IMF

WorldBank

Miscellaneous

Dummy variables

Data consolidation

Transformation

Calculation of secondary variables

Integrity check

Data inspection

Plots

Dependent variables

GDP Growth and HDI

Independent Variables

Economic Convergence

Beta

GDP per capita

Betasim GDP

GNI PPP per capita

Betasim GNI

Sigma

GDP per capita

GNI PPP per capita

Endogenous Beta

GDP per capita

GNI PPP per capita