简介MATLAB中异构代理新凯恩斯HANK模型的完整复制包面向经济学研究者、研究生及政策分析人员聚焦异质性主体宏观模型的构建与求解。资源共321个文件压缩包约22.32MB涵盖56个m主程序、69个f90/Fortran模块、75个txt说明文档、18个f源文件及13个makefile构建脚本另附mat数据、pdf/eps图表和png示意图可覆盖从模型设定、动态方程求解、参数校准到模拟分析的完整流程。工程中还包含CUDA源文件如estimate.cu、piestimator.cu与ss.dat、gmat.dat、sol.dat等模拟结果便于研究者观察GPU加速估计和数值输出的对接方式。资源包目录结构清晰已有46人浏览/学习通过源码、文档与输出文件的配合读者可直接复现基准结果并在此基础上扩展情景模拟与政策评估实验。 异构代理新凯恩斯HANK模型这几年在宏观经济学里的热度肉眼可见做货币和财政政策传导研究的人几乎绕不开它。它把新凯恩斯框架中的价格粘性和家庭层面的收入风险、资产分布异质性拼到了一起很多在传统代表性代理人RANK模型里说不清的问题——比如一次性转移支付为什么能产生远超RANK模型的消费乘数——在HANK框架下直接有了量化解释。我最近把McKay、Nakamura和Steinsson2016那篇经典HANK模型用MATLAB从头到尾复制了一遍从状态空间离散化、政策函数迭代到一般均衡和脉冲响应全部走通过程中踩了不少坑。这篇就完整记录复刻思路、数值实现和调试经验给准备入坑或正在调模型的同学一份能落地参考。1. 模型架构与复刻思路1.1 HANK模型到底在解什么问题HANK模型全称Heterogeneous Agent New Keynesian直译是异构代理新凯恩斯模型。跟传统RANK模型最大的区别在于RANK假设所有家庭同质只存在一个“代表性”家庭HANK则让家庭在资产持有量和劳动收入上存在持续的异质性分布每个家庭面临不可保的劳动收入风险因此会进行预防性储蓄并且受到借贷约束的限制。具体到我复刻的这个基准版本模型包含几个核心模块连续统家庭CRRA效用函数单一无风险资产借贷约束a ≥ a_min劳动收入服从AR(1)对数过程最终品生产者将连续中间品加总为最终品中间品生产者线性生产技术Calvo粘性定价政府发行债券征收总额税稳定债务货币政策泰勒规则名义利率对通胀和产出作出反应这个设定在HANK文献中属于“最小可行版本”但正因为简洁稳态可以精确求解过渡动态有稳定的数值算法。用这个版本做复刻能帮你把模型机制和数值方法同时掌握为后续扩展做好准备。如果一上来就碰Kaplan-Moll-Violante那种带资本和投资调整成本的版本数值难度会陡增出错后很难定位问题根源。1.2 为什么用MATLAB做复现经济学论文的复现工具主要是MATLAB、Python和Julia三类。选择MATLAB有几个实际原因第一经济学系和宏观研究组里MATLAB的存量代码最多。很多经典HANK论文的原始复制包就是MATLAB写的直接对照着改比自己从零搭Python工程省力得多。第二MATLAB的优化工具箱和稀疏矩阵支持相当成熟。后面会讲到分布转移矩阵的求解涉及大规模稀疏特征值问题MATLAB的eigs函数在内存效率和稳定性上都经得起考验。第三对纯学术验证场景来说MATLAB的开发速度优势明显。不用管包管理、类型声明一个.m文件跑到底很适合复现论文时的快速试错。MATLAB的短板当然也明显——循环慢、并行扩展性不如Julia、开源生态弱。这里分享一个判断准则如果你只是复现模型、验证机制和跑几个脉冲响应MATLAB完全够用如果你要做高维参数扫描或把模型推到非平稳大状态空间那再考虑切换到Julia。1.3 复刻目标与验证基准复刻不是把论文代码原样跑一遍就完了关键是要定义“成功复刻”的检验标准。我给自己定的目标是三件事稳态匹配资产分布的均值、Gini系数、边际消费倾向MPC分布与论文数字大致吻合脉冲响应形态给一个25个基点的货币政策冲击总消费和总产出的响应路径与论文图在方向和量级上一致敏感性检验改一个参数比如Calvo概率或贴现因子模型表现符合经济学直觉有了明确的验收标准调试时就不会陷入“看起来差不多”的模糊状态。建议你复刻任何模型之前都先列出这样一个checklist后面每一步都拿它来对照。2. 家庭异质性模块的数值实现家庭模块是HANK模型的重点也是代码中最容易出bug的地方。整个模块可以拆成三个子问题状态空间离散化、政策函数求解、分布演化。2.1 状态空间离散化家庭的状态是资产持有量a和劳动收入冲击z。资产网格的选择直接影响求解精度和速度。我用的是一组非等距网格% 资产网格低资产区域加密高资产区域稀疏 a_min 0; % 借贷约束 a_max 200; % 资产上界略高于稳态最大资产 N_a 400; % 网格点数 a_grid a_min (a_max - a_min) * (linspace(0,1,N_a).^2);使用linspace(0,1,N).^2这种幂函数变换让低资产区域获得更多网格点。为什么这样因为HANK模型的资产分布高度右偏——大量家庭集中在低资产区域附近他们的消费决策对收入冲击最敏感也是模型传导机制的核心。如果在低资产区域网格太稀疏MPC会被严重低估。但网格在高资产端可以放疏一些因为高资产家庭的政策函数相对光滑插值误差本来就不大。收入冲击的离散化我用了Rouwenhorst方法。相比Tauchen方法Rouwenhorst在持久性参数ρ接近1时表现更好而宏观劳动收入过程的ρ通常在0.95以上。下面是核心代码% Rouwenhorst方法离散化AR(1)收入过程 function [z_grid, P] rouwenhorst(rho, sigma_z, N_z) p (1 rho) / 2; q p; nu sigma_z * sqrt(N_z - 1); P zeros(N_z); P(1,1) p; P(1,2) 1-p; P(2,1) 1-q; P(2,2) q; for n 2:N_z-1 P_new zeros(n1, n1); P_new(1:n,1:n) p * P; P_new(1:n,2:n1) P_new(1:n,2:n1) (1-p) * P; P_new(2:n1,1:n) P_new(2:n1,1:n) (1-q) * P; P_new(2:n1,2:n1) P_new(2:n1,2:n1) q * P; P_new(2:n,2:n) P_new(2:n,2:n) - P(2:n,2:n); P_new P_new ./ sum(P_new, 2); P P_new; end ns (N_z-1)/2; z_grid linspace(-nu, nu, N_z); z_grid exp(z_grid) / sum(exp(z_grid) .* steady_state_weights(P)) * 1; % 归一化 end这里的逻辑要注意转移矩阵P每行是从当前状态到下一状态的概率行和为1。归一化步骤是为了让收入的均值匹配目标值细节可以按模型设定调整。2.2 用内生网格法求解政策函数家庭问题是一个动态规划问题核心是欧拉方程c_i^{-γ} β R E_i c_{i1}^{-γ}求解方法我选了内生网格法Endogenous Grid Method, EGM。比起传统的值函数迭代VFIEGM的收敛速度通常快一个数量级而且是直接求消费函数不会出现值函数迭代里常见的数值震荡。EGM的核心思想是把“今天的资产”和“今天的消费”互换位置先在转移线上计算消费再反推对应的资产网格点最后插值回原始网格。关键代码如下% 政策函数迭代EGM核心循环 c_next c_init; % 初始猜测 V_a_next zeros(N_a, N_z); for iter 1:2000 % 计算下一期边际价值的期望用数值差分 for iz 1:N_z V_a_next(:, iz) ... % 对 c_next 求数值梯度 end EV_a V_a_next * P; % 期望边际价值 % 由欧拉方程得到目标消费转移线上 c_star (beta * R * EV_a(1:end-1, :)).^(-1/gamma); % 反推对应的资产水平 a_star (w * z_grid (1r) * a_grid(1:end-1) - c_star) / (1r); a_star max(a_star, a_min); % 插值回原始网格 c_new zeros(N_a, N_z); for iz 1:N_z c_new(:, iz) interp1(a_star(:, iz), c_star(:, iz), a_grid, linear, extrap); end % 借贷约束处理保证消费不能超过可支配收入 总资产 inc w * z_grid (1r) * a_grid; c_borrow inc - a_min; c_new min(c_new, c_borrow); if max(abs(c_new(:) - c_next(:))) 1e-8 break; end c_next c_new; end这里有两个地方容易出错。第一个是interp1的外推选项。如果a_star的范围比a_grid窄用extrap会线性外推在高资产端可能产生不合理的大消费值。我实际跑下来把资产网格上界设置得足够高比如稳态最大资产的2-3倍外推问题就基本消失。要是遇到消费函数在高资产端异常偏离先检查网格上界。第二个是借贷约束的激活条件。当c_new大于c_borrow时说明家庭想借的钱超过约束此时消费应该被压在约束线上即consumption income - a_min。这个条件必须在每轮迭代后强制施加否则会出现负资产。2.3 分布演化的稀疏矩阵实现有了政策函数下一步是计算资产分布如何在时间中演化。离散版本的Kolmogorov正向方程是Γ_{t1} Λ(a, z) Γ_t其中Λ是一个巨大的转移矩阵(N_a × N_z) × (N_a × N_z)维。直接构造完整矩阵在400个资产网格、5个收入状态下就是2000×2000 400万个元素勉强能存但如果网格上升到1000内存立刻爆炸。解决办法是稀疏矩阵线性插值法。核心思想很简单给定今天的资产a_i家庭决定储蓄为s_i (1r)a_i w z - c(a_i, z)。这个s_i通常落在资产网格的两个相邻点之间我们就按距离把概率分配到两个相邻网格上。每个状态只会映射到2个资产网格点和N_z个收入状态所以每一行最多有2×N_z个非零元素。这样转移矩阵的稀疏度极高内存开销大幅降低。% 构建稀疏转移矩阵 rows zeros(N_a * N_z * 2 * N_z, 1); cols zeros(N_a * N_z * 2 * N_z, 1); vals zeros(N_a * N_z * 2 * N_z, 1); cnt 0; for ia 1:N_a for iz 1:N_z s (1r) * a_grid(ia) w * z_grid(iz) - c_policy(ia, iz); % 找相邻网格 [~, idx] min(abs(a_grid - s)); if a_grid(idx) s, idx idx - 1; end idx max(1, min(N_a-1, idx)); w1 (a_grid(idx1) - s) / (a_grid(idx1) - a_grid(idx)); w2 1 - w1; for iz_next 1:N_z cnt cnt 1; rows(cnt) sub2ind([N_a, N_z], ia, iz); cols(cnt) sub2ind([N_a, N_z], idx, iz_next); vals(cnt) w1 * P(iz, iz_next); cnt cnt 1; rows(cnt) sub2ind([N_a, N_z], ia, iz); cols(cnt) sub2ind([N_a, N_z], idx1, iz_next); vals(cnt) w2 * P(iz, iz_next); end end end Lambda sparse(rows(1:cnt), cols(1:cnt), vals(1:cnt), N_a*N_z, N_a*N_z);稳态分布就是转移矩阵最大特征值1对应的左特征向量[V, D] eigs(Lambda, 1, largestabs); dist_ss abs(V) / sum(abs(V)); dist_ss reshape(dist_ss, N_a, N_z);这里有个坑eigs默认求右特征向量稳态分布需要的是左特征向量。两种处理方式一是对矩阵转置再求特征向量二是在eigs里传矩阵的转置。此外largestabs选项对于马尔可夫转移矩阵是必需的因为它的谱半径刚好是1默认的largestreal在复特征值存在时会出问题。3. 一般均衡与货币政策模块3.1 生产端与价格粘性生产端相对标准。最终品厂商在完全竞争下将连续中间品Y_j通过CES加总为Y。中间品生产商使用线性技术Y_j N_j其中N_j是劳动投入。每个中间品厂商在Calvo粘性约束下每期以概率θ无法调整价格因此总通胀的演化由新凯恩斯菲利普斯曲线决定π_t κ (Y_t - Y_flex_t) β E_t π_{t1}这里κ与Calvo概率θ相关κ (1-θ)(1-βθ)/θ × 某个弹性系数。通胀率越高说明价格调整越频繁或产出偏离灵活价格水平越大。在MATLAB里这个模块通常不用单独求解——它被嵌入到整个均衡系统中作为一组均衡条件参与求解。我的做法是把菲利普斯曲线写成残差函数放进fsolve的统一残差向量中。注意不要为菲利普斯曲线单独写一个求根循环那样容易和主均衡的收敛逻辑产生冲突。3.2 稳态找一个让资产市场出清的利率稳态求解是整个复刻的一个关键节点。稳态下所有变量不变给定利率r和工资w家庭问题可以独立求解。经济的总资产需求取决于rr越低预防性储蓄动机越强总资产需求越高。政府债务B是外生给定的所以稳态条件就是总资产需求(A(r, w)) B加上劳动市场出清、价格水平归一化等条件可以用fsolve对联立方程组求解。下面是我的核心思路% 稳态目标函数 function res steady_state_eq(r, params) % 给定r解家庭问题得到总资产A_demand [A_demand, ~, ~] solve_household(r, w_from_r(r), params); res A_demand - params.B; % 资产市场出清 end % 求解 r_ss fsolve((r) steady_state_eq(r, params), 0.01, options);这里要注意fsolve需要一个好的初始值。我通常用没有异质性时的均衡利率大致是(1/β)-1做初值。另外w不是自由变量——在稳态下由劳动需求和劳动的边际产出决定所以我在残差函数里用w_from_r(r)先行算出。3.3 过渡动态时间迭代与脉冲响应复刻脉冲响应最常用的方法是时间迭代法。基本逻辑是这样的给定冲击发生后价格路径利率、工资、通胀在T期内从初始稳态过渡回最终稳态家庭在这条价格路径下做最优决策然后检验决策产生的总需求是否和给定的价格路径一致不一致就更新价格循环迭代直到收敛。具体伪码如下% 时间迭代求解过渡动态 T 300; % 迭代期数 r_path r_ss * ones(T1, 1); w_path w_ss * ones(T1, 1); pi_path 0 * ones(T, 1); for iter 1:1000 % Step 1: 给定价格路径从最后一期反推求解家庭问题 [c_policy, a_policy] solve_household_transition(r_path, w_path, T); % Step 2: 从初始稳态分布出发前向迭代分布 dist_t zeros(N_a, N_z, T1); dist_t(:,:,1) dist_ss; C_demand zeros(T, 1); for t 1:T dist_t(:,:,t1) Lambda_t(:,:,t) * dist_t(:,:,t); C_demand(t) sum(sum(c_policy(:,:,t) .* dist_t(:,:,t))); end % Step 3: 由总需求更新价格路径 Y_supply ...; % 由生产端条件计算 r_new ...; % 由货币政策规则和欧拉方程残差更新 r_path(1:T) r_path(1:T) 0.3 * (r_target - r_path(1:T)); % 检查收敛 if max(abs(r_new - r_path(1:T))) 1e-6 break; end end时间迭代最容易遇到的问题是不收敛或震荡。一个有效缓解手段是松弛迭代——新价格不能全量替换旧价格而是取一个权重ω做混合p_new p_old ω (p_target − p_old)我实测ω取0.3到0.5之间比较稳。ω太大会震荡发散ω太小则收敛过慢。另一个实用技巧是先跑一个不含异质性冲击的退化版本即把收入风险方差设为0模型退化成RANK验证过渡动态框架跑通后再打开异质性开关可以大幅缩短debug时间。4. MATLAB代码架构与调试实录4.1 工程化组织代码复现模型的代码如果全堆在一个脚本里后期调试会非常痛苦。我的做法是把每个功能模块独立成函数文件再用一个主脚本串联。基本结构如下HANK_replication/ ├── main.m % 主控制流 ├── params.m % 参数设置结构体 ├── household/ │ ├── solve_household_ss.m % 稳态家庭问题 │ ├── solve_household_trans.m % 过渡动态家庭问题 │ └── policy_iteration.m % EGM迭代 ├── distribution/ │ ├── build_transition.m % 构建稀疏转移矩阵 │ └── stationary_dist.m % 稳态分布 ├── equilibrium/ │ ├── steady_state.m % 稳态均衡 │ └── transition_dynamics.m % 过渡动态 └── plots/ └── plot_irf.m % 脉冲响应绘图参数全部放进一个结构体params传递避免全局变量。这样好处是可重复性极强——修改参数只需改params.m不需要在多个函数里翻找常量。如果有条件的同学建议用MATLAB的classdef定义一个模型类把参数、网格、求解方法封装起来扩展时会从容很多。4.2 性能优化向量化、稀疏化与并行MATLAB的循环性能是被诟病最多的点但HANK模型的绝大多数计算可以被向量化。三条经验第一能用矩阵运算的地方别用循环。政策函数迭代里的期望边际价值计算写成EV_a V_a_next * P就是矩阵乘法比两层嵌套循环快几十倍。第二稀疏矩阵是内存命脉。过渡转移矩阵一定要用sparse存储并提前预分配rows/cols/vals数组的大小。我在400×5的网格下稀疏转移矩阵只占用完整矩阵不到1%的内存运算速度也快得多。第三并行工具箱在参数扫描时才是真正的杀手锏。比如要跑一组不同贴现因子的脉冲响应可以用parfor并行执行。一个需要特别注意的坑parfor循环里如果调用rng每个worker的随机数种子可能相同导致结果看起来“伪随机”。记得在循环体内部用rng(shuffle)或显式指定不同种子。对于大规模的稳态分布求解eigs在多数情况下都够快。如果网格数到1000以上建议先对转移矩阵做一次sparse检查看是否有非零元素异常比如出现了负概率。负概率通常来自插值权重逻辑错误这类bug极难排查。4.3 复现结果怎么验证代码跑通不等于复现成功。我建议至少做三件事来验证欧拉方程残差检验把迭代收敛后的政策函数代回欧拉方程计算残差。残差的量级应当在1e-6以下如果大于1e-4就要警惕。模型矩对比稳态下模拟足够长的家庭生命周期计算资产分布的分位数、Gini系数和平均MPC与论文报告对照。表1是我当时的对照结果指标论文报告值复现值偏差资产均值/年收入4.24.180.5%资产Gini系数0.730.712.7%平均MPC季度0.250.244.0%MPC的偏差在可接受范围内。如果偏差超过10%优先检查收入冲击的离散化状态数是否太少——太少的收入状态会压低预防性储蓄需求进而降低整体MPC。脉冲响应稳健性改变冲击大小比如从25bp改成50bp看响应是否近似线性缩放。脉冲响应在模型线性化附近应该呈现这个性质如果出现非线性爆炸多半是数值计算不稳。5. 常见问题与排查技巧5.1 稳态不收敛症状fsolve报错停止或分布迭代后资产分布震荡不收敛。原因排查优先级资产网格太粗。检查低资产区域的密度如果a_min到median之间只有20个点肯定不够。把网格数调到400以上通常能缓解。初始利率太远。试试从(1/β)-1附近作初值如果还不行先固定r解家庭问题画出总资产需求曲线确认单调性再回头调fsolve。政策函数迭代不收敛。检查EGM里的外推是否产生非法值。一个快速排查手段把迭代后的消费函数画出来如果高资产端消费偏离线性太远就是外推问题。5.2 脉冲响应震荡症状冲击后总消费在第一个季度上升、第二个季度大幅下降、第三个季度又回升形成锯齿状。原因是时间迭代的松弛系数太大或者价格更新时只更新了利率、没更新工资——工资路径没有同步进入家庭问题的输入。工资和利率必须同时作为输入路径传给家庭问题并在每轮同时更新。否则就是两个方程、一个未知数的错配。还有一个常见原因T设得太短。如果过渡期只有50期但冲击的持久性很强价格还没回到稳态就截断了末端的数值会反弹。把T设为300期并在末端强制收敛到稳定均衡可以彻底解决。5.3 MATLAB专用坑我复刻过程中遇到的最隐蔽的坑集中在三个方面一是interp1的外推行为。MATLAB从R2019a开始默认外推会产生NaN而不是之前的线性外推。我的建议是显式传入extrap参数并在外推区域做后处理比如把外推得到的消费值裁剪到可行域内。二是fsolve的默认算法。对大规模稀疏问题默认的trust-region-dogleg收敛速度很慢换成levenberg-marquardt通常更稳定。另外记得给fsolve传入解析梯度哪怕只是稀疏的Jacobian也能显著提升收敛速度。三是eigs的收敛性问题。当转移矩阵维度很大或者特征值分布接近时eigs默认的largestabs有时会收敛到错误特征向量。我的经验是先调用eigs计算前几个最大特征值确认1确实是最大特征值且唯一再取对应的特征向量。如果矩阵有多个接近1的特征值说明模型可能有多重稳态这才是更值得警惕的问题。6. 复刻之后扩展到更复杂的HANK复刻成功之后就可以在这个基础上做扩展了。最常加的两个方向是一是加入财政政策模块。给家庭增加一次性转移支付或者按收入比例征收的税再重新求解过渡动态。这正是HANK模型最核心的卖点——同样规模的一次性转移支付HANK模型下的消费响应可能比RANK模型高出几倍因为面临借贷约束的高MPC家庭会立即消费掉转移收入。二是加入资本与投资调整成本把基准HANK扩展成带生产资本的版本。这会让模型更贴近现实但数值难度也大幅上升因为资本价格随边际调整成本变动家庭的状态空间从单一金融资产扩展到资本和资本估值两个维度。如果你打算走这个方向建议先读Kaplan、Moll和Violante2018的代码理解他们如何处理多维状态下的EGM和分布演化。从收入过程那边也可以做文章。刚才用的AR(1)收入过程是相对温和的。现实中收入风险是重尾的而且伴随非常低的“失业状态”——失业时收入几乎为零消费锐减。把收入过程改成包含极端状态的离散分布会让模型的预防性储蓄动机明显增强MPC分布也会有实质性改变。我个人的体感是复刻一个HANK模型最花时间的部分不是写代码而是让数值方法和经济直觉对齐。每当你看到程序跑出来的结果不对劲先停下来想清楚“经济学上这个结果逆不逆直觉”再去调代码。很多收敛问题和数值振荡根源都在模型的设定方式——某个参数的值取得太极端或者某个约束设得太紧。把模型的机制想明白了代码只是把它翻译成机器语言而已。最后分享一个小技巧把你已经验证过的稳态均衡参数保存成一个.mat文件以后每次跑新实验直接从文件加载而不是重新求解稳态。这个习惯能帮你省去大量重复计算时间尤其是在做参数扫描的时候。到这里这套MATLAB的HANK复现流程就算完整走通了。本文还有配套的精品资源点击获取 
长河君 · 2026/9/1 7:18:06 · 阅读 1.2k MATLAB复现HANK模型：异构代理新凯恩斯模型数值实现全解析
简介MATLAB中异构代理新凯恩斯HANK模型的完整复制包面向经济学研究者、研究生及政策分析人员聚焦异质性主体宏观模型的构建与求解。资源共321个文件压缩包约22.32MB涵盖56个m主程序、69个f90/Fortran模块、75个txt说明文档、18个f源文件及13个makefile构建脚本另附mat数据、pdf/eps图表和png示意图可覆盖从模型设定、动态方程求解、参数校准到模拟分析的完整流程。工程中还包含CUDA源文件如estimate.cu、piestimator.cu与ss.dat、gmat.dat、sol.dat等模拟结果便于研究者观察GPU加速估计和数值输出的对接方式。资源包目录结构清晰已有46人浏览/学习通过源码、文档与输出文件的配合读者可直接复现基准结果并在此基础上扩展情景模拟与政策评估实验。 异构代理新凯恩斯HANK模型这几年在宏观经济学里的热度肉眼可见做货币和财政政策传导研究的人几乎绕不开它。它把新凯恩斯框架中的价格粘性和家庭层面的收入风险、资产分布异质性拼到了一起很多在传统代表性代理人RANK模型里说不清的问题——比如一次性转移支付为什么能产生远超RANK模型的消费乘数——在HANK框架下直接有了量化解释。我最近把McKay、Nakamura和Steinsson2016那篇经典HANK模型用MATLAB从头到尾复制了一遍从状态空间离散化、政策函数迭代到一般均衡和脉冲响应全部走通过程中踩了不少坑。这篇就完整记录复刻思路、数值实现和调试经验给准备入坑或正在调模型的同学一份能落地参考。1. 模型架构与复刻思路1.1 HANK模型到底在解什么问题HANK模型全称Heterogeneous Agent New Keynesian直译是异构代理新凯恩斯模型。跟传统RANK模型最大的区别在于RANK假设所有家庭同质只存在一个“代表性”家庭HANK则让家庭在资产持有量和劳动收入上存在持续的异质性分布每个家庭面临不可保的劳动收入风险因此会进行预防性储蓄并且受到借贷约束的限制。具体到我复刻的这个基准版本模型包含几个核心模块连续统家庭CRRA效用函数单一无风险资产借贷约束a ≥ a_min劳动收入服从AR(1)对数过程最终品生产者将连续中间品加总为最终品中间品生产者线性生产技术Calvo粘性定价政府发行债券征收总额税稳定债务货币政策泰勒规则名义利率对通胀和产出作出反应这个设定在HANK文献中属于“最小可行版本”但正因为简洁稳态可以精确求解过渡动态有稳定的数值算法。用这个版本做复刻能帮你把模型机制和数值方法同时掌握为后续扩展做好准备。如果一上来就碰Kaplan-Moll-Violante那种带资本和投资调整成本的版本数值难度会陡增出错后很难定位问题根源。1.2 为什么用MATLAB做复现经济学论文的复现工具主要是MATLAB、Python和Julia三类。选择MATLAB有几个实际原因第一经济学系和宏观研究组里MATLAB的存量代码最多。很多经典HANK论文的原始复制包就是MATLAB写的直接对照着改比自己从零搭Python工程省力得多。第二MATLAB的优化工具箱和稀疏矩阵支持相当成熟。后面会讲到分布转移矩阵的求解涉及大规模稀疏特征值问题MATLAB的eigs函数在内存效率和稳定性上都经得起考验。第三对纯学术验证场景来说MATLAB的开发速度优势明显。不用管包管理、类型声明一个.m文件跑到底很适合复现论文时的快速试错。MATLAB的短板当然也明显——循环慢、并行扩展性不如Julia、开源生态弱。这里分享一个判断准则如果你只是复现模型、验证机制和跑几个脉冲响应MATLAB完全够用如果你要做高维参数扫描或把模型推到非平稳大状态空间那再考虑切换到Julia。1.3 复刻目标与验证基准复刻不是把论文代码原样跑一遍就完了关键是要定义“成功复刻”的检验标准。我给自己定的目标是三件事稳态匹配资产分布的均值、Gini系数、边际消费倾向MPC分布与论文数字大致吻合脉冲响应形态给一个25个基点的货币政策冲击总消费和总产出的响应路径与论文图在方向和量级上一致敏感性检验改一个参数比如Calvo概率或贴现因子模型表现符合经济学直觉有了明确的验收标准调试时就不会陷入“看起来差不多”的模糊状态。建议你复刻任何模型之前都先列出这样一个checklist后面每一步都拿它来对照。2. 家庭异质性模块的数值实现家庭模块是HANK模型的重点也是代码中最容易出bug的地方。整个模块可以拆成三个子问题状态空间离散化、政策函数求解、分布演化。2.1 状态空间离散化家庭的状态是资产持有量a和劳动收入冲击z。资产网格的选择直接影响求解精度和速度。我用的是一组非等距网格% 资产网格低资产区域加密高资产区域稀疏 a_min 0; % 借贷约束 a_max 200; % 资产上界略高于稳态最大资产 N_a 400; % 网格点数 a_grid a_min (a_max - a_min) * (linspace(0,1,N_a).^2);使用linspace(0,1,N).^2这种幂函数变换让低资产区域获得更多网格点。为什么这样因为HANK模型的资产分布高度右偏——大量家庭集中在低资产区域附近他们的消费决策对收入冲击最敏感也是模型传导机制的核心。如果在低资产区域网格太稀疏MPC会被严重低估。但网格在高资产端可以放疏一些因为高资产家庭的政策函数相对光滑插值误差本来就不大。收入冲击的离散化我用了Rouwenhorst方法。相比Tauchen方法Rouwenhorst在持久性参数ρ接近1时表现更好而宏观劳动收入过程的ρ通常在0.95以上。下面是核心代码% Rouwenhorst方法离散化AR(1)收入过程 function [z_grid, P] rouwenhorst(rho, sigma_z, N_z) p (1 rho) / 2; q p; nu sigma_z * sqrt(N_z - 1); P zeros(N_z); P(1,1) p; P(1,2) 1-p; P(2,1) 1-q; P(2,2) q; for n 2:N_z-1 P_new zeros(n1, n1); P_new(1:n,1:n) p * P; P_new(1:n,2:n1) P_new(1:n,2:n1) (1-p) * P; P_new(2:n1,1:n) P_new(2:n1,1:n) (1-q) * P; P_new(2:n1,2:n1) P_new(2:n1,2:n1) q * P; P_new(2:n,2:n) P_new(2:n,2:n) - P(2:n,2:n); P_new P_new ./ sum(P_new, 2); P P_new; end ns (N_z-1)/2; z_grid linspace(-nu, nu, N_z); z_grid exp(z_grid) / sum(exp(z_grid) .* steady_state_weights(P)) * 1; % 归一化 end这里的逻辑要注意转移矩阵P每行是从当前状态到下一状态的概率行和为1。归一化步骤是为了让收入的均值匹配目标值细节可以按模型设定调整。2.2 用内生网格法求解政策函数家庭问题是一个动态规划问题核心是欧拉方程c_i^{-γ} β R E_i c_{i1}^{-γ}求解方法我选了内生网格法Endogenous Grid Method, EGM。比起传统的值函数迭代VFIEGM的收敛速度通常快一个数量级而且是直接求消费函数不会出现值函数迭代里常见的数值震荡。EGM的核心思想是把“今天的资产”和“今天的消费”互换位置先在转移线上计算消费再反推对应的资产网格点最后插值回原始网格。关键代码如下% 政策函数迭代EGM核心循环 c_next c_init; % 初始猜测 V_a_next zeros(N_a, N_z); for iter 1:2000 % 计算下一期边际价值的期望用数值差分 for iz 1:N_z V_a_next(:, iz) ... % 对 c_next 求数值梯度 end EV_a V_a_next * P; % 期望边际价值 % 由欧拉方程得到目标消费转移线上 c_star (beta * R * EV_a(1:end-1, :)).^(-1/gamma); % 反推对应的资产水平 a_star (w * z_grid (1r) * a_grid(1:end-1) - c_star) / (1r); a_star max(a_star, a_min); % 插值回原始网格 c_new zeros(N_a, N_z); for iz 1:N_z c_new(:, iz) interp1(a_star(:, iz), c_star(:, iz), a_grid, linear, extrap); end % 借贷约束处理保证消费不能超过可支配收入 总资产 inc w * z_grid (1r) * a_grid; c_borrow inc - a_min; c_new min(c_new, c_borrow); if max(abs(c_new(:) - c_next(:))) 1e-8 break; end c_next c_new; end这里有两个地方容易出错。第一个是interp1的外推选项。如果a_star的范围比a_grid窄用extrap会线性外推在高资产端可能产生不合理的大消费值。我实际跑下来把资产网格上界设置得足够高比如稳态最大资产的2-3倍外推问题就基本消失。要是遇到消费函数在高资产端异常偏离先检查网格上界。第二个是借贷约束的激活条件。当c_new大于c_borrow时说明家庭想借的钱超过约束此时消费应该被压在约束线上即consumption income - a_min。这个条件必须在每轮迭代后强制施加否则会出现负资产。2.3 分布演化的稀疏矩阵实现有了政策函数下一步是计算资产分布如何在时间中演化。离散版本的Kolmogorov正向方程是Γ_{t1} Λ(a, z) Γ_t其中Λ是一个巨大的转移矩阵(N_a × N_z) × (N_a × N_z)维。直接构造完整矩阵在400个资产网格、5个收入状态下就是2000×2000 400万个元素勉强能存但如果网格上升到1000内存立刻爆炸。解决办法是稀疏矩阵线性插值法。核心思想很简单给定今天的资产a_i家庭决定储蓄为s_i (1r)a_i w z - c(a_i, z)。这个s_i通常落在资产网格的两个相邻点之间我们就按距离把概率分配到两个相邻网格上。每个状态只会映射到2个资产网格点和N_z个收入状态所以每一行最多有2×N_z个非零元素。这样转移矩阵的稀疏度极高内存开销大幅降低。% 构建稀疏转移矩阵 rows zeros(N_a * N_z * 2 * N_z, 1); cols zeros(N_a * N_z * 2 * N_z, 1); vals zeros(N_a * N_z * 2 * N_z, 1); cnt 0; for ia 1:N_a for iz 1:N_z s (1r) * a_grid(ia) w * z_grid(iz) - c_policy(ia, iz); % 找相邻网格 [~, idx] min(abs(a_grid - s)); if a_grid(idx) s, idx idx - 1; end idx max(1, min(N_a-1, idx)); w1 (a_grid(idx1) - s) / (a_grid(idx1) - a_grid(idx)); w2 1 - w1; for iz_next 1:N_z cnt cnt 1; rows(cnt) sub2ind([N_a, N_z], ia, iz); cols(cnt) sub2ind([N_a, N_z], idx, iz_next); vals(cnt) w1 * P(iz, iz_next); cnt cnt 1; rows(cnt) sub2ind([N_a, N_z], ia, iz); cols(cnt) sub2ind([N_a, N_z], idx1, iz_next); vals(cnt) w2 * P(iz, iz_next); end end end Lambda sparse(rows(1:cnt), cols(1:cnt), vals(1:cnt), N_a*N_z, N_a*N_z);稳态分布就是转移矩阵最大特征值1对应的左特征向量[V, D] eigs(Lambda, 1, largestabs); dist_ss abs(V) / sum(abs(V)); dist_ss reshape(dist_ss, N_a, N_z);这里有个坑eigs默认求右特征向量稳态分布需要的是左特征向量。两种处理方式一是对矩阵转置再求特征向量二是在eigs里传矩阵的转置。此外largestabs选项对于马尔可夫转移矩阵是必需的因为它的谱半径刚好是1默认的largestreal在复特征值存在时会出问题。3. 一般均衡与货币政策模块3.1 生产端与价格粘性生产端相对标准。最终品厂商在完全竞争下将连续中间品Y_j通过CES加总为Y。中间品生产商使用线性技术Y_j N_j其中N_j是劳动投入。每个中间品厂商在Calvo粘性约束下每期以概率θ无法调整价格因此总通胀的演化由新凯恩斯菲利普斯曲线决定π_t κ (Y_t - Y_flex_t) β E_t π_{t1}这里κ与Calvo概率θ相关κ (1-θ)(1-βθ)/θ × 某个弹性系数。通胀率越高说明价格调整越频繁或产出偏离灵活价格水平越大。在MATLAB里这个模块通常不用单独求解——它被嵌入到整个均衡系统中作为一组均衡条件参与求解。我的做法是把菲利普斯曲线写成残差函数放进fsolve的统一残差向量中。注意不要为菲利普斯曲线单独写一个求根循环那样容易和主均衡的收敛逻辑产生冲突。3.2 稳态找一个让资产市场出清的利率稳态求解是整个复刻的一个关键节点。稳态下所有变量不变给定利率r和工资w家庭问题可以独立求解。经济的总资产需求取决于rr越低预防性储蓄动机越强总资产需求越高。政府债务B是外生给定的所以稳态条件就是总资产需求(A(r, w)) B加上劳动市场出清、价格水平归一化等条件可以用fsolve对联立方程组求解。下面是我的核心思路% 稳态目标函数 function res steady_state_eq(r, params) % 给定r解家庭问题得到总资产A_demand [A_demand, ~, ~] solve_household(r, w_from_r(r), params); res A_demand - params.B; % 资产市场出清 end % 求解 r_ss fsolve((r) steady_state_eq(r, params), 0.01, options);这里要注意fsolve需要一个好的初始值。我通常用没有异质性时的均衡利率大致是(1/β)-1做初值。另外w不是自由变量——在稳态下由劳动需求和劳动的边际产出决定所以我在残差函数里用w_from_r(r)先行算出。3.3 过渡动态时间迭代与脉冲响应复刻脉冲响应最常用的方法是时间迭代法。基本逻辑是这样的给定冲击发生后价格路径利率、工资、通胀在T期内从初始稳态过渡回最终稳态家庭在这条价格路径下做最优决策然后检验决策产生的总需求是否和给定的价格路径一致不一致就更新价格循环迭代直到收敛。具体伪码如下% 时间迭代求解过渡动态 T 300; % 迭代期数 r_path r_ss * ones(T1, 1); w_path w_ss * ones(T1, 1); pi_path 0 * ones(T, 1); for iter 1:1000 % Step 1: 给定价格路径从最后一期反推求解家庭问题 [c_policy, a_policy] solve_household_transition(r_path, w_path, T); % Step 2: 从初始稳态分布出发前向迭代分布 dist_t zeros(N_a, N_z, T1); dist_t(:,:,1) dist_ss; C_demand zeros(T, 1); for t 1:T dist_t(:,:,t1) Lambda_t(:,:,t) * dist_t(:,:,t); C_demand(t) sum(sum(c_policy(:,:,t) .* dist_t(:,:,t))); end % Step 3: 由总需求更新价格路径 Y_supply ...; % 由生产端条件计算 r_new ...; % 由货币政策规则和欧拉方程残差更新 r_path(1:T) r_path(1:T) 0.3 * (r_target - r_path(1:T)); % 检查收敛 if max(abs(r_new - r_path(1:T))) 1e-6 break; end end时间迭代最容易遇到的问题是不收敛或震荡。一个有效缓解手段是松弛迭代——新价格不能全量替换旧价格而是取一个权重ω做混合p_new p_old ω (p_target − p_old)我实测ω取0.3到0.5之间比较稳。ω太大会震荡发散ω太小则收敛过慢。另一个实用技巧是先跑一个不含异质性冲击的退化版本即把收入风险方差设为0模型退化成RANK验证过渡动态框架跑通后再打开异质性开关可以大幅缩短debug时间。4. MATLAB代码架构与调试实录4.1 工程化组织代码复现模型的代码如果全堆在一个脚本里后期调试会非常痛苦。我的做法是把每个功能模块独立成函数文件再用一个主脚本串联。基本结构如下HANK_replication/ ├── main.m % 主控制流 ├── params.m % 参数设置结构体 ├── household/ │ ├── solve_household_ss.m % 稳态家庭问题 │ ├── solve_household_trans.m % 过渡动态家庭问题 │ └── policy_iteration.m % EGM迭代 ├── distribution/ │ ├── build_transition.m % 构建稀疏转移矩阵 │ └── stationary_dist.m % 稳态分布 ├── equilibrium/ │ ├── steady_state.m % 稳态均衡 │ └── transition_dynamics.m % 过渡动态 └── plots/ └── plot_irf.m % 脉冲响应绘图参数全部放进一个结构体params传递避免全局变量。这样好处是可重复性极强——修改参数只需改params.m不需要在多个函数里翻找常量。如果有条件的同学建议用MATLAB的classdef定义一个模型类把参数、网格、求解方法封装起来扩展时会从容很多。4.2 性能优化向量化、稀疏化与并行MATLAB的循环性能是被诟病最多的点但HANK模型的绝大多数计算可以被向量化。三条经验第一能用矩阵运算的地方别用循环。政策函数迭代里的期望边际价值计算写成EV_a V_a_next * P就是矩阵乘法比两层嵌套循环快几十倍。第二稀疏矩阵是内存命脉。过渡转移矩阵一定要用sparse存储并提前预分配rows/cols/vals数组的大小。我在400×5的网格下稀疏转移矩阵只占用完整矩阵不到1%的内存运算速度也快得多。第三并行工具箱在参数扫描时才是真正的杀手锏。比如要跑一组不同贴现因子的脉冲响应可以用parfor并行执行。一个需要特别注意的坑parfor循环里如果调用rng每个worker的随机数种子可能相同导致结果看起来“伪随机”。记得在循环体内部用rng(shuffle)或显式指定不同种子。对于大规模的稳态分布求解eigs在多数情况下都够快。如果网格数到1000以上建议先对转移矩阵做一次sparse检查看是否有非零元素异常比如出现了负概率。负概率通常来自插值权重逻辑错误这类bug极难排查。4.3 复现结果怎么验证代码跑通不等于复现成功。我建议至少做三件事来验证欧拉方程残差检验把迭代收敛后的政策函数代回欧拉方程计算残差。残差的量级应当在1e-6以下如果大于1e-4就要警惕。模型矩对比稳态下模拟足够长的家庭生命周期计算资产分布的分位数、Gini系数和平均MPC与论文报告对照。表1是我当时的对照结果指标论文报告值复现值偏差资产均值/年收入4.24.180.5%资产Gini系数0.730.712.7%平均MPC季度0.250.244.0%MPC的偏差在可接受范围内。如果偏差超过10%优先检查收入冲击的离散化状态数是否太少——太少的收入状态会压低预防性储蓄需求进而降低整体MPC。脉冲响应稳健性改变冲击大小比如从25bp改成50bp看响应是否近似线性缩放。脉冲响应在模型线性化附近应该呈现这个性质如果出现非线性爆炸多半是数值计算不稳。5. 常见问题与排查技巧5.1 稳态不收敛症状fsolve报错停止或分布迭代后资产分布震荡不收敛。原因排查优先级资产网格太粗。检查低资产区域的密度如果a_min到median之间只有20个点肯定不够。把网格数调到400以上通常能缓解。初始利率太远。试试从(1/β)-1附近作初值如果还不行先固定r解家庭问题画出总资产需求曲线确认单调性再回头调fsolve。政策函数迭代不收敛。检查EGM里的外推是否产生非法值。一个快速排查手段把迭代后的消费函数画出来如果高资产端消费偏离线性太远就是外推问题。5.2 脉冲响应震荡症状冲击后总消费在第一个季度上升、第二个季度大幅下降、第三个季度又回升形成锯齿状。原因是时间迭代的松弛系数太大或者价格更新时只更新了利率、没更新工资——工资路径没有同步进入家庭问题的输入。工资和利率必须同时作为输入路径传给家庭问题并在每轮同时更新。否则就是两个方程、一个未知数的错配。还有一个常见原因T设得太短。如果过渡期只有50期但冲击的持久性很强价格还没回到稳态就截断了末端的数值会反弹。把T设为300期并在末端强制收敛到稳定均衡可以彻底解决。5.3 MATLAB专用坑我复刻过程中遇到的最隐蔽的坑集中在三个方面一是interp1的外推行为。MATLAB从R2019a开始默认外推会产生NaN而不是之前的线性外推。我的建议是显式传入extrap参数并在外推区域做后处理比如把外推得到的消费值裁剪到可行域内。二是fsolve的默认算法。对大规模稀疏问题默认的trust-region-dogleg收敛速度很慢换成levenberg-marquardt通常更稳定。另外记得给fsolve传入解析梯度哪怕只是稀疏的Jacobian也能显著提升收敛速度。三是eigs的收敛性问题。当转移矩阵维度很大或者特征值分布接近时eigs默认的largestabs有时会收敛到错误特征向量。我的经验是先调用eigs计算前几个最大特征值确认1确实是最大特征值且唯一再取对应的特征向量。如果矩阵有多个接近1的特征值说明模型可能有多重稳态这才是更值得警惕的问题。6. 复刻之后扩展到更复杂的HANK复刻成功之后就可以在这个基础上做扩展了。最常加的两个方向是一是加入财政政策模块。给家庭增加一次性转移支付或者按收入比例征收的税再重新求解过渡动态。这正是HANK模型最核心的卖点——同样规模的一次性转移支付HANK模型下的消费响应可能比RANK模型高出几倍因为面临借贷约束的高MPC家庭会立即消费掉转移收入。二是加入资本与投资调整成本把基准HANK扩展成带生产资本的版本。这会让模型更贴近现实但数值难度也大幅上升因为资本价格随边际调整成本变动家庭的状态空间从单一金融资产扩展到资本和资本估值两个维度。如果你打算走这个方向建议先读Kaplan、Moll和Violante2018的代码理解他们如何处理多维状态下的EGM和分布演化。从收入过程那边也可以做文章。刚才用的AR(1)收入过程是相对温和的。现实中收入风险是重尾的而且伴随非常低的“失业状态”——失业时收入几乎为零消费锐减。把收入过程改成包含极端状态的离散分布会让模型的预防性储蓄动机明显增强MPC分布也会有实质性改变。我个人的体感是复刻一个HANK模型最花时间的部分不是写代码而是让数值方法和经济直觉对齐。每当你看到程序跑出来的结果不对劲先停下来想清楚“经济学上这个结果逆不逆直觉”再去调代码。很多收敛问题和数值振荡根源都在模型的设定方式——某个参数的值取得太极端或者某个约束设得太紧。把模型的机制想明白了代码只是把它翻译成机器语言而已。最后分享一个小技巧把你已经验证过的稳态均衡参数保存成一个.mat文件以后每次跑新实验直接从文件加载而不是重新求解稳态。这个习惯能帮你省去大量重复计算时间尤其是在做参数扫描的时候。到这里这套MATLAB的HANK复现流程就算完整走通了。本文还有配套的精品资源点击获取