Calibration - continued

library(data.table)
library(fst)

wiot <- read_fst("output/wiot.fst", as.data.table = TRUE)
trade_flows <- read_fst("output/trade_flows.fst", as.data.table = TRUE)

Estimating \(\alpha_j^k\)

This is the share of expenditures in industry \(k\) for importer country \(j\). This must be computed from the destination values, that is, industry \(k\) is found in the columns of WIOD’s input-output matrix. It’s the industry where the importer country is spending inputs from all countries and all other industries. We compute \(\alpha_j^k\) as the proportion of country \(j\) expenditures in this industry to total expenditures.

#' share of expenditures in k industry for in_country j. This must be computed
#' from in_ind!! It's the industry k where the importer country is spending
alpha_jk <- wiot[, by = c("in_country", "in_ind"),
                        .(alpha_jk = sum(value))]
alpha_jk[, by = "in_country",
         alpha_jk := alpha_jk / sum(alpha_jk)]

#' Compatibilize industry name with the rest of data
setnames(alpha_jk, "in_ind", "out_ind")
setkey(alpha_jk, in_country, out_ind)
#' Check alphas sum to one
alpha_jk[, by = in_country, sum(alpha_jk)]
##     in_country V1
##  1:        AUS  1
##  2:        AUT  1
##  3:        BEL  1
##  4:        BGR  1
##  5:        BRA  1
##  6:        CAN  1
##  7:        CHE  1
##  8:        CHN  1
##  9:        CYP  1
## 10:        CZE  1
## 11:        DEU  1
## 12:        DNK  1
## 13:        ESP  1
## 14:        EST  1
## 15:        FIN  1
## 16:        FRA  1
## 17:        GBR  1
## 18:        GRC  1
## 19:        HRV  1
## 20:        HUN  1
## 21:        IDN  1
## 22:        IND  1
## 23:        IRL  1
## 24:        ITA  1
## 25:        JPN  1
## 26:        KOR  1
## 27:        LTU  1
## 28:        LUX  1
## 29:        LVA  1
## 30:        MEX  1
## 31:        MLT  1
## 32:        NLD  1
## 33:        NOR  1
## 34:        POL  1
## 35:        PRT  1
## 36:        ROU  1
## 37:        ROW  1
## 38:        RUS  1
## 39:        SVK  1
## 40:        SVN  1
## 41:        SWE  1
## 42:        TUR  1
## 43:        TWN  1
## 44:        USA  1
##     in_country V1

Estimating \(\pi_{ij}^k\)

The trade shares are defined as \(\pi_{ij}^k\equiv x_{ij}^k / \sum_{i'} x_{i'j}^k\), where the summation index represents all countries, \(i'\in \mathcal{I}\). This can be done straight from the trade_flows data.table previously computed.

pi_ijk <- trade_flows[, by = c("in_country", "out_ind"),
                      .(out_country,
                        pi_ijk = value / sum(value))]
#' Compute $\pi_{jj}^k$ since its useful to estimate trade costs
pi_jjk <- pi_ijk[out_country == in_country]
pi_ijk <- merge(pi_ijk, pi_jjk[, -c("out_country")], 
                 by = c("in_country", "out_ind"), 
                 all.x = TRUE)
setnames(pi_ijk, c("pi_ijk.x", "pi_ijk.y"), c("pi_ijk", "pi_jjk"))
setkeyv(pi_ijk, c("out_country", "in_country", "out_ind"))
rm(pi_jjk)
#' Check $\pi_{ij}^k$ sum to one over in_country and out_ind
pi_ijk[, by = c("in_country", "out_ind"), sum(pi_ijk)]
##       in_country out_ind V1
##    1:        AUS       1  1
##    2:        AUS       2  1
##    3:        AUS       3  1
##    4:        AUS       4  1
##    5:        AUS       5  1
##   ---                      
## 2416:        USA      51  1
## 2417:        USA      52  1
## 2418:        USA      53  1
## 2419:        USA      54  1
## 2420:        USA      55  1

Matching wages to guarantee balanced trade

We assume trade shares in the baseline economy (in the model) are equal to the ones observed in the data (as well as consumption shares over sectors per country), and use this to find the wages that guarantee balanced trade. So, we need to compute what goes into assumption A5 of Costinot, Donaldson, and Komunjer (2012).

\[ \begin{equation} \sum_j\sum_k \pi_{ij}^k\alpha_j^k\gamma_j = \gamma_i \end{equation} \]

After some manipulation of the above equation we arrive at a linear system of the form:

\[ \begin{equation} \Lambda\cdot\gamma = \gamma \end{equation} \]

where \(\gamma=[\gamma_1, \gamma_2, \ldots, \gamma_I]'\) and \(\lambda_{ij}=\sum_k \pi_{ij}^k \alpha_{j}^k\) are the entries of the matrix \(\Lambda\).

We have already estimated, \(\alpha_{j}^k\) and \(\pi_{ij}^k\) from data. With those values we can compute \(\lambda_{ij}\) and assemble the \(\Lambda\) matrix. We are left to solve the system \(\Lambda\cdot\gamma=\gamma\). Beware the constraint that \(\mathbf{1}\cdot\gamma=1\), that is, we need to normalize the resulting eigenvector. This system is the definition of the eigenvector \(\gamma\) associated to the \(\Lambda\)’s matrix egeinvalue equal to one! Therefore, all we need to find \(\gamma\) is this eigenvector and normalize it, such that its entries sum to one. (Note this is not the same as the eigenvector’s norm equal to one!).

lambda_ij <- copy(pi_ijk)
lambda_ij[alpha_jk, on = c("in_country", "out_ind"),
                    lambda_ij := pi_ijk*alpha_jk]
lambda_ij <- lambda_ij[, by = c("out_country", "in_country"),
                       .(lambda_ij = sum(lambda_ij))][
                         order(out_country)]
#' Create the $\Lambda$ matrix
Lambda <- dcast(lambda_ij, out_country~in_country, value.var = "lambda_ij")
Lambda <- as.matrix(Lambda[, -1])
#' Eigen value = 1 and corresponding eigenvector 
# idx <- which.min(abs(eigen(Lambda)$values - 1))
# gamma <- eigen(Lambda)$vectors[, idx] / sum(eigen(Lambda)$vectors[, idx])
# gamma_i <- data.table(out_country = colnames(Lambda), gamma_i = as.numeric(gamma))

#' Other method to compute gamma_i, iterate!
gamma_old <- rep(1/nrow(Lambda), nrow(Lambda))
for (i in seq_len(1000)) {
  gamma_new <- tcrossprod(Lambda, t(gamma_old))
  if (max(abs(gamma_old - gamma_new)) < 1e-5) 
    break
  
  gamma_old <- gamma_new
}

gamma_i <- data.table(out_country = colnames(Lambda), 
                      gamma_i = as.numeric(gamma_new))
#' Check the sum of $\gamma_i=1$, minimum and maximum values
gamma_i[, .(sum = sum(gamma_i), min = min(gamma_i), max = max(gamma_i))]
##    sum          min       max
## 1:   1 0.0001495473 0.2279988

Once we have the solution \(\gamma\), the vector of coutries shares on world income, it is possible to calculate the wages that guaratee assumption A5 as \(w_i=\gamma_i/L_i \cdot\sum_{i'}w_{i'}L_{i'}\), where \(L_i\) represents the labor force in country \(i\) and I normalize USA’s wage to 1.

#' Get Socio-Economic data
employed <- read_fst("output/employed.fst", as.data.table = TRUE)
#' Now merge with gamma
gamma_i <- merge(gamma_i, employed, by = "out_country")
#' Normalize USA wages to 1 and compute the world's wage bill $\sum_i w_iL_i$
wld_wage <- gamma_i[out_country == "USA", employed_i / gamma_i]
gamma_i[, wage_i := gamma_i * wld_wage / employed_i]
setkey(gamma_i, out_country)
#' Merge everything back to trade_flows, the main database
trade_flows <- trade_flows[alpha_jk, on = c("in_country", "out_ind"),
                            nomatch = 0][
                              pi_ijk, on = c("out_country", "in_country", "out_ind"),
                              nomatch = 0][
                                gamma_i, on = c("out_country"), nomatch = 0][
                                  gamma_i, 
                                  on = c("in_country == out_country"),
                                  nomatch = 0]
setnames(trade_flows, 
         c("i.gamma_i", "i.employed_i", "i.wage_i"),
         c("gamma_j", "employed_j", "wage_j"))
#' Write trade_flows to fst file
write_fst(trade_flows, "output/trade_flows.fst")

You should check we have the right values of \(\gamma\). Equation (3) below must hold (i.e., balanced trade).

\[ \sum_{j=1}^{I} \sum_{k=1}^{K} \pi_{i j}^{k} \alpha_{j}^{k} \gamma_{j}=\gamma_{i}, \] ## Python code

# share of expenditures in k industry for in_country j. This must be computed
# from in_ind!! It's the industry k where the importer country is spending
alpha_jk = (wiot
            .groupby(['in_country', 'in_ind'], as_index=False)['value']
            .sum())
alpha_jk['alpha_jk'] = (alpha_jk
                        .groupby('in_country', as_index=False)['value']
                        .transform(lambda x: x/sum(x)))
alpha_jk.drop('value', axis=1, inplace=True)
# Compatibilize industry name with the rest of data
alpha_jk.rename(columns={'in_ind': 'out_ind'}, inplace=True)
# Check alphas sum to one
alpha_jk.groupby('in_country')['alpha_jk'].sum().describe()

# Estimating pi_ijk
pi_ijk = trade_flows[['in_country', 'out_ind', 'out_country']]
pi_ijk['pi_ijk'] = (trade_flows
                    .groupby(['in_country', 'out_ind'], as_index=False)['value']
                    .transform(lambda x: x/sum(x)))
# Compute pi_jjk since its useful to estimate trade costs
pi_jjk = pi_ijk[pi_ijk['out_country'] == pi_ijk['in_country']]
pi_ijk = pi_ijk.merge(pi_jjk.drop('out_country', axis=1),
                      on=['in_country', 'out_ind'],
                      how='left')
pi_ijk.rename(columns={'pi_ijk_x': 'pi_ijk', 'pi_ijk_y': 'pi_jjk'},
              inplace=True)
# ' Check $\pi_{ij}^k$ sum to one over in_country and out_ind
pi_ijk.groupby(['in_country', 'out_ind'])['pi_ijk'].sum().describe()
# Clean environment
del pi_jjk

# Matching wages to guarantee balanced trade
lambda_ij = pi_ijk.merge(alpha_jk, on=jk_cols, how='left')
lambda_ij['lambda_ij'] = lambda_ij['pi_ijk']*lambda_ij['alpha_jk']
lambda_ij = (lambda_ij
             .groupby(ij_cols, as_index=False)['lambda_ij']
             .sum())
# Creates the Lambda matrix
Lambda_mat = lambda_ij.pivot(index='out_country',
                             columns='in_country',
                             values='lambda_ij')
# Vector of gamma guesses
n_countries = Lambda_mat.shape[0]
gamma_old = pd.Series(np.ones(n_countries) / n_countries,
                      index=Lambda_mat.index)
# Iteration to find eigenvalue
for i in range(1000):
    gamma_new = Lambda_mat @ gamma_old
    if max(abs(gamma_new - gamma_old)) < 1e-5:
        break

    gamma_old = gamma_new.copy()

gamma_i = gamma_new.to_frame(name='gamma_i').reset_index()
# Merge Socio-Economic data
gamma_i = gamma_i.merge(employed, on='out_country', how='left')
# Normalize USA wages to 1 and compute the world's wage bill \sum_i w_i L_i
usa_idx = gamma_i['out_country'] == 'USA'
wld_wage = gamma_i.loc[usa_idx, 'employed_i'] / gamma_i.loc[usa_idx, 'gamma_i']
gamma_i['wage_i'] = (gamma_i['gamma_i'] * wld_wage.values
                     / gamma_i['employed_i'])
# Check gamma_i dataframe
gamma_i
# Merge everything back to trade_flows, the main database
trade_flows = (trade_flows
               .merge(alpha_jk, on=jk_cols, how='left')
               .merge(pi_ijk, on=ijk_cols, how='left')
               .merge(gamma_i, on=['out_country'], how='left')
               .merge(gamma_i, left_on='in_country', right_on='out_country',
                      how='left'))
trade_flows.drop('out_country_y', axis=1, inplace=True)
trade_flows.rename(columns={'out_country_x': 'out_country',
                            'gamma_i_x': 'gamma_i',
                            'employed_i_x': 'employed_i',
                            'wage_i_x': 'wage_i',
                            'gamma_i_y': 'gamma_j',
                            'employed_i_y': 'employed_j',
                            'wage_i_y': 'wage_j'},
                   inplace=True)

References

Costinot, Arnaud, Dave Donaldson, and Ivana Komunjer. 2012. “What Goods Do Countries Trade? A Quantitative Exploration of Ricardo’s Ideas.” The Review of Economic Studies 79 (2): 581–608.