贸易成本结构估计：引力、Head-Ries 指数与 EK 反演
本页自包含：结构引力方程推导 → 多边阻力项（MR）→ Head-Ries 指数从贸易份额反推 → EK 反演法 → 边境效应 → PPML(ppmlhdfe) 估计 Stata 代码 → EK 反演与贸易成本矩阵可视化 Python 代码 → 数据来源 → 关税/非关税分解 + 论文案例。
01 经济环境与假设
贸易成本是两地之间、把一件商品从卖方运到买方手中所需的全部代价——不仅是运费，还包括关税、非关税壁垒、语言文化差异、信息摩擦、合同执行成本，甚至"看不见"的家乡偏好。结构引力给出了贸易流与贸易成本的理论映射，本页要解决的问题是反过来：给定观测到的贸易流，如何把贸易成本 $\tau_{ni}$ 测出来？
核心假设清单（每条附设定理由 / 参数含义 / 经济直觉）：
- 假设 A1：Armington / EK 偏好下的引力框架。两国 $n$（进口）、$i$（出口），产品有来源国偏好（Armington）或生产率有 Fréchet 分布（EK），最终都导出同一形式的引力方程。理由：引力方程不是纯经验关系，而是有微观结构的——这是它能做"结构估计"而非"简化估计"的根本。参数含义：$\sigma>1$ 是 Armington 替代弹性，$\theta$ 是 EK 贸易弹性。直觉：产品越是同质（$\sigma$ 小 / $\theta$ 大），贸易流对成本越敏感。
- 假设 A2：iceberg 贸易成本 $\tau_{ni}\ge 1$。1 单位商品从 $i$ 出发只有 $1/\tau_{ni}$ 到达 $n$，损失部分"融解"在路上。理由：省去运输部门建模，且与引力/EK 闭式兼容。参数含义：$\tau_{ni}=1$ 代表无成本（仅理论极限），$\tau_{nn}=1$（国内无跨境成本）。直觉：关税 $t_{ni}$ 只是 $\tau_{ni}$ 的一小部分。
- 假设 A3：贸易成本可分解为可观测 + 不可观测。$\tau_{ni}=\text{距离}^{ρ}\cdot(1+t_{ni})\cdot\exp(\text{其他壁垒})$。理由：距离、关税可观测，是政策分析的抓手；剩下的非关税壁垒、边境效应是残差。参数含义：$\rho$ 是距离弹性，$t_{ni}$ 是从价关税。直觉：中美贸易成本 60% 以上不是关税，而是关税以外的一切。
- 假设 A4：结构估计 vs 简化估计。简化估计只跑 $\ln X_{ni}=\beta_0+\beta_1\ln\text{dist}+\epsilon$，得到的是条件相关；结构估计先推出理论引力方程，再估计出"与福利相关"的弹性（如 $\sigma-1$ 或 $\theta$）。理由：只有结构弹性才能做 ACR 福利反事实。直觉：$\ln$ 距离系数 −1 到底意味着什么福利？只有结构模型告诉你。
符号与参数总表
| 符号 | 含义 | 典型取值/来源 |
|---|---|---|
| $X_{ni}$ | $n$ 从 $i$ 的进口额 | Comtrade/BACI |
| $Y_i$ | $i$ 国总产出 | PWT/WB WDI |
| $Y$ | 世界总产出 | $\sum_i Y_i$ |
| $\tau_{ni}$ | 双边贸易成本（iceberg） | 本页估计对象 |
| $\Pi_i,P_n$ | outward / inward 多边阻力 | 结构引力核心 |
| $\sigma$ | Armington 替代弹性 | $4\sim 9$ |
| $\theta$ | EK 贸易弹性 | $4\sim 5$ |
| $t_{ni}$ | 从价关税（可观测） | WITS/TRAINS |
| $\pi_{ni}$ | $n$ 从 $i$ 进口份额 | $X_{ni}/X_n$ |
02 结构引力方程推导
由 Armington 偏好（CES 复合商品）+ 完全竞争 + 贸易成本 $p_{ni}=p_i\tau_{ni}$，$n$ 地对 $i$ 地产品的支出份额为：
设定理由：这是 Anderson-van Wincoop (2003) 的核心贡献——把早期引力方程里被忽略的"多边阻力"显式写出来。早期工作（Tinbergen 1962）只写 $X_{ni}=G Y_iY_n/\tau_{ni}$，漏掉 $\Pi_i,P_n$，导致遗漏变量偏误。
经济直觉：加拿大把货卖给美国明尼苏达，成本取决于美加距离，但也取决于"美国从日本、德国买同样的货有多便宜"——这就是 $P_n$。同样的明尼苏达—曼尼托巴距离，相对于"明尼苏达—加州"的距离，才决定了边境效应有多大。这正是 McCallum (1995) 美加边境之谜能被解开的钥匙。
$\Pi_i$ 与 $P_n$ 由下面的定义式自洽决定（结构引力的"一般均衡"闭环）：
03 对数线性化与 MR 项处理
把 (TC2) 两边取对数：
关键问题：$\Pi_i$ 只随出口国变、$P_n$ 只随进口国变。因此在面板回归里，它们正好被出口国固定效应 $\mu_i$ 和进口国固定效应 $\nu_n$吸收：
这就是为什么今天的引力回归必须加出口国/进口国（或时间-出口国、时间-进口国）固定效应——不是为了"稳健性"，而是为了恰好控制住多边阻力，否则距离系数有偏。这是 Anderson-van Wincoop 2003 之后引力文献的铁律。但注意：对数化把 $X_{ni}=0$ 的观测全部丢掉了——小国之间的零贸易流恰恰是贸易成本最高的样本，丢掉它们会系统性低估成本。这正是下一节 PPML 的理由。
当误差项与解释变量相关（异方差）时，$\ln X$ 的 OLS 估计不是 $E[X]$ 的矩，会有偏。Silva & Tenreyro (2006) 证明：在结构引力下，应该直接对 $X_{ni}$ 用泊松伪极大似然（PPML），把 $X_{ni}$ 的水平值当因变量，零贸易流保留。这是 03 引力与 PPML 页的核心，本页 Stata 代码直接用 ppmlhdfe。
04 Head-Ries 指数：从贸易份额反推
PPML 只能估计"成本的弹性"（距离系数、关税系数），不能给出双边贸易成本的绝对水平。Head & Ries (2001) 与 Novy (2013) 提出一个天才的恒等式：利用"双向贸易流 + 国内销售"，可以在不估计 MR、不假设模型参数的情况下，直接把 $\tau_{ni}$ 反推出来。
由结构引力，$n$ 从 $i$ 进口 $X_{ni}$，$i$ 从 $n$ 进口 $X_{in}$，国内销售 $X_{nn}$（$n$ 卖给自己）、$X_{ii}$（$i$ 卖给自己）。四式相除，$\Pi_i,P_n$ 和 $Y_i,Y_n$ 全部消掉：
Novy (2013) 的福利形式：在已知替代弹性 $\sigma$ 时，可把 HR 指数换算成"从价等值贸易成本"（ad-valorem tariff equivalent）：
经济直觉：如果 $n$ 卖给自己的份额 $X_{nn}/Y_n$ 很高，说明跨境很难；如果 $n$ 和 $i$ 互相卖给对方的份额很高，说明跨境容易。HR 指数把"相对国内贸易的跨境贸易"这个比例直接转成成本。它是衡量"国内市场整合度"的利器——中国省间贸易成本就是用省间铁路货运量 + 省总产出构造 $X_{nn},X_{ni}$ 后算 HR。
(TC4) 里 $X_{ni}$ 与 $X_{in}$ 相乘，天然假设" $n$ 到 $i$ 的成本 ≈ $i$ 到 $n$ 的成本"。这在对称壁垒（距离、共同语言）下成立；但关税、非对称边境管制会让 $\tau_{ni}\ne\tau_{in}$，HR 只能给出"均值"。非对称场景要用 EK 反演（下一节）或结构引力估计。
05 EK 反演法
在 EK 模型里，贸易份额为 $\pi_{ni}=\dfrac{T_i(w_i\tau_{ni})^{-\theta}}{\sum_k T_k(w_k\tau_{nk})^{-\theta}}$。取 $n$ 对自己的份额 $\pi_{nn}=\dfrac{T_n(w_n\tau_{nn})^{-\theta}}{\sum_k(\cdot)}$（$\tau_{nn}=1$），两式相除，分母相同消掉：
经济直觉：如果 $n$ 从 $i$ 进口份额远低于 $n$ 自给份额，同时 $i$ 的生产率 $T_i$ 并不比 $n$ 差、工资 $w_i$ 也不更高，那唯一的解释就是" $n$ 到 $i$ 的贸易成本很高"。EK 反演把这一差额归因到 $\tau_{ni}$。它的代价是需要先估出生产率 $T_i$ 与工资 $w_i$——而这两者本身又是模型内生的，需要用 inversion 从均衡反推。
06 边境效应（McCallum 之谜）
McCallum (1995, AER) 估计加拿大各省与美国各州的贸易流，发现：即使控制了 GDP、距离，加拿大省际贸易流量是加拿大省—美国州贸易流量的22 倍。这就是著名的"边境效应"（home bias / border effect）——一道无形的边境，把国际贸易压缩到远低于预测水平。
在引力方程里加一个"共同国家"虚拟变量 $\text{samecountry}_{ni}$：
经济直觉：小国（加拿大）的边境效应被放大，因为它替代选项少——加拿大省际贸易受边境影响的相对变化大；大国（美国）省—州贸易受边境影响小。这就是为什么 McCallum 用加拿大数据会得到 22 倍，而用美国数据会小得多。理解 MR 之后，边境效应不再是谜，而是结构模型里的一个矩。
07 估计方法与代码
7.1 Stata：ppmlhdfe 估计结构引力
当代引力估计的金标准是 ppmlhdfe（Correia, Guimarães, Zylkin 2020）：PPML 估计 + 高维固定效应，自动吸收出口国/进口国/时间固定效应，且保留零贸易流。
* ============================================================
* 结构引力 PPML 估计 + 贸易成本弹性 (Anderson-van Wincoop 2003)
* 数据: 双边贸易流面板, 每对(i,n,t)一行
* 需要: ssc install ppmlhdfe, replace
* ============================================================
clear all
set more off
* --- 0. 载入双边贸易数据 ---
* flow_ni : n从i进口额(水平值, 不是对数!)
use "gravity_panel.dta", clear
* 构造时间×出口国 / 时间×进口国固定效应(吸收MR)
egen exp_yr = group(exp_iso year) // 出口国-时间 FE
egen imp_yr = group(imp_iso year) // 进口国-时间 FE
* --- 1. PPML 基准: 出口国-时间 + 进口国-时间 FE ---
* 这是 AER 2003 之后的标准: 严格控制多边阻力
ppmlhdfe flow_ni ln_dist ln_tariff contig comlang rta, ///
absorb(exp_yr imp_yr) cluster(pair_id)
* 系数解释:
* ln_dist 系数 = -(sigma-1)*rho (贸易成本的距离弹性)
* ln_tariff 系数 = -(sigma-1) (关税弹性, 直接得 sigma-1)
* rta 系数 >0: 区域贸易协定降低成本
* --- 2. 边境效应: 加 "同国内部" 虚拟变量 ---
gen same_country = (exp_iso == imp_iso) // 国内销售=1
ppmlhdfe flow_ni ln_dist same_country contig comlang, ///
absorb(exp_yr imp_yr) cluster(pair_id)
* exp(_b[same_country]) 即国内贸易相对跨境的倍数(边境效应)
* --- 3. 用估出的 sigma-1 换算贸易成本弹性 ---
* 若 ln_tariff 系数为 -1.5, 则 sigma-1 = 1.5, sigma = 2.5
* 贸易成本 tau 每升 10%, 贸易流下降 1.5*10% = 15%
* --- 4. 稳健性: 对比 OLS 对数化(故意保留零流问题) ---
gen ln_flow = ln(flow_ni)
reg ln_flow ln_dist contig comlang i.exp_iso i.imp_iso i.year
* 注意: 这条 OLS 会自动丢掉 flow_ni=0 的行, 系数有偏
7.2 Python：Head-Ries 指数 + EK 反演 + 成本矩阵
# ============================================================
# 贸易成本反推: Head-Ries 指数 + EK 反演法
# 输入: 双边贸易矩阵 X (n从i进口, 对角=国内销售)
# 生产率 T_i, 工资 w_i, 贸易弹性 theta
# 输出: HR成本矩阵 tau_HR, EK反演矩阵 tau_EK
# ============================================================
import numpy as np
def head_ries(X):
"""Head-Ries 指数: tau_HR[n,i] = sqrt(X[n,i]*X[i,n] / (X[n,n]*X[i,i]))
X[n,i] = n 从 i 进口 (含对角元 X[n,n] = n 自给)"""
X = np.asarray(X, dtype=float)
n = X.shape[0]
diag = np.diag(X) # X[n,n]
num = X * X.T # X[n,i]*X[i,n]
den = diag[:, None] * diag[None, :] # X[n,n]*X[i,i]
hr = np.sqrt(num / den)
np.fill_diagonal(hr, 1.0) # 国内成本=1
return hr
def novy_tariff_equiv(X, sigma=5.0):
"""Novy (2013): 从价等值贸易成本 = (X_nn X_ii / (X_ni X_in))^(1/(2(sigma-1))) - 1"""
hr = head_ries(X) # = sqrt(X_ni X_in / (X_nn X_ii))
# 注意: hr = sqrt(bilater/domestic), Novy 用 domestic/bilater = 1/hr^2
novy = (1.0 / hr**2) ** (1.0 / (2 * (sigma - 1))) - 1.0
np.fill_diagonal(novy, 0.0)
return novy
def ek_inversion(X, T, w, theta):
"""EK 反演: tau_ni = ((pi_ni/pi_nn)/(T_i/T_n))^(-1/theta) * w_i/w_n
X[n,i] = n 从 i 进口; T[i]=i地生产率; w[i]=i地工资"""
X = np.asarray(X, float)
row_sum = X.sum(axis=1, keepdims=True) # n地总进口(含国内)
pi = X / row_sum # pi[n,i] = n从i进口份额
pi_nn = np.diag(pi)[:, None] # n 自给份额 (n,n)
ratio_pi = pi / pi_nn # pi_ni / pi_nn
T = np.asarray(T, float); w = np.asarray(w, float)
ratio_T = T[None, :] / T[:, None] # T_i / T_n (n行,i列)
ratio_w = w[None, :] / w[:, None] # w_i / w_n
tau = (ratio_pi / ratio_T) ** (-1.0 / theta) * ratio_w
np.fill_diagonal(tau, 1.0)
return tau
# ------------------------------------------------------------
# 用模拟贸易矩阵演示 (实际替换为 BACI/海关数据)
# ------------------------------------------------------------
if __name__ == "__main__":
rng = np.random.default_rng(0)
N = 5
# 模拟: 国内销售大, 跨境贸易小; 构造一个 (n从i进口) 矩阵
X = rng.uniform(0.5, 2.0, (N, N))
np.fill_diagonal(X, rng.uniform(8.0, 12.0, N)) # 国内销售远大于跨境
X = X + X.T # 近似对称贸易
hr = head_ries(X)
print("Head-Ries 指数矩阵 (越接近0跨境越难):")
print(np.round(hr, 3))
print("HR 行平均(排除对角):",
np.round(np.nanmean(np.where(np.eye(N) == 0, hr, np.nan)), 3))
novy = novy_tariff_equiv(X, sigma=5.0)
print("\nNovy 从价等值成本 (sigma=5):")
print(np.round(novy, 3))
# EK 反演
T_prod = rng.uniform(0.8, 1.3, N)
w = rng.uniform(0.9, 1.1, N)
tau_ek = ek_inversion(X, T_prod, w, theta=4.0)
print("\nEK 反演双边成本 (theta=4):")
print(np.round(tau_ek, 3))
print("EK 反演对角=1:", np.allclose(np.diag(tau_ek), 1.0))
- 对角元是命门：HR 与 EK 反演都需要 $X_{nn}$（国内销售）。BACI 只有跨境贸易，必须用"总产出 − 出口"近似 $X_{nn}$；
- 对称性：HR 自动用 $X_{ni}X_{in}$，要求数据双向都有；单向数据只能用 EK 反演；
- Stata vs Python 分工：回归弹性用
ppmlhdfe（Stata 生态成熟），反推水平值用 Python（矩阵运算快）。
08 数据来源
| 数据 | 来源 | 用途 |
|---|---|---|
| 双边贸易流 $X_{ni}$ | UN Comtrade；Cepii BACI（已修正 CIF/FOB 不一致）；WITS | 引力回归、HR 指数 |
| 引力变量（距离/语言/边界/RTA） | Cepii GeoDist（Mayer-Zignago）；CEPII Gravity 数据库 | $\ln\tau_{ni}$ 的可观测部分 |
| 关税 $t_{ni}$ | WTO TAO；UNCTAD TRAINS；WITS（整合后） | 关税 vs 非关税分解 |
| GDP / 总产出 | Penn World Table；IMF WEO；World Bank WDI | $Y_i,Y_n,Y$ |
| 国内销售 $X_{nn}$ | 总产出 − 总出口（工业普查/国民账户）；WIOD 对角块 | HR 指数分母 |
| 中国微观 | 海关数据库（HS8 位企业-产品-目的地，2000-2016）；中国工业企业库 | 企业层面引力、加工贸易 |
| 省间贸易 $X_{ni}^{\text{省}}$ | 省间铁路货运量；引力模型反推（如 Tombe-Zhu 2019） | 中国国内市场整合 |
实操提醒：(i) Comtrade 的进口用 CIF、出口用 FOB，$X_{ni}\ne X_{in}$，必须用 BACI 做镜像调和；(ii) 关税用 HS 编码，贸易流用 SITC/ISIC，必须做 HS↔ISIC↔SITC 多版本匹配（WITS 自带 concordance）；(iii) 中国省间贸易没有官方统计，常用"铁路货运量 + 引力反推"近似，会有测量误差。
09 模拟与应用：成本分解与趋势
9.1 贸易成本矩阵的描述性分析
算出 $\tau_{ni}$ 后，按 (i) 距离分桶、(ii) 是否同 RTA、(iii) 是否同语言、(iv) 是否同国家（边境效应）分组，看平均成本如何变化。典型发现：距离每翻一倍，贸易成本上升约 20–30%；同 RTA 使成本下降 5–15%；同语言使成本下降约 5%。
9.2 关税 vs 非关税壁垒的贡献分解
从 PPML 回归可得：可观测关税只解释双边总贸易成本的一小部分（通常 < 20%），剩下 80% 是非关税壁垒（技术性贸易壁垒 TBT、检验检疫、许可证、文化距离、信息摩擦）。分解方法：把 $\ln\tau_{ni}=\rho\ln\text{dist}+\ln(1+t_{ni})+\boldsymbol\gamma'\mathbf{Z}+\tilde\tau_{ni}$，关税项贡献 $\ln(1+t_{ni})$，残差 $\tilde\tau_{ni}$ 即非关税。这正是"中美贸易战加了 25% 关税，但实际成本上升远不止"的量化来源。
9.3 贸易成本随时间下降趋势
用 HR 指数逐年算 $\tau_{ni,t}$，1950–2020 年全球平均贸易成本下降约 30–50%（Novy 2013 估计），但近 10 年下降停滞甚至回升——这是"慢全球化"（slowbalization）和贸易战的量化证据。中国内部省间贸易成本自 1995 年以来持续下降，但地区间差距仍大，统一大市场建设的空间就在这里。
10 论文案例、常见错误与难度
10.1 英文经典论文
10.2 中文顶刊与中国应用
10.3 常见错误
ln(flow) 直接把 flow=0 的样本丢掉。小国之间零贸易恰恰是高成本观测，丢掉会严重低估贸易成本弹性。必须用 PPML 保留水平值。
不加出口国/进口国固定效应（或时间×国固定效应），距离与边境系数都有偏——这正是 McCallum 22 倍之谜的来源。没加 MR 的引力回归在 2024 年不能投稿。
HR 用 $X_{ni}X_{in}$，本质假设对称。若研究非对称关税或单边制裁，HR 会把双边成本"平均化"，应改用 EK 反演。
关税是 HS6/HS8 位，贸易流常是 SITC/ISIC。直接按 HS 年份合并会因编码版本变动（HS1996→2002→2007→2012→2017→2022）错配。必须用 WITS/BACI 的官方 concordance 并按版本映射。
HR 与 EK 反演都需要对角元 $X_{nn}$。BACI 只含跨境贸易，直接算 HR 会除零。必须用"总产出 − 出口"或 WIOD 对角块补国内销售。
10.4 难度与前置
难度：★★★☆☆。需要掌握：结构引力与 PPML（03 页）、EK 模型（02 页）、ACR 福利（04 页）。本页是 03 页的深化——03 教你跑 PPML，本页教你把 PPML 系数 + HR + EK 反演拼成一套完整的"贸易成本工具箱"。与 GVC 增加值互补：贸易成本决定"流多少"，增加值分解决定"流的里面有多少是真本国造的"。