Preamble

This paper seeks to characterize the evolving structure of postwar global capitalist society using

The PWT contains internationally-comparable annual data on real GDP output-method \(Y\), real GDP expenditure-method \(Y_e\), population \(N\), employment \(E\), average hours of labor per employed person \(l\), labor share of income \(\omega\), consumption share of income \(c\), end-of-period capital stock \(K\), capital services \(K^s\), and the internal rate of return of capital \(\tilde{r}\) for 168 countries. The period under consideration ranges from 1954 to 2017 with an emphasis on the 1971-2017 subperiod.

The data are interpreted naively, i.e. under the strong assumption that the variables in the data set match exactly the corresponding Marxist categories. Admittedly, this assumption implies a measurement error of unknown size. With the resources at hand, the cost of transforming the data set of these dimensions to make its variables conform to the theoretical categories seems to be much higher than its potential benefits.

The mapping of the primary theoretical and empirical categories follows:

Symbol Marxist variable PWT variable
\(N\) population population
\(E\) employed workers engaged persons
\(L\) direct labor time hours worked by engaged persons
\(Y\) value added real GDP output method
\(Y_e\) value added real GDP expenditure method
\(W\) variable capital labor income
\(C\) variable capital consumption spending
\(K_s\) constant capital capital services
\(K\) capital stock fixed capital stock

The largest discrepancies between Marxist and empirical categories are related to (1) the likely overestimation of the annual value added as a host of ostensibly nonproductive activities are deemed to add value, and (2) the likely underestimation of the labor share of income (value added), as national accounting conventions include the hefty salaries and labor benefits of managers as labor income.

However, much of the focus of this paper is on the empirical time path of these categories. The time-series characteriztion of the structural variables of interest (productivity, exploitation, and profitability measures) is affected by these discrepancies only inasmuch as such discrepancies grow or decay significantly over time. As long as the magnitude of the discrepancies between theoretical and empirical categories remains relatively constant, the time-series implications of the analysis are barely affected.

To allow full replicability of results, data sets and R code are stored in this repository:

https://github.com/jhuato/deeper_structures

library(readr)
x  <- read_csv("C:/Users/jhuat/Downloads/pwt91.csv")
## Parsed with column specification:
## cols(
##   .default = col_double(),
##   countrycode = col_character(),
##   country = col_character(),
##   currency_unit = col_character(),
##   i_cig = col_character(),
##   i_xm = col_character(),
##   i_xr = col_character(),
##   i_outlier = col_character(),
##   i_irr = col_character()
## )
## See spec(...) for full column specifications.
x <- x[ which(x$year > 1953 ), ]
x$c <- x$csh_c
x$c[x$c<0 | x$c>1] <- NA 
x$Y <- x$rgdpo 
x$Ye <- x$rgdpe
x$om <- x$labsh
x$l <- x$avh
x$Ks <- x$rkna
x$W <- x$labsh * x$Y
x$C <- x$c * x$Ye
x$K <- x$rnna
x$N <- x$pop
x$E <- x$emp
x$L <- x$avh * x$E
x <- x[ , c(1, 4, 53:64)]
str(x)
## Classes 'tbl_df', 'tbl' and 'data.frame':    11648 obs. of  14 variables:
##  $ countrycode: chr  "ABW" "ABW" "ABW" "ABW" ...
##  $ year       : num  1954 1955 1956 1957 1958 ...
##  $ c          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ Y          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ Ye         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ om         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ l          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ Ks         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ W          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ C          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ K          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ N          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ E          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ L          : num  NA NA NA NA NA NA NA NA NA NA ...
summary(x)
##  countrycode             year            c                Y           
##  Length:11648       Min.   :1954   Min.   :0.0256   Min.   :      20  
##  Class :character   1st Qu.:1970   1st Qu.:0.5312   1st Qu.:    6338  
##  Mode  :character   Median :1986   Median :0.6294   Median :   27156  
##                     Mean   :1986   Mean   :0.6233   Mean   :  273105  
##                     3rd Qu.:2001   3rd Qu.:0.7360   3rd Qu.:  139871  
##                     Max.   :2017   Max.   :0.9998   Max.   :18383838  
##                                    NA's   :2153     NA's   :1902      
##        Ye                 om              l              Ks       
##  Min.   :      18   Min.   :0.090   Min.   :1354   Min.   :0.008  
##  1st Qu.:    6157   1st Qu.:0.445   1st Qu.:1791   1st Qu.:0.210  
##  Median :   27255   Median :0.545   Median :1963   Median :0.478  
##  Mean   :  276014   Mean   :0.533   Mean   :1978   Mean   :0.534  
##  3rd Qu.:  140987   3rd Qu.:0.624   3rd Qu.:2141   3rd Qu.:0.829  
##  Max.   :18396068   Max.   :0.900   Max.   :2911   Max.   :3.034  
##  NA's   :1902       NA's   :3832    NA's   :8385   NA's   :4601   
##        W                  C                  K                  N            
##  Min.   :      18   Min.   :      11   Min.   :      16   Min.   :   0.0044  
##  1st Qu.:    4928   1st Qu.:    4206   1st Qu.:   20010   1st Qu.:   1.5902  
##  Median :   19159   Median :   17247   Median :  103732   Median :   6.0550  
##  Mean   :  195376   Mean   :  163964   Mean   : 1145903   Mean   :  30.8136  
##  3rd Qu.:  101824   3rd Qu.:   80432   3rd Qu.:  615449   3rd Qu.:  19.7770  
##  Max.   :10720879   Max.   :13076158   Max.   :94903728   Max.   :1409.5175  
##  NA's   :4267       NA's   :2153       NA's   :1928       NA's   :1902       
##        E                  L          
##  Min.   :  0.0012   Min.   :     89  
##  1st Qu.:  0.9321   1st Qu.:   5387  
##  Median :  3.0132   Median :  12094  
##  Mean   : 14.8520   Mean   :  62007  
##  3rd Qu.:  8.5834   3rd Qu.:  44736  
##  Max.   :792.5753   Max.   :1723336  
##  NA's   :3023       NA's   :8385

The variables (hereon called ``parameters’’) that capture the deeper structures of a capitalist society are, respectively,

  1. Labor productivity, measured by
  1. Exploitation or balance of class forces, measured by
  1. Profitability of capital, measured by

This computes the annual parameters for each country:

x$yN <- x$Y/x$N
x$yE <- x$Y/x$E
x$yL <- x$Y/x$L
x$wN <- (1 - x$om)*x$Y/x$N
x$wE <- (1 - x$om)*x$Y/x$E
x$wL <- (1 - x$om)*x$Y/x$L
x$s <- (1 - x$om)/x$om
x$ss <- (1 - x$c)/x$c
x$h <- x$Ks/x$W
x$rho <- x$Y/x$K
x$r <- (1 - x$om)*x$Y/x$K
x$ka <- (1 - x$c)*x$Y/x$K
str(x)
## Classes 'tbl_df', 'tbl' and 'data.frame':    11648 obs. of  26 variables:
##  $ countrycode: chr  "ABW" "ABW" "ABW" "ABW" ...
##  $ year       : num  1954 1955 1956 1957 1958 ...
##  $ c          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ Y          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ Ye         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ om         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ l          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ Ks         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ W          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ C          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ K          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ N          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ E          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ L          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ yN         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ yE         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ yL         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ wN         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ wE         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ wL         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ s          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ ss         : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ h          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ rho        : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ r          : num  NA NA NA NA NA NA NA NA NA NA ...
##  $ ka         : num  NA NA NA NA NA NA NA NA NA NA ...

These structures are defined globally, by agggregation of

N <- tapply(x$N, x$year, sum, na.rm=TRUE) 
E <- tapply(x$E, x$year, sum, na.rm=TRUE)
L <- tapply(x$L, x$year, sum, na.rm=TRUE) 
Y <- tapply(x$Y, x$year, sum, na.rm=TRUE) 
Ye <- tapply(x$Ye, x$year, sum, na.rm=TRUE) 
W <- tapply(x$W, x$year, sum, na.rm=TRUE) 
C <- tapply(x$C, x$year, sum, na.rm=TRUE) 
K <- tapply(x$K, x$year, sum, na.rm=TRUE)
Ks <- tapply(x$Ks, x$year, sum, na.rm=TRUE)
yN <- Y/N
yE <- Y/E
yL <- Y/L
wN <- W/N
wE <- W/E
wL <- W/L
s <- (Y-W)/W
ss <- (Y - C)/C
h <- Ks/W
rho <- Y/K
r <- (Y-W)/K
ka <- (Y-C)/K
year <- c(1954:2017)
g <- data.frame(year, N, E, L, Y, Ye, W, C, K, yN, yE, yL, wN, wE, wL, 
                s, ss, h, rho, r, ka)
str(g)
## 'data.frame':    64 obs. of  21 variables:
##  $ year: int  1954 1955 1956 1957 1958 1959 1960 1961 1962 1963 ...
##  $ N   : num  2047 2120 2161 2205 2249 ...
##  $ E   : num  846 874 892 910 935 ...
##  $ L   : num  674430 688711 699984 706030 704023 ...
##  $ Y   : num  7574162 8217579 8535440 8851155 8996442 ...
##  $ Ye  : num  7667910 8314487 8643475 8952750 9124558 ...
##  $ W   : num  4620027 4954604 5190471 5366630 5436475 ...
##  $ C   : num  4961254 5339159 5528277 5722745 5824452 ...
##  $ K   : num  28542956 30570574 32027560 33548059 35012958 ...
##  $ yN  : num  3699 3875 3949 4014 4001 ...
##  $ yE  : num  8955 9401 9570 9731 9617 ...
##  $ yL  : num  11.2 11.9 12.2 12.5 12.8 ...
##  $ wN  : num  2257 2337 2401 2434 2418 ...
##  $ wE  : num  5462 5668 5819 5900 5812 ...
##  $ wL  : num  6.85 7.19 7.42 7.6 7.72 ...
##  $ s   : num  0.639 0.659 0.644 0.649 0.655 ...
##  $ ss  : num  0.527 0.539 0.544 0.547 0.545 ...
##  $ h   : num  9.61e-07 1.01e-06 1.02e-06 1.11e-06 1.29e-06 ...
##  $ rho : num  0.265 0.269 0.267 0.264 0.257 ...
##  $ r   : num  0.103 0.107 0.104 0.104 0.102 ...
##  $ ka  : num  0.0915 0.0942 0.0939 0.0933 0.0906 ...

Also, these structures are defined internationally by taking account, in due turn, of each country’s

The analysis seeks to empirically characterize

The last section reports on the results of VARs, SVARs, VECM, and SVECM estimations, in an attempt to capture structural relations among

Global results

The following plots show the different measures of labor productivity (per capita, per worker, and per hour) in contrast with the corresponding measures of the real wage (per capita, per worker, and per hour):

byN <- coef(lm(log(g$yN) ~ g$year))[2]
bwN <- coef(lm(log(g$wN) ~ g$year))[2]
plot(g$year, g$yN, type="l", ylim=c(2000, 15000), 
     main="Global per-capita productivity and wage", 
     xlab="Year", ylab="2011 USD")
lines(g$year, g$wN, lty=2)
mtext(bquote(hat(y[N]) ==.(round(byN, 4))~","~ hat(w[N]) ==.(round(bwN, 4))))
legend("topleft", bty="n", lty=c(1, 2), c(expression(paste(y[N]), paste(w[N]))))

byE <- coef(lm(log(g$yE) ~ g$year))[2]
bwE <- coef(lm(log(g$wE) ~ g$year))[2]
plot(g$year, g$yE, type="l", ylim=c(5000, 35000), 
     main="Global per-worker productivity and wage", 
     xlab="Year", ylab="2011 USD")
lines(g$year, g$wE, lty=2)
mtext(bquote(hat(y[E]) ==.(round(byE, 4))~","~ hat(w[E]) ==.(round(bwE, 4))))
legend("topleft", bty="n", lty=c(1, 2), c(expression(paste(y[E]), paste(w[E]))))

byL <- coef(lm(log(g$yL[-c(1:16)]) ~ g$year[-c(1:16)]))[2]
bwL <- coef(lm(log(g$wL[-c(1:16)]) ~ g$year[-c(1:16)]))[2]
plot(g$year[-c(1:16)], g$yL[-c(1:16)], type="l", ylim=c(0, 20), 
     main="Global per-hour productivity and wage", 
     xlab="Year", ylab="2011 USD")
lines(g$year[-c(1:16)], g$wL[-c(1:16)], lty=2)
mtext(bquote(hat(y[L])  ==.(round(byL, 4))~","~ hat(w[L])  ==.(round(bwL, 4))))
legend("topleft", bty="n", lty=c(1, 2), c(expression(paste(y[L]), paste(w[L]))))

These plots show, respectively, the exploitation and expenditure-exploitation rates:

bs <- coef(lm(log(g$s) ~ g$year))[2]
bss <- coef(lm(log(g$ss) ~ g$year))[2]
plot(g$year, g$s, type="l", ylim=c(0.4, 1),
     main="Global exploitation rates", 
     xlab="Year", ylab=" ")
lines(g$year, g$ss, lty=2)
mtext(bquote(hat(sigma)  ==.(round(bs, 4))~","~ hat(sigma[e]) ==.(round(bss, 4))))
legend("topleft", bty="n", lty=c(1, 2), c(expression(paste(sigma), paste(sigma[e] ))))

The flow-flow capital consumption plot follows:

bh <- coef(lm(log(g$h) ~ g$year))[2]
plot(g$year, g$h, type="l", 
     main="Global flow-flow capital composition", 
     xlab="Year", ylab="h")
mtext(bquote(hat(h) ==.(round(bh, 4))))

The output-capital ratio plot follows:

brho <- coef(lm(log(g$rho) ~ g$year))[2]
plot(g$year, g$rho, type="l", 
     main="Global output-capital ratio", 
     xlab="Year", ylab=expression(paste(rho)))
mtext(bquote(hat(rho) ==.(round(brho, 4))))

The following plots show the profit and accumulation rates:

br <- coef(lm(log(g$r) ~ g$year))[2]
bka <- coef(lm(log(g$ka) ~ g$year))[2]
plot(g$year, g$r, type="l", ylim=c(0.08, .13),
     main="Global profit rates", 
     xlab="Year", ylab=" ")
lines(g$year, g$ka, lty=2)
mtext(bquote(hat(r)  ==.(round(br, 4))~","~ hat(kappa) ==.(round(bka, 4))))
legend("topleft", bty="n", lty=c(1, 2), c(expression(paste(r), paste(kappa))))

International results

The following density plots show the international distribution of the structural parameters in selected years: 1975, 1995, and 2015:

yN75<-log(na.omit(subset(x$yN, x$year==1975)))
yN95<-log(na.omit(subset(x$yN, x$year==1995)))
yN15<-log(na.omit(subset(x$yN, x$year==2015)))

yE75<-log(na.omit(subset(x$yE, x$year==1975)))
yE95<-log(na.omit(subset(x$yE, x$year==1995)))
yE15<-log(na.omit(subset(x$yE, x$year==2015)))

yL75<-log(na.omit(subset(x$yL, x$year==1975)))
yL95<-log(na.omit(subset(x$yL, x$year==1995)))
yL15<-log(na.omit(subset(x$yL, x$year==2015)))

wN75<-log(na.omit(subset(x$wN, x$year==1975)))
wN95<-log(na.omit(subset(x$wN, x$year==1995)))
wN15<-log(na.omit(subset(x$wN, x$year==2015)))

wE75<-log(na.omit(subset(x$wE, x$year==1975)))
wE95<-log(na.omit(subset(x$wE, x$year==1995)))
wE15<-log(na.omit(subset(x$wE, x$year==2015)))

wL75<-log(na.omit(subset(x$wL, x$year==1975)))
wL95<-log(na.omit(subset(x$wL, x$year==1995)))
wL15<-log(na.omit(subset(x$wL, x$year==2015)))

s75<-log(na.omit(subset(x$s, x$year==1975)))
s95<-log(na.omit(subset(x$s, x$year==1995)))
s15<-log(na.omit(subset(x$s, x$year==2015)))

ss75<-log(na.omit(subset(x$ss, x$year==1975)))
ss95<-log(na.omit(subset(x$ss, x$year==1995)))
ss15<-log(na.omit(subset(x$ss, x$year==2015)))

h75<-log(na.omit(subset(x$h, x$year==1975)))
h95<-log(na.omit(subset(x$h, x$year==1995)))
h15<-log(na.omit(subset(x$h, x$year==2015)))

rho75<-log(na.omit(subset(x$rho, x$year==1975)))
rho95<-log(na.omit(subset(x$rho, x$year==1995)))
rho15<-log(na.omit(subset(x$rho, x$year==2015)))

r75<-log(na.omit(subset(x$r, x$year==1975)))
r95<-log(na.omit(subset(x$r, x$year==1995)))
r15<-log(na.omit(subset(x$r, x$year==2015)))

ka75<-log(na.omit(subset(x$ka, x$year==1975)))
ka95<-log(na.omit(subset(x$ka, x$year==1995)))
ka15<-log(na.omit(subset(x$ka, x$year==2015)))

# Multi plot function:
# The input of the function MUST be a numeric list
plot.multi.dens <- function(s)
{
        junk.x = NULL
        junk.y = NULL
        for(i in 1:length(s))
        {
                junk.x = c(junk.x, density(s[[i]])$x)
                junk.y = c(junk.y, density(s[[i]])$y)
        }
        xr <- range(junk.x)
        yr <- range(junk.y)
        plot(density(s[[1]]), xlim = xr, ylim = yr, main = "", xlab = "")
        for(i in 1:length(s))
        {
                lines(density(s[[i]]), xlim = xr, ylim = yr, lty=i)
        }
}

library(Hmisc)
## Warning: package 'Hmisc' was built under R version 3.6.3
## Loading required package: lattice
## Loading required package: survival
## Warning: package 'survival' was built under R version 3.6.3
## Loading required package: Formula
## Loading required package: ggplot2
## Warning: package 'ggplot2' was built under R version 3.6.3
## 
## Attaching package: 'Hmisc'
## The following objects are masked from 'package:base':
## 
##     format.pval, units
# Density plots:
plot.multi.dens(list(yN75,yN95,yN15))
title(main = expression(paste(log(y[N]), ": Per-capita productivity")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(yE75,yE95,yE15))
title(main = expression(paste(log(y[E]), ": Per-worker productivity")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(yL75,yL95,yL15))
title(main = expression(paste(log(y[L]), ": Per-hour productivity")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(wN75,wN95,wN15))
title(main = expression(paste(log(w[N]), ": Per-capita wage")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(wE75,wE95,wE15))
title(main = expression(paste(log(w[E]), ": Per-worker wage")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(wL75,wL95,wL15))
title(main = expression(paste(log(w[L]), ": Per-hour wage")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(s75,s95,s15))
title(main = expression(paste(log(sigma), ": Exploitation rate")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(ss75,ss95,ss15))
title(main = expression(paste(log(sigma[e]), ": Expenditure-exploitation rate")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(h75,h95,h15))
title(main = expression(paste(log(h), ": Capital composition")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(rho75,rho95,rho15))
title(main = expression(paste(log(rho), ": Output-capital ratio")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(r75,r95,r15))
title(main = expression(paste(log(r), ": Profit rate")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

plot.multi.dens(list(ka75,ka95,ka15))
title(main = expression(paste(log(kappa), ": Accumulation rate")))
legend("topleft", bty="n", legend=c("1975","1995","2015"), lty =(1:3))

Unweighted aggregation: one country, one observation unit

Under current international law, national states are endowed with formal political sovereignity. History shows that the weight of the juridical figure is relative in terms of effective international relations. As a first approximation to the characterization of the deep international structures of postwar global capitalist society, one must consider the time-series features of plain unweighted aggregates (means, standard errors of means, medians, etc.) of the primary variables (those with people and value flows as measurement units), and their result structural ratios.

The following plots show the evolution of these unweighted aggregates, as well as their resulting structural parameters, over the 1971-2017 period. Let a primary variable for a country \(i\) and a year \(t\) be denoted as \(X_{it}\). Then the mean of this variable across all countries in year \(t\) is denoted as \(\bar{X}_t.\) The time mean of this variable is then denoted as \(\bar{X}^T \equiv \sum_{i=1}^T\) where \(T\) is the length of the time-series sample. Similarly, the cross-country median of the variable is denoted as \(\tilde{X}_t\) and its time mean as \(\tilde{X}^T\). In general, the superscript \(T\) denotes a time mean of the corresponding variable. The geometric mean of the annual growth rate of any variable \(X\) over the period is denoted as \(\hat{X}\). Confidence intervals around the cross-country mean of a variable (in year \(t\), not indicated) are constructed as \([\bar{X} \pm 2 \ \text{se}(\bar{X})]\), where \(\text{se}(\bar{X}) = \bar{X}/\sqrt{n}\), where \(n\) is the length of the country sample in that particular year.

library(plotrix)
## Warning: package 'plotrix' was built under R version 3.6.3
iyN <- tapply(x$yN, x$year, mean, na.rm=TRUE)
summary(iyN)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    4254   10089   13441   12187   14547   19152
biyN <- coef(lm(log(iyN[c(20:64)]) ~ year[c(20:64)]))[2]
seiyN <- tapply(x$yN, x$year, std.error, na.rm=TRUE)
miyN <- tapply(x$yN, x$year, median, na.rm=TRUE)
ciH <- iyN + 2*seiyN
ciL <- iyN - 2*seiyN
plot(year[c(20:64)], iyN[c(20:64)], type="l", lty=0, 
     main=expression(paste("Unweighted mean (95% CI) and median of ", y[N])), 
     xlab="Year", ylab="2011 USD/capita", ylim=c(4000, 25000))
polygon(c(year[c(20:64)], rev(year[c(20:64)])), 
        c(ciL[c(20:64)], rev(ciH[c(20:64)])), 
        col="gray90", border=NA)
lines(year[c(20:64)], iyN[c(20:64)], lty=1)
lines(year[c(20:64)], miyN[c(20:64)], lty=2)
legend("top", bty="n", 
       legend=c(expression(paste(bar(y[N]))), expression(paste(tilde(y[N])))), 
       lty =c(1,2))
mtext(bquote(bar(y[N])^T  == .(round(mean(iyN[c(20:64)]), 1))~ "
             "~se(bar(y[N]))^T ==.(round(mean(seiyN[c(20:64)]), 1))~" 
             "~tilde(y[N])^T  == .(round(mean(miyN[c(20:64)]), 1))~" 
             "~hat(bar(y[N]))  ==.(round(biyN, 3))))

iyE <- tapply(x$yE, x$year, mean, na.rm=TRUE)
summary(iyE)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   12623   22748   29410   27679   32053   42156
biyE <- coef(lm(log(iyE[c(20:64)]) ~ year[c(20:64)]))[2]
seiyE <- tapply(x$yE, x$year, std.error, na.rm=TRUE)
miyE <- tapply(x$yE, x$year, median, na.rm=TRUE)
ciH <- iyE + 2*seiyE
ciL <- iyE - 2*seiyE
plot(year[c(20:64)], iyE[c(20:64)], type="l", lty=0, 
     main=expression(paste("Unweighted mean (95% CI) of ", y[E])), 
     xlab="Year", ylab="2011 USD/worker", ylim=c(0, 55000))
polygon(c(year[c(20:64)], rev(year[c(20:64)])), 
        c(ciL[c(20:64)], rev(ciH[c(20:64)])), 
        col="gray90", border=NA)
lines(year[c(20:64)], iyE[c(20:64)], lty=1)
lines(year[c(20:64)], miyN[c(20:64)], lty=2)
legend("top", bty="n", 
       legend=c(expression(paste(bar(y[E]))), expression(paste(tilde(y[E])))), 
       lty =c(1,2))
mtext(bquote(bar(y[E])^T  == .(round(mean(iyE[c(20:64)]), 1))~ "
             "~se(bar(y[E]))^T ==.(round(mean(seiyE[c(20:64)]), 1))~" 
             "~tilde(y[E])^T  == .(round(mean(miyE[c(20:64)]), 1))~" 
             "~hat(bar(y[E]))  ==.(round(biyE, 3))))

iyL <- tapply(x$yL, x$year, mean, na.rm=TRUE)
summary(iyL)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   7.668  12.762  17.079  18.917  23.592  34.646
biyL <- coef(lm(log(iyL[c(20:64)]) ~ year[c(20:64)]))[2]
seiyL <- tapply(x$yL, x$year, std.error, na.rm=TRUE)
miyL <- tapply(x$yL, x$year, median, na.rm=TRUE)
ciH <- iyL + 2*seiyL
ciL <- iyL - 2*seiyL
plot(year[c(20:64)], iyL[c(20:64)], type="l", lty=0, 
     main=expression(paste("Unweighted mean (95% CI) of ", y[L])), 
     xlab="Year", ylab="2011 USD/hour", ylim=c(10, 40))
polygon(c(year[c(20:64)], rev(year[c(20:64)])), 
        c(ciL[c(20:64)], rev(ciH[c(20:64)])), 
        col="gray90", border=NA)
lines(year[c(20:64)], iyL[c(20:64)], lty=1)
lines(year[c(20:64)], miyL[c(20:64)], lty=2)
legend("top", bty="n", 
       legend=c(expression(paste(bar(y[L]))), expression(paste(tilde(y[L])))), 
       lty =c(1,2))
mtext(bquote(bar(y[L])^T  == .(round(mean(iyL[c(20:64)]), 1))~ "
             "~se(bar(y[L]))^T ==.(round(mean(seiyL[c(20:64)]), 1))~" 
             "~tilde(y[L])^T  == .(round(mean(miyL[c(20:64)]), 1))~" 
             "~hat(bar(y[L]))  ==.(round(biyL, 3))))

iwN <- tapply(x$wN, x$year, mean, na.rm=TRUE)
summary(iwN)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1962    4884    6488    6301    7534   11237
biwN <- coef(lm(log(iwN[c(20:64)]) ~ year[c(20:64)]))[2]
seiwN <- tapply(x$wN, x$year, std.error, na.rm=TRUE)
miwN <- tapply(x$wN, x$year, median, na.rm=TRUE)
ciH <- iwN + 2*seiwN
ciL <- iwN - 2*seiwN
plot(year[c(20:64)], iwN[c(20:64)], type="l", lty=0, 
     main=expression(paste("Unweighted mean (95% CI) of ", w[N])), 
     xlab="Year", ylab="2011 USD/capita", ylim=c(1500, 15000))
polygon(c(year[c(20:64)], rev(year[c(20:64)])), 
        c(ciL[c(20:64)], rev(ciH[c(20:64)])), 
        col="gray90", border=NA)
lines(year[c(20:64)], iwN[c(20:64)], lty=1)
lines(year[c(20:64)], miwN[c(20:64)], lty=2)
# lines(year[c(20:64)], iwN[c(20:64)] + 2*seiwN[c(20:64)], lty=2)
# lines(year[c(20:64)], iwN[c(20:64)] - 2*seiwN[c(20:64)], lty=2)
legend("top", bty="n", 
       legend=c(expression(paste(bar(w[N]))), expression(paste(tilde(w[N])))), 
       lty =c(1,2))
mtext(bquote(bar(w[N])^T  == .(round(mean(iwN[c(20:64)]), 1))~ "
             "~se(bar(w[N]))^T ==.(round(mean(seiwN[c(20:64)]), 1))~" 
             "~tilde(w[N])^T  == .(round(mean(miwN[c(20:64)]), 1))~" 
             "~hat(bar(w[N]))  ==.(round(biwN, 3))))


  1. Note that \(l\) is the average per-worker hours worked.↩︎

  2. Note that \(\omega\) is the labor share of total income and \(c\) is the consumption share of total expenditure.↩︎