详细解读间接推断

间接推断（Indirect Inference）详解

第一章：间接推断的核心思想

1.1 基本直觉

前面讲了三种主流方法：

• NFXP：完全求解动态规划，直接在原模型上做 MLE
• CCP：非参数估计 CCP，然后反演
• GMM/SMM：匹配矩条件

间接推断的根本思想完全不同：

何不用一个更简单、易于估计的"辅助模型"来间接推断原始模型的参数？

具体逻辑：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11

原始模型 θ（我们关心的参数）
    ↓
    生成虚拟数据
    ↓
拟合"辅助模型"
    ↓
辅助模型的参数 β
    ↓
存在映射：β = β(θ)
    ↓
反演：θ = β^{-1}(β)

比喻：

假设我想知道一个黑盒的内部结构（原始模型的参数$\theta$），但直接打开黑盒太困难。

我可以：

1. 向黑盒输入各种信号（模拟数据）
2. 观察黑盒的输出（虚拟数据的统计特征）
3. 用一个简单的模型拟合这个输入-输出关系（辅助模型）
4. 通过辅助模型的参数，反推黑盒的参数

1.2 名字的来源：“间接"vs"直接”

直接推断：

在原模型上直接做推断（NFXP、MLE 都是）

1

原模型 → 参数估计

间接推断：

用辅助模型作"中介"，间接推断原模型

1

原模型 → 辅助模型 → 参数估计

虽然多了一步，但如果辅助模型比原模型简单得多，总成本可能更低。

1.3 间接推断的三个要素

要素 1：原始模型（Structural Model）

定义动态规划问题，参数为$\theta \in \Theta$。

产生的模型含义：CCP $P_\theta(x)$，状态转移$F_\theta$，等等。

要素 2：虚拟数据生成（Simulation）

给定$\theta$，模拟虚拟数据$\{x_s,t, d_s,t\}$（$s=1,\ldots,N_{sim}$; $t=1,\ldots,T$）。

要素 3：辅助模型（Auxiliary Model）

用某个容易估计的模型（通常线性、logit 等）拟合虚拟数据。

辅助模型有参数$\beta$，通常用 OLS、MLE 等标准方法估计。

关键映射：

存在一个函数$\beta(\theta)$，使得用参数$\theta$的原模型生成的虚拟数据，用辅助模型拟合后，得到$\beta \approx \beta(\theta)$。

反演：

从真实数据用辅助模型估计$\hat{\beta}$，然后反演得到$\hat{\theta} = \beta^{-1}(\hat{\beta})$。

第二章：间接推断的数学框架

2.1 正向映射：从$\theta$到$\beta$

定义： 给定参数$\theta$，定义函数

其中：

• $Z^{sim}(\theta)$ = 用参数$\theta$生成的虚拟数据
• $Q_N(\beta; \cdot)$ = 辅助模型的目标函数（如 OLS 的残差平方和）

具体例子：

假设原始模型是动态汽车维修，参数$\theta = (\rho, RC)$。

辅助模型是简单的 logit：

其中$p_t = Pr(d_t=1|x_t)$。

映射过程：

1. 给定$\theta = (\rho, RC) = (0.001, 13)$
2. 用这个参数模拟 10000 辆车，每辆 10 年的记录
3. 对虚拟数据做 logit 回归：$d_t \sim \log(x_t)$
4. 得到回归系数$\beta(\theta) = (\alpha(\theta), \beta(\theta))$

比如，$\beta(0.001, 13) = (0.5, 0.02)$。

这个$\beta(0.001, 13)$就是正向映射的结果。

2.2 反向映射：从观测数据估计$\theta$

步骤 1：对真实数据做辅助模型估计

从真实数据$\{Z_1, \ldots, Z_N\}$，用辅助模型估计：

例如，用 logit 回归，得到$\hat{\beta} = (\hat{\alpha}, \hat{\beta})$。

步骤 2：反演参数

求解方程：

得到$\hat{\theta}_{indirect} = \beta^{-1}(\hat{\beta})$。

问题： 这个反演方程怎么解？

通常没有闭式解，需要数值方法：

1

搜索 θ，使得 β(θ) 尽可能接近 β̂

2.3 搜索框架：从$\beta$反演$\theta$

目标函数：

其中$W$是权重矩阵（通常是单位矩阵或$\hat{\beta}$的协方差矩阵逆）。

算法流程：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20

初始化：θ_0 = 初始猜测

迭代：for iter = 1, 2, ...

  β(θ) 的计算：
    模拟虚拟数据：Z^sim ~ 模型(θ_current)
    用辅助模型拟合虚拟数据：β(θ_current) = Aux_Reg(Z^sim)

  目标函数：
    J(θ) = [β(θ) - β̂]' W [β(θ) - β̂]

  梯度：
    ∇J = 2(∂β/∂θ)' W [β(θ) - β̂]

  参数更新：
    θ_new = θ_current - step_size × ∇J

  收敛检验：
    if ||θ_new - θ_current|| < ε:
      return θ_new

关键点： 需要计算$\partial \beta / \partial \theta$（梯度）

这可以用数值差分：

对每个参数，需要两次完整的模拟和辅助模型回归。

2.4 渐近分布

定理（Smith, 1990）： 在适当条件下，

其中：

更具体地：

其中：

• $\frac{\partial \beta}{\partial \theta}$ = 辅助模型参数对原模型参数的敏感性（Jacobian）
• $V_{\beta}$ = 辅助模型估计量的渐近方差

解释：

间接推断的方差取决于：

1. 映射$\beta(\theta)$的陡峭程度（$\partial \beta / \partial \theta$大 → 方差小）
2. 辅助模型参数的精度（$V_{\beta}$小 → 方差小）

比喻：

如果$\beta(\theta)$对$\theta$的变化很敏感（斜率大），则从$\hat{\beta}$的小变化就能推断出$\theta$的变化，所以最终的方差小。

反之，如果$\beta(\theta)$很平缓，则$\theta$的不确定性很大。

第三章：辅助模型的选择

3.1 辅助模型的性质要求

要求 1：容易估计

应该能用标准方法（OLS、logit、非线性回归等）快速估计。

要求 2：能识别参数

从$\hat{\beta}$应该能唯一反演$\hat{\theta}$。

这要求映射$\beta(\theta)$是一一对应（至少在参数空间的相关部分）。

要求 3：信息量充足

辅助模型应该能"捕捉"原模型参数对数据特征的主要影响。

如果辅助模型对某些参数无敏感性，就无法识别那些参数。

要求 4：计算效率

虽然单次辅助模型估计快，但总的计算量（多次虚拟数据生成+多次辅助模型拟合）仍可观。

3.2 常见的辅助模型选择

选择 1：Logit/Probit 回归（离散选择模型）

最常见于动态离散选择模型的间接推断。

形式：

其中$\Lambda$是 logit 或 probit 函数。

优点：

• 与原模型（CCP）形式相近，容易解释
• logit 有闭式 MLE
• 可轻松加入多项式项、交互项

缺点：

• 仍需多次 logit 估计（虽然比原模型快）

选择 2：线性回归（状态转移）

如果关心状态转移的特性（如里程增长的均值和方差）。

形式：

优点：

• OLS 估计非常快
• 有闭式解

缺点：

• 可能过于简化状态转移
• 对某些非线性效应的捕捉不足

选择 3：配置模型（Distribution Moments）

用样本矩作为辅助参数。

形式：

优点：

• 计算非常快（只需计算统计量）
• 接近 SMM 的思想

缺点：

• 参数$\beta$的个数有限
• 可能遗漏某些重要的模式

选择 4：混合模型

综合多个简单模型的信息。

例子：

用 logit 回归的系数+状态分布的矩作为$\beta$。

$\beta = (\text{logit coefficients}, \text{moments})$

优点： 信息最充分

缺点： 参数维度增加，反演更复杂

3.3 参数匹配的有效性

关键问题： 从$\hat{\beta}$能否确实恢复$\theta_0$？

答案取决于映射的可逆性

定义信息标准：Jacobian 矩阵

如果$\text{rank}(J(\theta_0)) = k$（参数维度），则映射在$\theta_0$附近可逆。

检验映射的有效性：

1
2
3
4

在几个θ值处计算β(θ)，绘图：
  - 如果β(θ)单调或有明确的映射 → 可逆
  - 如果β(θ)平坦或振荡 → 可能不可逆
  - 如果β(θ)有多条通向同一β̂ → 多重根

3.4 辅助模型的设计艺术

设计原则 1：保留关键特征

辅助模型应该对原模型的主要参数敏感。

比如，如果关键参数是折现因子$\beta$，辅助模型应该能捕捉不同$\beta$值下 CCP 的差异（如陡峭程度）。

设计原则 2：计算复杂度 vs 信息量的权衡

• 太简单的辅助模型：快，但可能识别不了所有参数
• 太复杂的辅助模型：信息充分，但计算成本高，失去了简洁性

设计原则 3：经济学的可解释性

辅助模型的参数应该有清晰的经济学含义，便于诊断。

比如，logit 的系数可以解释为边际效应。

第四章：间接推断的实现步骤

4.1 完整的估计流程

第一阶段：前期准备

步骤 1.1：数据整理与描述统计

从真实数据计算需要的统计量：

• 按不同状态分组的选择频率
• 状态分布
• 自相关系数
• 等等

1

计算 β̂：真实数据上的辅助模型估计

步骤 1.2：原模型的参数空间定义

确定参数的合理范围：

• 成本参数应为正
• 折现因子应在(0, 1)
• 等等

1

定义 Θ = [θ_min, θ_max]

步骤 1.3：辅助模型的选择与实现

选定辅助模型（如 logit），准备好估计代码。

1

定义 β = f_aux(data)

第二阶段：计算正向映射$\beta(\theta)$

步骤 2.1：设定网格点

在参数空间内选择一些测试点$\{\theta^{(1)}, \theta^{(2)}, \ldots, \theta^{(R)}\}$。

步骤 2.2：对每个网格点，计算$\beta$值

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

for r = 1 to R:

  # 模拟虚拟数据
  Z_sim = simulate(θ^(r), N_sim=5000, T=20)

  # 用辅助模型拟合
  β^(r) = aux_model(Z_sim)

  # 存储结果
  β_grid(r, :) = β^(r)

结果： 得到映射表$\{\theta^{(r)}, \beta^{(r)}\}_{r=1}^R$

这可以用插值来近似连续的$\beta(\theta)$函数。

第三阶段：反演参数

步骤 3.1：优化搜索

1
2
3
4
5
6

最小化：||β(θ) - β̂||^2

使用：
  - 初值：某个合理的猜测（如分析值）
  - 搜索方向：梯度下降、拟牛顿法等
  - 梯度计算：数值差分（多次调用β(θ)）

步骤 3.2：收敛与验证

检查：

• 目标函数值是否接近最小
• 参数是否在合理范围内
• 梯度是否接近零

1
2
3
4

if ||∇J(θ̂)|| < ε and J(θ̂) < tolerance:
  接受 θ̂
else:
  诊断问题

4.2 计算$\beta(\theta)$的实现细节

问题： 每次计算$\beta(\theta)$都需要：

1. 模拟 N_sim×T 个数据点
2. 拟合辅助模型

这仍然有成本。

优化 1：虚拟数据的规模

选择$N_{sim} \times T$使得$\beta(\theta)$的估计充分精确。

通常$N_{sim} \times T = 5000$到 50000 已足够。

对比 NFXP：每次需要求解完整的贝尔曼方程。

间接推断快很多。

优化 2：缓存与插值

在固定的网格点上计算$\beta(\theta)$一次，然后用插值法（如三次样条）来近似$\beta(\theta)$在任意点的值。

1
2
3
4
5
6

β_grid = compute_grid(θ_grid_points)
β_interp = spline_fit(θ_grid_points, β_grid)

# 在优化中使用插值
for iter:
  β(θ) ≈ β_interp(θ)  # 快速查询，无需再模拟

这样可以加速优化，代价是引入插值误差（通常很小）。

优化 3：梯度的近似

如果用数值梯度$\partial \beta / \partial \theta$，需要$2k$次函数评估（$k$=参数数）。

替代 1：有限差分的稀疏计算

• 只对参数中的"关键"维度计算数值梯度
• 其他维度假设梯度为零或很小

替代 2：自动微分

• 如果能对虚拟数据生成和辅助模型拟合代码进行自动微分
• 可以精确计算梯度，成本更低

4.3 辅助模型估计的标准化

问题： 如果虚拟数据不同，辅助模型的估计量会有随机波动。

这会导致$\beta(\theta)$的估计有噪声，影响反演的精确性。

解决方案 1：增加虚拟样本规模

用更多虚拟数据（$N_{sim}$更大），使估计更精确。

但这增加计算成本。

解决方案 2：固定随机种子

在给定$\theta$时，用固定的随机种子生成虚拟数据。

这样$\beta(\theta)$是$\theta$的确定函数（无随机波动）。

1
2
3
4

function β(θ):
  set_seed(12345)  # 固定种子
  Z_sim = simulate(θ)
  return aux_model(Z_sim)

这样$\beta(\theta)$是光滑的确定性函数，便于优化。

解决方案 3：多次模拟取平均

对每个$\theta$，用不同的随机种子模拟多次，对$\beta$值取平均。

1

β(θ) = (1/M) Σ_m aux_model(simulate(θ, seed=m))

这样也得到确定函数，且噪声更小。

第五章：间接推断的优势与劣势

5.1 相比 NFXP 的优势

5.2 相比 GMM/SMM 的优势与劣势

间接推断 vs GMM（标准 GMM，非模拟）

间接推断 vs SMM（模拟矩方法）

5.3 间接推断的主要优势

优势 1：灵活使用现有工具

辅助模型通常是标准的计量经济模型（logit、线性回归等）。

有成熟的软件包和理论支持。

例子： 可以直接用 R 的glm()函数或 Stata 的logit命令。

优势 2：避免完整的动态规划求解

不需要：

• 离散化状态空间
• 求解贝尔曼方程
• 计算复杂的梯度

只需：

• 模拟虚拟数据（相对直接）
• 拟合辅助模型（标准方法）

优势 3：相对清晰的诊断

如果估计失败或结果奇怪，可以逐步诊断：

• 虚拟数据的特征是否合理？
• 辅助模型的拟合是否好？
• 映射$\beta(\theta)$是否可逆？
• 反演的目标函数是否有问题？

优势 4：易于扩展和修改

改变辅助模型很容易：

• 从 logit 改成 probit
• 加入多项式项或交互项
• 改变被拟合的统计量

无需重新设计原始模型或动态规划框架。

5.4 间接推断的主要劣势

劣势 1：虚拟数据生成仍需原模型

虽然避免了求解贝尔曼方程的完整优化，但仍需：

• 计算 CCP（从某处）
• 模拟状态转移

这仍然需要对原模型有完整的 specification。

劣势 2：映射可能不可逆

如果映射$\beta(\theta)$不是一一对应（flat 或多值），反演会失败。

例子：

假设$\beta(\theta_1) = \beta(\theta_2)$对于$\theta_1 \neq \theta_2$。

则从$\hat{\beta}$无法唯一确定$\theta$。

这种情况无法事前完全避免，需要试验。

劣势 3：模拟噪声的传导

虽然相对轻微，但虚拟数据的有限规模导致$\beta(\theta)$有估计误差。

这会影响最终参数估计的精确性。

与 SMM 类似，需要选择合适的$N_{sim}$。

劣势 4：反演步骤引入复杂性

虽然单次$\beta(\theta)$计算简单，但反演需要多次计算，涉及优化问题。

如果映射$\beta(\theta)$有多个局部最小值（非凸），优化困难。

劣势 5：渐近效率通常次优

间接推断估计量的方差取决于：

1. 原模型参数到$\beta$的映射（$\partial \beta / \partial \theta$）
2. 辅助模型的估计精度

如果映射"不陡"（即$\partial \beta / \partial \theta$小），则方差很大。

相比之下，MLE 在模型正确时是最有效的。

第六章：具体应用案例

6.1 案例 1：汽车维修的间接推断

原始模型： 动态维修决策

参数：$\theta = (\rho, RC)$

• $\rho$：维护成本系数
• $RC$：维修成本

辅助模型： Logit 回归

参数：$\beta = (\alpha, \beta_1, \beta_2)$

映射过程：

对不同的$(\rho, RC)$值，模拟虚拟数据，用 logit 回归：

关键观察：

• 增加$\rho$（成本增加）→ $\alpha$、$\beta_1$减小（维修率下降）
• 增加$RC$（维修成本高）→ $\alpha$、$\beta_1$增加（虽然会延迟维修，但在状态高时更可能维修）

观察数据：

从真实的汽车维修数据，拟合同样的 logit 模型：

得到$\hat{\beta} = (0.42, 0.017, 0.00008)$

反演：

搜索$(\rho, RC)$使得$\beta(\rho, RC) = (0.42, 0.017, 0.00008)$

比如，发现$(\rho, RC) = (0.001, 13)$时，$\beta$最接近观测值。

6.2 案例 2：设备替换的间接推断

原始模型： 企业的设备替换决策

参数：$\theta = (\theta_1, \theta_2)$

• $\theta_1$：运营成本的年龄敏感性
• $\theta_2$：折现因子$\beta$

辅助模型： 企业设备年龄的分布+替换率

映射过程：

对不同$(\theta_1, \theta_2)$值，模拟企业运营过程，计算辅助参数。

优点：

• 三个矩条件可以识别两个参数
• 辅助参数有清晰的经济学含义
• 计算远比 NFXP 快

6.3 案例 3：多维状态的间接推断

原始模型： 两维状态（主设备年龄+从设备年龄）

参数：$\theta = (\theta_1, \ldots, \theta_5)$（5 个参数）

为什么不用 NFXP？

五维离散化会导致$50 \times 50 = 2500$个状态，每个参数的梯度计算成本$O(k N^3) = O(5 \times 2500^3)$不可行。

为什么不用 GMM？

矩条件的设计困难（多维状态的统计特征复杂）。

用间接推断：

• 辅助模型：两个二元 logit（对主、从设备分别做替换决策的 logit）
• 参数：$\beta = (\alpha_1, \beta_{1,1}, \beta_{1,2}, \alpha_2, \beta_{2,1}, \beta_{2,2})$（6 个）

每次迭代：

1. 模拟虚拟企业（包含两维状态转移）
2. 拟合两个 logit（快速）
3. 反演 5 个参数

计算成本：远低于 NFXP。

第七章：间接推断的关键技术细节

7.1 虚拟数据生成的精细设计

问题 1：初始条件

虚拟代理人从哪个初始状态开始？

选择 1：稳态分布

从模型的稳态分布中抽样初始状态。

优点： 虚拟数据的分布与模型的长期均衡一致

缺点： 计算稳态分布可能困难（需要长模拟）

选择 2：数据的初始分布

从观测数据的初始状态分布中抽样。

优点： 保证虚拟数据与真实数据的初始条件相同

缺点： 可能不反映模型的真实稳态

选择 3：burn-in 期

初始运行$T_{burn}$期（不用于拟合辅助模型），然后再收集数据。

优点： 减少初始条件的影响

缺点： 额外的计算（虽然小）

问题 2：虚拟代理人的数量与时期

权衡：

• $N_{sim}$和$T$越大 →$\beta(\theta)$估计越精确 → 参数反演越精确
• 但$N_{sim} \times T$太大 → 计算成本高

经验法则：

• 小模型（1 维状态）：$N_{sim} \times T = 5,000$
• 中等模型（2 维状态）：$N_{sim} \times T = 20,000$
• 大模型（3 维状态）：$N_{sim} \times T = 50,000$

通常$N_{sim} = N$（真实数据大小），$T = 平均观测时间$。

问题 3：虚拟数据的缺失值与不平衡

真实数据可能有：

• 缺失值
• 不平衡的面板（有人在第 5 年退出）

虚拟数据应该镜像这些特征，以保证$\beta(\theta)$能代表真实数据的拟合。

做法：

1
2

虚拟数据生成后，随机删除一些观测，
使得虚拟数据的不平衡特征与真实数据匹配

7.2 映射反演的数值方法

方法 1：直接搜索（Direct Search）

不计算导数，直接搜索$\theta$。

1
2
3

目标函数：J(θ) = ||β(θ) - β̂||²

算法：单纯形法（Nelder-Mead）

优点： 无需导数，稳健

缺点： 收敛慢，需要更多函数评估

方法 2：数值梯度+拟牛顿法

计算有限差分梯度$\partial J / \partial \theta$。

1
2
3
4
5
6
7

∇J(θ) = 2(∂β/∂θ)' [β(θ) - β̂]

其中 ∂β/∂θ 通过数值差分：
  ∂β_i/∂θ_j ≈ [β(θ + εe_j) - β(θ - εe_j)] / (2ε)

每个参数维度需要2次β(θ)计算
总共2k次（k = 参数数）

用 BFGS 等拟牛顿法优化。

优点： 收敛快（二阶方法）

缺点： 需要$2k$次虚拟数据生成和回归

方法 3：插值+梯度混合

在稀疏网格点上计算$\beta(\theta)$，用插值近似。

在优化中用插值的$\beta$和其导数。

1
2
3
4
5
6
7

β_grid = compute_on_grid(θ_grid_points, N_grid)  # 一次性，离线
β_interp(θ) = interpolate(β_grid, θ)              # 快速
∂β_interp/∂θ = interpolate_gradient(β_grid, θ)   # 快速

在线优化：
  minimize J(θ) = ||β_interp(θ) - β̂||²
  用梯度 ∇J = 2(∂β_interp/∂θ)' [β_interp(θ) - β̂]

优点： 平衡了精度和效率

缺点： 引入插值误差（通常很小）

方法 4：平行计算

对不同的$\theta$值并行计算$\beta(\theta)$。

现代多核计算机可以显著加速。

7.3 梯度的数值稳定性

问题： 数值梯度的步长$\epsilon$如何选择？

步长太小：

• 舍入误差导致导数估计不准
• $\beta(\theta+\epsilon) - \beta(\theta)$可能被噪声淹没

步长太大：

• 泰勒展开不准
• 高阶项影响

推荐： $\epsilon = 10^{-6} \times ||\theta||$

即步长与参数的大小相关。

验证梯度：

计算数值梯度和解析梯度（如果有）并比较。

或用不同的$\epsilon$计算梯度，检查一致性。

7.4 收敛诊断

诊断 1：目标函数的减小

$J(\theta)$应该单调减小（可能有小的波动）。

1
2

if J(θ_new) > J(θ_old):
  可能步长太大或梯度计算有误

诊断 2：梯度的大小

在最优点，梯度应接近零。

1
2

if ||∇J(θ̂)|| > 0.01:
  可能未收敛，继续迭代或调整参数

诊断 3：参数的变化

连续迭代的参数变化应减小。

1
2
3

参数变化趋势：
  iter 1-10: 0.1, 0.08, 0.06, 0.04, ...（递减好）
  iter 1-10: 0.1, 0.05, 0.10, 0.03, ...（振荡，可能问题）

诊断 4：多起点搜索

从不同初值开始，检查是否收敛到同一个点。

如果多个初值给出不同的结果，说明有多个局部最小值。

1
2
3
4
5
6
7
8

for init in [init1, init2, init3, ...]:
  θ_hat = optimize(J, init)
  print(θ_hat)

if 所有θ_hat相同:
  全局最优解可能已找到
else:
  存在多个局部最小值，需要进一步调查

第八章：间接推断与其他方法的深入对比

8.1 信息内容的对比

NFXP（完全信息 MLE）：

使用的信息：完整的数据分布（或更准确说，数据的完整对数似然）

GMM（矩方法）：

使用的信息：选定的矩条件

SMM（模拟矩）：

使用的信息：通过虚拟数据匹配的矩

间接推断：

使用的信息：虚拟数据与真实数据在辅助模型上的相似性

信息量排序：

NFXP ≥ GMM（如果矩条件来自 score） > 间接推断 > GMM（任意矩）

直观：NFXP 用尽了数据的所有信息（完整似然）；GMM/SMM 用部分信息；间接推断通过辅助模型"压缩"了信息。

8.2 计算复杂性的对比

假设：

• 参数$k=5$
• 状态$N_x=100$
• 真实数据$N=10,000$
• 虚拟数据（若需）$N_{sim}=10,000$

NFXP：

• 每次评估：贝尔曼求解 $O(I_{bell} \times N_x) \approx 100 \times 100 = 10,000$

  • 梯度计算：$O(k \times N_x^2) \approx 5 \times 10,000 = 50,000$
  • 总：~60,000 ops

• 优化迭代：~500 次

• 总计：~3000 万 ops

间接推断：

• 每次评估$\beta(\theta)$：

  • 模拟：$O(N_{sim} \times T)$（数据生成，很快）
  • 回归：$O(N_{sim} \times T)$（logit 拟合，快）
  • 总：~200,000 ops

• 梯度计算：$2k$次$\beta(\theta)$，即~400,000 ops

• 优化迭代：~100 次

• 总计：~4000 万 ops

或者（如果用网格插值）：

• 离线网格计算：网格点数$\times 200,000$（假设 50 个网格点）
• 在线优化：用插值，只需~100 次迭代，每次非常快
• 总计：~1000 万 ops

结论： 间接推断和 NFXP 在这个例子中计算量相当，但：

• NFXP 的计算更"集中"（求解 DP）
• 间接推断的计算更"分散"（多次模拟和回归）

在维数更高时（如$N_x=500$），NFXP 的成本呈立方增长，间接推断线性增长。

8.3 稳健性的对比

模型误设的影响：

假设真实的 CCP 形式不是 logit 而是 probit，但估计时用 logit 辅助模型。

NFXP 的表现：

所有参数都会有系统性偏差。

无法检测这个误设（没有过度识别检验）。

间接推断的表现：

logit 回归仍能拟合虚拟数据（只是参数估计有偏）。

关键是：虚拟数据的分布特征（如 CCP 的形状）是否与真实数据匹配。

如果真实 CCP 的形状与 logit 的预测差异很大，反演会失败（或得到奇怪的参数值）。

可能有某种"诊断"（如检查是否能在多个初值下收敛）。

GMM/SMM 的表现：

通过 J 检验可以直接检测模型不匹配。

8.4 样本量与渐近性质

小样本（$N$很小）：

• NFXP：MLE 在小样本中相对稳健，但数值优化可能困难
• 间接推断：虚拟数据足够大（固定），所以$\beta(\theta)$估计精确；不确定性主要来自真实数据的$\hat{\beta}$
• GMM/SMM：样本矩不精确，可能影响优化

大样本（$N$很大）：

• NFXP：$\sqrt{N}$-相容，方差$\sim N^{-1}$
• 间接推断：也是$\sqrt{N}$-相容，方差$\sim N^{-1}$（来自$\hat{\beta}$）
• GMM/SMM：相同的渐近性质

第九章：间接推断的高级话题

9.1 映射可逆性的深入分析

定义： 映射$\beta: \Theta \to \mathcal{B}$是可逆的，如果存在连续函数$\theta(\beta)$使得$\theta(\beta(\theta)) = \theta$。

充分条件： Jacobian 矩阵$J(\theta) = \partial \beta / \partial \theta$在所有$\theta \in \Theta$有满秩（列秩=参数数$k$）。

检验方法 1：数值检验

在参数空间的网格点上计算$J(\theta)$的秩。

1
2
3
4
5

for θ in grid:
  J(θ) = numerical_jacobian(β, θ)
  rank(θ) = rank(J(θ))
  if rank(θ) < k:
    警告：θ处映射可能不可逆

检验方法 2：映射的可视化

对简化的 2 参数模型，绘制$\beta(\theta)$的轨迹。

1
2
3

绘图：(β_1(θ), β_2(θ)) 当θ扫过参数空间
  - 如果轨迹覆盖β空间，且无自交 → 可逆
  - 如果轨迹相交或自折 → 不可逆

问题情形 1：Flat 映射

$\beta(\theta)$随$\theta$变化很小。

这导致：

• Jacobian 的列接近平行
• 秩接近 0
• 反演：小的$\hat{\beta}$波动导致大的$\theta$波动（方差爆炸）

解决方案：

• 改变辅助模型（选择对参数更敏感的模型）
• 增加辅助参数的维数

问题情形 2：多值映射

$\beta(\theta_1) = \beta(\theta_2)$对于$\theta_1 \neq \theta_2$。

反演时无法唯一确定$\theta$。

解决方案：

• 改变辅助模型
• 在参数空间的受限区域内工作（如假设某些参数已知）

9.2 辅助模型的多层策略

思路： 不用单一的辅助模型，而是多个分层的辅助模型。

第一层：粗粒度辅助模型

快速计算，识别参数的大致范围。

例子： 简单 logit，线性在状态上。

第二层：精细辅助模型

加入非线性项、交互项等。

例子： Logit，二次多项式在状态上，加入滞后项。

第三层：最细致辅助模型

接近原模型的结构。

例子： 完整的动态 logit（考虑决策的动态相关性）。

实施流程：

1
2
3
4
5
6
7
8

θ_0 = 初值

for level = 1 to 3:
  用level层的辅助模型
  β(θ)_level = 计算正向映射
  θ_level = 反演，用θ_{level-1}作为初值

return θ_3

优点：

• 粗层快速缩小搜索范围
• 细层从好初值开始，容易收敛

缺点：

• 实施复杂

9.3 动态辅助模型

前面的辅助模型是"静态"的（不考虑时间序列结构）。

改进：动态辅助模型

考虑虚拟代理人的历史决策。

例子： 动态 logit，包含滞后因变量

其中$d_{t-1}$是前期决策，$\mathbb{E}[V_{t+1}]$是对未来价值的期望。

优点：

• 更好地捕捉动态模型的特征
• 可能更好地识别参数（如折现因子）

缺点：

• 辅助模型更复杂，计算更慢
• Jacobian 可能更复杂

9.4 部分识别（Partial Identification）

如果映射$\beta(\theta)$不可逆，我们仍可以得到部分识别的结果。

思路： 确定$\theta$的一个集合（而非点估计），使得对应的$\beta(\theta)$都与$\hat{\beta}$在某个容差内。

实施：

扫描参数空间，找所有使目标函数$J(\theta) < \epsilon$的$\theta$。

1
2
3
4
5
6

Θ_partial = []
for θ in parameter_grid:
  if J(θ) < ε:
    Θ_partial.append(θ)

绘制Θ_partial的凸包或其他描述

这给出参数的"可信集"，比点估计更诚实（承认了不确定性）。

第十章：实施例程与代码框架

10.1 概念伪代码框架

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144

class IndirectInferenceEstimator:

  def __init__(self, data, structural_model, aux_model):
    self.data = data
    self.struct_model = structural_model  # 原始模型
    self.aux_model = aux_model            # 辅助模型
    self.N = len(data)

  # ===== 步骤1：计算真实数据上的辅助参数 =====
  def estimate_aux_beta_real(self):
    """对真实数据拟合辅助模型"""
    β_real = self.aux_model.estimate(self.data)
    return β_real

  # ===== 步骤2：计算虚拟数据上的辅助参数（正向映射） =====
  def compute_beta_sim(self, theta):
    """
    给定theta，计算β(theta)
    1. 模拟虚拟数据
    2. 拟合辅助模型
    3. 返回参数
    """
    # 模拟虚拟数据
    data_sim = self.struct_model.simulate(
      theta,
      N_sim=self.N,
      T=self.average_T
    )

    # 拟合辅助模型
    beta_sim = self.aux_model.estimate(data_sim)

    return beta_sim

  # ===== 步骤3：定义反演的目标函数 =====
  def objective(self, theta, W=None):
    """
    目标函数：||β(θ) - β̂||²_W
    """
    beta_sim = self.compute_beta_sim(theta)
    diff = beta_sim - self.beta_real

    if W is None:
      W = np.eye(len(diff))

    return diff @ W @ diff

  # ===== 步骤4：梯度计算（数值） =====
  def gradient(self, theta, eps=1e-6):
    """
    数值梯度：∇J(θ)
    使用中心差分
    """
    k = len(theta)
    grad = np.zeros(k)

    for i in range(k):
      theta_plus = theta.copy()
      theta_plus[i] += eps * (1 + abs(theta[i]))

      theta_minus = theta.copy()
      theta_minus[i] -= eps * (1 + abs(theta[i]))

      J_plus = self.objective(theta_plus)
      J_minus = self.objective(theta_minus)

      grad[i] = (J_plus - J_minus) / (2 * eps * (1 + abs(theta[i])))

    return grad

  # ===== 步骤5：参数估计（优化） =====
  def estimate(self, theta_init, method='BFGS'):
    """
    主要的反演过程
    """
    # 先计算真实数据的辅助参数
    self.beta_real = self.estimate_aux_beta_real()
    print(f"β̂ (from real data) = {self.beta_real}")

    # 优化搜索
    if method == 'BFGS':
      from scipy.optimize import minimize

      result = minimize(
        fun=self.objective,
        x0=theta_init,
        method='BFGS',
        jac=self.gradient,
        options={'disp': True, 'gtol': 1e-5}
      )

      theta_hat = result.x

    elif method == 'nelder-mead':
      from scipy.optimize import minimize

      result = minimize(
        fun=self.objective,
        x0=theta_init,
        method='Nelder-Mead',
        options={'disp': True}
      )

      theta_hat = result.x

    # 最优点的β值
    beta_opt = self.compute_beta_sim(theta_hat)

    return {
      'theta_hat': theta_hat,
      'beta_hat': beta_opt,
      'beta_real': self.beta_real,
      'J_final': self.objective(theta_hat),
      'success': result.success
    }

  # ===== 诊断方法 =====
  def check_mapping(self, theta_grid):
    """
    绘制β(θ)的映射，检查可逆性
    """
    betas = []
    for theta in theta_grid:
      beta = self.compute_beta_sim(theta)
      betas.append(beta)

    betas = np.array(betas)

    # 检查是否单调或有自交
    return betas

  def diagnostic_report(self, result):
    """
    输出诊断报告
    """
    print("="*60)
    print("INDIRECT INFERENCE DIAGNOSTIC REPORT")
    print("="*60)
    print(f"θ̂ = {result['theta_hat']}")
    print(f"β̂ (real data) = {result['beta_real']}")
    print(f"β(θ̂) (sim data) = {result['beta_hat']}")
    print(f"Difference: {result['beta_hat'] - result['beta_real']}")
    print(f"Final J = {result['J_final']:.6f}")
    print(f"Convergence: {result['success']}")

10.2 使用示例（伪代码）

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67

# ===== 定义模型 =====

class StructuralModel:
  def __init__(self, theta):
    self.theta = theta

  def simulate(self, theta, N_sim, T):
    """模拟虚拟数据"""
    x = np.zeros((N_sim, T))
    d = np.zeros((N_sim, T))

    # 初始状态
    x[:, 0] = np.random.normal(100, 50, N_sim)

    for t in range(T-1):
      # CCP（可以来自某个公式或查表）
      P = self.compute_CCP(x[:, t], theta)

      # 决策
      d[:, t] = np.random.binomial(1, P)

      # 状态转移
      x[:, t+1] = x[:, t] + (1 - d[:, t]) * 1 + np.random.normal(0, 5, N_sim)

    return {'x': x, 'd': d}

class AuxiliaryModel:
  def estimate(self, data):
    """拟合logit回归"""
    from sklearn.linear_model import LogisticRegression

    X = data['x'].flatten()
    y = data['d'].flatten()

    # 加入二次项
    X_aug = np.column_stack([np.ones_like(X), X, X**2])

    # Logit回归
    log_reg = LogisticRegression(fit_intercept=False)
    log_reg.fit(X_aug, y)

    return log_reg.coef_[0]

# ===== 估计 =====

# 读入数据
data_real = read_csv('vehicle_data.csv')

# 定义模型
struct = StructuralModel(None)
aux = AuxiliaryModel()

# 创建估计器
estimator = IndirectInferenceEstimator(data_real, struct, aux)

# 初值
theta_init = np.array([0.001, 13.0])

# 估计
result = estimator.estimate(theta_init, method='BFGS')

# 诊断
estimator.diagnostic_report(result)

print(f"\n最终参数估计：")
print(f"  维护成本系数 ρ = {result['theta_hat'][0]:.6f}")
print(f"  维修成本 RC = {result['theta_hat'][1]:.2f}")

第十一章：间接推断的高阶应用与拓展

11.1 “间接推断的间接推断”

思路： 有时候，即使辅助模型的估计也很困难。

可以再套一层间接推断！

例子： 原模型很复杂 → 第一辅助模型（相对简单）→ 第二辅助模型（非常简单）

1
2
3
4
5
6
7

原模型 θ
  ↓ [模拟+第一辅助]
第一辅助参数 β₁
  ↓ [虚拟数据的虚拟数据+第二辅助]
第二辅助参数 β₂
  ↓ [反演]
估计的θ

这很少见，但在极端复杂的模型中可能有用。

11.2 组合方法：CCP+间接推断

思路： 不直接计算原模型的 CCP，而是先用 CCP 半参数法估计经验 CCP，然后在间接推断中使用。

1
2
3
4
5
6
7
8

真实数据
  ↓ [非参数CCP估计]
经验CCP
  ↓
原模型（使用经验CCP）
  ↓ [虚拟数据生成]
  + [间接推断]
→ 参数估计

优点：

• 避免了假设特定的 CCP 形式（如 logit）
• 可能比直接用 logit CCP 更稳健

缺点：

• 多了一个非参数估计步骤
• 计算复杂度增加

11.3 贝叶斯间接推断

思路： 不仅估计参数，还给出参数的后验分布。

框架：

1
2
3
4

后验 ∝ 似然 × 先验

其中似然来自间接推断的"拟合度"：
  p(θ | 数据) ∝ p(β̂ | β(θ)) × p(θ)

实施： 用 MCMC 方法（如 Metropolis-Hastings）抽样后验。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12

def log_posterior(theta, beta_real, prior):
  # 似然：β(θ)与β̂的匹配度
  beta_sim = compute_beta_sim(theta)
  log_lik = -0.5 * (beta_sim - beta_real).T @ (beta_sim - beta_real)

  # 先验
  log_prior = prior.log_pdf(theta)

  return log_lik + log_prior

# MCMC抽样
samples = mcmc_sampler(log_posterior, n_iterations=10000)

优点：

• 自然地处理不确定性
• 可以加入先验信息
• 完整的后验分布，而非点估计

缺点：

• MCMC 计算成本高
• 需要选择先验

11.4 多矩目标（多个辅助模型）

有时候，不是拟合单一的辅助模型，而是多个。

例子：

• 辅助模型 1：CCP 的 logit 拟合
• 辅助模型 2：状态分布的拟合
• 辅助模型 3：自相关结构的拟合

参数：

目标函数：

其中$\lambda_j$是权重，平衡不同目标的重要性。

第十二章：实践指南与决策

12.1 何时使用间接推断

优先使用间接推断的情况：

1. ✓ 参数空间有限（$k \leq 10$）

   • 间接推断的反演问题规模可控

2. ✓ 虚拟数据生成容易

   • 原模型虽然动态复杂，但可以模拟

3. ✓ 存在"好的"辅助模型

   • 容易快速估计，且对参数敏感

4. ✓ 计算资源受限

   • 不想完整求解贝尔曼方程（NFXP 的成本）

5. ✓ 希望灵活调整模型

   • 易于改变辅助模型，快速试验

6. ✓ 多维状态（$\geq 3$维）

   • NFXP 不可行，GMM 矩条件难以设计

不太适合的情况：

1. ✗ 模型 specification 非常不确定

   • GMM 或非参数方法可能更稳健

2. ✗ 需要完整的似然信息

   • NFXP 的最优效率无法达到

3. ✗ 映射$\beta(\theta)$不可逆

   • 反演会失败，需要重新设计

4. ✗ 参数维度很高（$k > 15$）

   • 反演的优化问题高维，困难

12.2 辅助模型的选择决策树

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19

问题特征分析

├─ 关心什么参数？
│  ├─ 成本参数 → CCP型辅助（logit）
│  ├─ 折现因子 → 动态特征（如自相关）
│  └─ 状态转移参数 → 状态分布/转移型辅助
│
├─ 数据的主要特征是什么？
│  ├─ 选择频率的差异大 → logit拟合好
│  ├─ 状态分布有特殊形状 → 矩条件好
│  └─ 时间序列结构复杂 → 动态辅助好
│
├─ 计算能力如何？
│  ├─ 充足 → 可用复杂辅助模型
│  └─ 有限 → 简单模型（线性回归、简单logit）
│
└─ 是否需要标准误？
   ├─ 是 → 选易于导数计算的辅助（如logit）
   └─ 否 → 更自由的选择

12.3 实施检查清单

• 原模型能否模拟？ 代码是否能生成虚拟数据
• 辅助模型能否快速估计？ 与原模型相比，至少快 10 倍
• 映射是否存在？ 在几个$\theta$值处计算$\beta(\theta)$，检查是否有意义
• 映射是否可逆？ 检查 Jacobian 秩，或绘图检查
• 初值是否合理？ 从粗略分析或其他方法得到
• 虚拟样本规模充分？ $N_{sim} \times T$足够使$\beta(\theta)$精确
• 目标函数可优化？ 从多个初值开始，是否收敛到同一点
• 诊断统计量好？ 梯度接近零，参数变化合理
• 结果可信？ 参数值是否与直觉或文献一致
• 敏感性分析？ 对虚拟样本规模、初值、辅助模型的依赖

12.4 问题排查指南

问题 1：优化不收敛

症状：目标函数$J$不减小，或振荡

原因可能：

• 初值太差
• 学习率（步长）太大
• 映射$\beta(\theta)$噪声大（虚拟数据规模太小）
• 映射不可逆（多个$\theta$给相同$\beta$）

解决：

• 尝试不同初值
• 减小步长或用更保守的优化方法（Nelder-Mead）
• 增加虚拟数据规模
• 检查映射的可逆性

问题 2：参数估计值不合理

症状：$\hat{\theta}$在参数空间的边界，或与直觉不符

原因可能：

• 辅助模型信息不足，无法识别该参数
• 映射中该参数的方向不清晰（梯度小）
• 数据与模型有本质不匹配

解决：

• 改变辅助模型，增加对该参数的敏感性
• 检查梯度$\partial \beta / \partial \theta$，看该参数的影响
• 审视模型假设和数据特征的匹配度

问题 3：标准误很大

症状：估计的参数有很大不确定性

原因：

• 映射"不陡"（$\partial \beta / \partial \theta$小）
• 辅助模型对参数不敏感
• 真实数据噪声大

解决：

• 选择更敏感的辅助模型
• 如果可能，增加数据量
• 接受这个不确定性（反映模型的识别性质）

第十三章：与其他方法的最终对比

13.1 五种方法的全面对比表

13.2 方法选择的决策流程

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

START

计算资源充足且时间充裕？
├─ YES → 是否完全相信模型的specification？
│        ├─ YES → NFXP（最优效率）
│        └─ NO → CCP（更稳健）
│
└─ NO → 状态维数多少？
         ├─ 1-2维 → GMM或SMM
         ├─ 3维 → SMM或间接推断
         └─ 4+维 → 间接推断或SMM（需特殊设计）

END

13.3 各方法的"黄金准则"

NFXP：

• 黄金准则：如果能用，就用它（最优渐近效率）
• 但前提是计算资源充足、模型 specification 清晰

CCP：

• 黄金准则：计算与稳健性的最佳折衷
• 适合绝大多数应用

GMM/SMM：

• 黄金准则：高维、复杂模型的首选
• 灵活、稳健、相对快速

间接推断：

• 黄金准则：当存在"好的"辅助模型时的绝佳选择
• 计算效率与实施灵活性的结合

第十四章：总结与实践建议

14.1 间接推断的精髓总结

核心思想：

用简单模型间接推断复杂模型的参数。

通过虚拟数据建立两个模型之间的映射，再通过反演参数。

三大支柱：

1. 虚拟数据生成：用原模型模拟数据（相对容易，无需优化）
2. 辅助模型拟合：用简单模型拟合虚拟和真实数据（快速标准）
3. 参数反演：反演映射，得到参数估计（需优化，但低维）

14.2 实践建议

第一阶段：问题理解与设计

1. 充分理解原始模型的动态
2. 识别关键参数和要估计的内容
3. 探索哪些简单模型可能作为辅助
4. 对几个参数值手工计算$\beta(\theta)$，检查映射的合理性

第二阶段：原型实现

1. 实现虚拟数据生成代码
2. 实现辅助模型拟合代码
3. 手工计算一些$\beta(\theta)$值，建立"映射表"
4. 可视化映射，检查可逆性

第三阶段：正式估计

1. 对真实数据估计$\hat{\beta}$
2. 设置优化问题（目标函数、初值、方法）
3. 运行优化，多个初值验证
4. 进行诊断和敏感性分析

第四阶段：结果评估与报告

1. 检查参数的合理性（符号、大小、与文献的对比）
2. 报告标准误或置信区间
3. 进行 sensitivity analysis（对虚拟样本规模、初值、辅助模型的依赖）
4. 讨论局限和假设

14.3 常见陷阱与规避方法

14.4 与文献的对接

关键论文：

• Smith, Richard A. (1990). “Estimating nonlinear time-series models using simulated vector autoregressions.” Journal of Applied Econometrics, 5(1), 63-84.

  • 原始论文，奠定理论基础

• Gourieroux, C., Monfort, A., & Renault, E. (1993). “Indirect inference.” Journal of Applied Econometrics, 8(S1), S85-S118.

  • 重要推广和应用

• Hyslop, D. R., & Imbens, G. W. (2001). “Bias from classical and bootstrap standard errors of within-group estimators.” Econometrica, 69(2), 481-487.

  • 渐近性质

• 在产业组织和劳动经济学的应用：Bajari, Benkard, Levin 等人的系列论文

最后：快速参考卡

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30

┌────────────────────────────────────┐
│ 间接推断的快速检查清单              │
├────────────────────────────────────┤
│                                    │
│ 【前期准备】                       │
│ ☐ 原模型能模拟？                   │
│ ☐ 选好辅助模型？                   │
│ ☐ 检查β(θ)的映射？                 │
│                                    │
│ 【估计过程】                       │
│ ☐ 计算β̂（真实数据）？              │
│ ☐ 设置优化问题？                   │
│ ☐ 多个初值验证？                   │
│                                    │
│ 【诊断报告】                       │
│ ☐ 梯度接近零？                     │
│ ☐ 参数在合理范围？                 │
│ ☐ 标准误可接受？                   │
│ ☐ 敏感性分析完成？                 │
│                                    │
└────────────────────────────────────┘

【五行总结】

间接推断：
  用简单模型拟合虚拟数据，
  建立原模型参数与辅助参数的映射，
  反演得到真实参数估计。

  灵活、实用、通常有效。

附录：数学符号速查