ARTICLE DETAIL

资讯详情

深耕编程入门与网站建设的一线实战洞察。

基于MATLAB的GARCH与Realized GARCH模型实现与对比

基于MATLAB的GARCH与Realized GARCH模型实现与对比 简介本资源是一套面向金融工程与量化分析学习者的Realized GARCH波动率建模MATLAB实现代码适用于具备基础时间序列知识和MATLAB编程能力的高年级本科生、研究生及初级量化从业者用于解决传统GARCH模型对日内波动信息利用不足的问题。压缩包共4个文件全部为.m脚本总计2KB其中主程序hw_code.m负责整体流程调度f1.m、f2.m、f3.m分别承担数据预处理、Realized GARCH核心估计与模型诊断功能结构清晰、模块分工明确便于理解高频波动率融入条件方差建模的技术路径。目前已有230人学习下载可直接运行复现Hansen Lunde2005提出的Realized GARCH框架涵盖日内实际波动率计算、非线性似然估计、参数稳定性检验等关键环节是掌握现代波动率建模实践的重要教学级参考代码。 上周帮人看一个 hw.zip差点没把我绕晕。里边的目录乱七八糟但核心东西就两件用 MATLAB 估计 GARCH(1,1)还有 Realized GARCH。这种作业在金融计量课程里太常见了但很多同学卡在“到底怎么让代码跑起来”这一步要么是数据读不进去要么是优化器报 NaN要么是搞不清 Realized GARCH 和普通 GARCH 的区别。今天把这套东西从头到尾拆一遍包括模型公式、MATLAB 实现、优化器选择以及我实际踩过的几个坑。如果你正在写类似的作业或者刚开始接触已实现波动率模型这篇文章应该能帮你省下不少折腾的时间。先说一下背景。hw.zip 这类文件我拿到的版本是一个典型的金融计量课程作业包里面没有太多工程化代码更多是“把模型跑通、把图画出来、把结果解释清楚”。所以下面所有内容我都围绕这个目标来讲。你不会看到复杂的 C 或 Python 工程架构只有 MATLAB 脚本、目标函数、优化器以及怎么避开那些让新手崩溃的细节。1. 项目整体设计与思路拆解1.1 hw.zip里装了什么从文件名反推项目结构先别急着解压跑代码我拿到任何作业包第一件事是看文件命名和目录结构。hw.zip 这个名字看起来随意但里面的文件通常能反映整套作业的逻辑。我这次看到的是这么一组东西readme.txt课程说明、数据来源、应提交哪些图表data.csv 或 data.xlsx日收益率序列、已实现方差序列estimate_garch.mGARCH 模型估计主脚本estimate_realized_garch.mRealized GARCH 估计主脚本garch_likelihood.m、realized_garch_likelihood.m对数似然函数realized_variance.m用高频数据计算已实现方差的函数。这种结构本身就在暗示一件事先估计一个普通 GARCH 作为基准再引入高频信息估计 Realized GARCH最后对比波动率预测效果。理解了这层关系你就不会把两个模型的代码混在一起改到崩溃。从作业角度看这个 hw.zip 要解决的核心问题可以提炼成三句话第一用 MATLAB 自己写 GARCH 类模型的似然函数而不是只调用现成的garch函数第二理解已实现测度Realized Measure如何在条件方差方程里起作用第三能对估计结果做简单的预测和解释。所以下面我按这三个目标展开。1.2 为什么用Realized GARCH而不是普通GARCH传统 GARCH(1,1) 的流行程度不用多说几乎所有金融计量课都会讲。它用日收益率的平方去更新条件方差但日收益率平方是个噪声非常大的“信号”一个异常大的日收益会让方差估计跳一下然后又慢慢衰减。真实波动率的持续性往往被这种噪声掩盖。Realized GARCH 的思路是用日内高频数据构造一个更精确的已实现方差Realized Variance, RV并把它直接放进条件方差方程。你可以简单理解成普通 GARCH 是“用昨天的坏消息和昨天的波动预测今天”Realized GARCH 是“用昨天更真实的波动观测来修正预测”。Hansen、Huang 和 Shek 在 2012 年提出这个模型后很多实证研究都发现它在样本外预测上明显优于传统 GARCH。这个作业之所以要同时写两个模型本质上是让学生亲眼看到加入已实现测度之后波动率持续性参数比如 β会下降因为大部分波动信息已经由 RV 直接提供了条件方差方程不再需要把“昨天的方差”拖得那么久。这个对比是写报告时的核心亮点。1.3 为什么选MATLAB实现有些同学会问这种模型用 R 的rugarch或者 Python 的arch库不是更方便吗确实方便但很多金融计量课程仍然指定 MATLAB原因有三点。一是 MATLAB 的 Optimization Toolbox 里fmincon、fminunc非常成熟处理小规模约束优化很稳定二是 Econometrics Toolbox 自带garch、estimate、forecast但 Realized GARCH 没有现成函数必须自己写这逼着学生去理解模型细节三是 MATLAB 画图、出表、写报告在学术圈有很长的使用传统特别是老一代教师非常习惯。我做项目时也不会盲目排斥用 MATLAB。虽然这种语言在工程上被吐槽不少但做这种几百个观测值的小规模波动率模型它完全够用而且调试起来很直观。你只需要注意一点不要把所有代码堆在一个脚本里最好把似然函数、约束函数、数据预处理分开这样后面排查问题会轻松很多。2. 核心细节解析与实操要点2.1 GARCH(1,1)模型与似然函数先复习一下最标准的 GARCH(1,1)r_t μ ε_tε_t σ_t z_tz_t ~ N(0,1)σ_t^2 ω α ε_{t-1}^2 β σ_{t-1}^2约束条件一般是 ω 0α ≥ 0β ≥ 0且 α β 1保证过程平稳。这里的 ε_{t-1} 是 t-1 期的收益率去均值后的残差也就是“新闻冲击”。在 MATLAB 里自己写估计最核心的部分是对数似然函数。假设扰动项服从正态分布单个观测的对数似然是log L_t -0.5 log(2π) - 0.5 log(σ_t^2) - 0.5 ε_t^2 / σ_t^2把所有时期的对数似然加总取负数就得到我们需要最小化的目标函数。为什么取负数因为 MATLAB 的优化器默认是求最小值。递推的时候要注意初始值。σ_1^2 不能是 0否则对数似然直接 NaN。我一般用样本方差作为初始条件方差ε_1 用第一个样本点的去均值残差。你可以在代码里加一个“预热”阶段比如用前三期平均但样本方差通常就够用。2.2 Realized GARCH模型的结构Realized GARCH 比普通 GARCH 多了一条“测量方程”。我用对数形式这是 Hansen 论文里的标准写法r_t μ ε_tε_t σ_t z_tlog σ_t^2 ω β log σ_{t-1}^2 γ log x_{t-1}log x_t ξ φ log σ_t^2 τ(z_t) u_t其中 τ(z_t) τ_1 z_t τ_2 (z_t^2 - 1)u_t ~ N(0, σ_u^2)。第一个方程是收益率方程第二个是条件方差方程第三个是测量方程它把已实现方差 x_t 和潜在的条件方差 σ_t^2 联系起来。τ(z_t) 的引入是为了捕捉杠杆效应当 z_t 为负也就是坏消息时τ(z_t) 会系统性地让已实现方差变大。为什么用 log 形式因为 x_t 一定大于 0log 变换后不需要额外加参数约束来保证方差为正数值上稳定很多。这条非常重要很多同学的 Realized GARCH 怎么调都不收敛就是因为用了线性形式 σ_t^2 ω β σ_{t-1}^2 γ x_{t-1}然后被方差的非负约束折磨得不行。估计 Realized GARCH 的似然函数由两部分组成收益率部分的密度加上已实现方差测量方程部分的密度。也就是说每一期都有两个“观测”参与估计收益 r_t 和已实现方差 x_t。这比普通 GARCH 信息量更大参数也更多。2.3 已实现方差RV的计算已实现方差的定义很直接把第 t 个交易日按固定时间间隔分成 M 段计算每段对数收益率的平方和RV_t Σ_{j1}^{M} r_{t,j}^2这里的 r_{t,j} 是日内第 j 个收益率。最常用的是 5 分钟频率。频率太高微观结构噪声会污染结果频率太低又会丢失日内波动信息。5 分钟是经验平衡点但实际做作业的时候老师一般会直接给你一个算好的 RV 序列或者给一份分钟数据让你自己算。如果要用 MATLAB 自己算有几个细节必须注意第一只计算交易时段内的收益率不要把隔夜收益直接混进当天的 RV否则会引入隔夜信息第二分钟数据里经常有空缺需要先按时间戳对齐第三如果某天数据太少比如只有 10% 的观测那这一天的 RV 可能严重低估建议在 Excel 或 MATLAB 里标记缺失不要硬算。3. 实操过程与核心环节实现3.1 数据准备与预处理拿到手的数据通常是一个 Excel 或 CSV里面至少三列日期、日收益率、已实现方差。我习惯先readtable读进来再做合法性检查% load_data.m data readtable(data.csv); dates data.Date; ret data.ret; % 日收益率小数形式比如 0.0123 表示 1.23% rv data.rv; % 已实现方差正数 % 清理缺失、非有限值 valid isfinite(ret) isfinite(rv) (rv 0); dates dates(valid); ret ret(valid); rv rv(valid); % 基础统计确认数据没有异常 disp([样本量: , num2str(length(ret))]); disp([收益率均值: , num2str(mean(ret))]); disp([RV 最小值: , num2str(min(rv))]);这一步最容易犯的错误是没检查rv 0。如果 RV 序列里有 0 或者负数后面取对数就会得到-Inf整个优化直接崩溃。宁可删掉几个观测也不要让模型在非法数据上跑。如果你的作业只给了日收益率没有 RV那就要先自己算一个示例 RV 序列。最粗糙的做法是用日收益率的平方代替但那样做 Realized GARCH 的意义就不大了。更合理的是找分钟数据调用realized_variance.m计算。3.2 MATLAB实现GARCH(1,1)估计我用的是fmincon因为要处理 α β 1 这个约束。先写目标函数% garch_likelihood.m function nll garch_likelihood(theta, r) mu theta(1); omega theta(2); alpha theta(3); beta theta(4); T length(r); e r - mu; sigma2 zeros(T, 1); sigma2(1) var(r); % 初始方差用样本方差 for t 2:T sigma2(t) omega alpha * e(t-1)^2 beta * sigma2(t-1); end % 负对数似然正态分布 nll sum(0.5 * log(2*pi) 0.5 * log(sigma2) 0.5 * e.^2 ./ sigma2); end然后写一个非线性约束函数% garch_constraint.m function [c, ceq] garch_constraint(theta) c theta(3) theta(4) - 1; % 要求 alpha beta 1 ceq []; end主脚本里设置初始值和边界% run_garch.m r ret; % 从数据文件读入 theta0 [0; 0.05 * var(r); 0.10; 0.85]; lb [-Inf; 0; 0; 0]; ub [ Inf; Inf; 1; 1]; options optimoptions(fmincon, ... Display, iter, ... Algorithm, interior-point, ... MaxIterations, 1000); [theta_hat, nll_hat] fmincon((th) garch_likelihood(th, r), theta0, ... [], [], [], [], lb, ub, garch_constraint, options); % 提取参数 mu_hat theta_hat(1); omega_hat theta_hat(2); alpha_hat theta_hat(3); beta_hat theta_hat(4); fprintf(mu %.6f\n, mu_hat); fprintf(omega %.6f\n, omega_hat); fprintf(alpha %.6f\n, alpha_hat); fprintf(beta %.6f\n, beta_hat);这里初始值非常重要ω 设置为 0.05 × 样本方差α 设成 0.1β 设成 0.85这是金融日收益率序列的典型经验值。如果你随便设成 1 和 1优化器很容易撞到约束边界然后报错。我试过用fminunc跑但如果 α β 的约束被违反递推出来的 σ_t^2 会发散成非常大的数最终对数似然变成 NaN优化器直接放弃。所以建议老老实实用带约束的fmincon。3.3 MATLAB实现Realized GARCH估计Realized GARCH 的目标函数要复杂一些因为每一期的对数似然由两部分组成。我按对数形式写% realized_garch_likelihood.m function nll realized_garch_likelihood(theta, r, x) mu theta(1); omega theta(2); beta theta(3); gamma theta(4); xi theta(5); phi theta(6); tau1 theta(7); tau2 theta(8); log_sigma_u theta(9); sigma_u2 exp(log_sigma_u); % 保证测量方程扰动方差为正 T length(r); e r - mu; log_sigma2 zeros(T, 1); log_x log(x); % 初始条件方差 log_sigma2(1) log(var(r)); % 递推条件方差 for t 2:T log_sigma2(t) omega beta * log_sigma2(t-1) gamma * log_x(t-1); end z e ./ sqrt(exp(log_sigma2)); % 收益率部分对数似然 ll_r -0.5 * log(2*pi) - 0.5 * log_sigma2 - 0.5 * z.^2; % 测量方程残差 u log_x - xi - phi * log_sigma2 - tau1 * z - tau2 * (z.^2 - 1); % 测量方程对数似然 ll_x -0.5 * log(2*pi) - 0.5 * log(sigma_u2) - 0.5 * u.^2 / sigma_u2; nll -sum(ll_r ll_x); end主脚本里同样设置参数边界% run_realized_garch.m x rv; % 已实现方差正数 % 可以先跑一遍普通GARCH用普通GARCH的mu/omega/beta作为初值 theta0 [mu_hat; 0.05; 0.85; 0.10; 0.1; 1.0; 0; 0.1; log(0.5)]; lb [-Inf; -Inf; 0; 0; -Inf; 0; -Inf; -Inf; log(0.0001)]; ub [ Inf; Inf; 1; 1; Inf; 2; Inf; Inf; log(10)]; options optimoptions(fmincon, ... Display, iter, ... Algorithm, interior-point, ... MaxIterations, 2000); [theta_rg, nll_rg] fmincon((th) realized_garch_likelihood(th, r, x), ... theta0, [], [], [], [], lb, ub, [], options);参数初值有个技巧先用普通 GARCH 估计出 ω、β再套到这里来。γ 初始 0.1 表示已实现方差的滞后影响。φ 初始 1.0因为测量方程里log x_t和log σ_t^2应该接近同尺度。log_sigma_u初始log(0.5)大体对应测量方程残差标准差μ_u 的方差不能设成 0否则测量方程似然会爆炸。有一点要注意我这里的实际使用了“对数方差”递推所以 ω、β、γ 的解释和线性 GARCH 不一样。写报告时不要拿普通 GARCH 的持续性概念直接套上去但你可以看 β 和 γ 的相对大小来讨论“已实现信息”的作用。3.4 结果输出与波动率预测估计完成之后通常还要把条件方差序列画出来并做一期预测。GARCH(1,1) 的下一期条件方差预测是σ_{T1}^2 ω α ε_T^2 β σ_T^2Realized GARCH 的下一期预测则需要测量方程帮助更新log σ_{T1}^2 ω β log σ_T^2 γ log x_T但因为log x_T是已知的所以可以直接代入。也可以用测量方程预测 x_{T1} 的均值但作业一般只要求预测波动率。画图代码很简单% plot_volatility.m sigma2_garch garch_volatility(theta_hat, ret); % 函数里返回条件方差序列 sigma2_rg realized_garch_volatility(theta_rg, ret, rv); figure; plot(dates, sqrt(sigma2_garch), LineWidth, 1); hold on; plot(dates, sqrt(sigma2_rg), LineWidth, 1); hold off; legend(GARCH(1,1), Realized GARCH); title(Conditional Volatility Comparison); ylabel(Volatility);如果你用的是 MATLAB 2022b 之后版本图表的中文和标题处理会更方便但注意保存图片时用exportgraphics而不是旧的print避免字体丢失。4. 常见问题与排查技巧实录4.1 优化不收敛与NaN问题这是最让人崩溃的。症状通常是fmincon跑了几步之后直接退出提示 “Objective function returned NaN/Inf”。原因大概率是递推过程中出现了非正方差或者对数似然函数的某个部分出现了 NaN。排查思路很简单先在目标函数里加一行检查把异常时的参数值打出来。if any(~isfinite(nll)) || any(sigma2 0) fprintf(非有限值: theta ); disp(theta); keyboard; % 进入调试模式 end根据我的经验常见触发点有三个初始方差设成 0 或负数α、β 初值太大递推几期就爆炸数据里有缺失值导致mean(r)或var(r)返回 NaN。解决办法就是数据先清洗初始值按经验值设给参数加边界如果还不行考虑把目标函数里的 θ 做一个变换比如用log(alpha)代替alpha这样优化器永远不会生成负数。4.2 RV序列处理细节RV 序列看起来简单实际坑很多。我第一次跑 Realized GARCH 的时候结果里面的 φ 估计是负的怎么看都不合理。后来发现原始数据里的 RV 是“已实现波动率”标准差形式不是“已实现方差”我直接取了 log导致模型全乱。这里要注意模型里的 x_t 是已实现方差单位是收益率平方。如果数据文件里给的是已实现波动率Realized Volatility需要先平方再代入模型。另外RV 序列里如果有个别极端值比如 0.01 这类正常水平的几十倍也会让优化结果不稳定。可以尝试对极端值做 winsorize或者用中位数过滤。还有一点计算 RV 时用多少频率不是随便定的。你可以画一个波动率签名图volatility signature plot横轴是采样频率纵轴是平均 RV看哪个频率附近曲线开始平稳。5 分钟到 10 分钟之间通常比较稳定。如果作业没有要求直接用 5 分钟即可。4.3 MATLAB版本与工具箱问题如果你用的 MATLAB 版本比较老比如 R2020a 之前有些函数可能不存在。最常见的是readtable对 CSV 的读取差异以及optimoptions的写法。老版本用optimset但代码逻辑是一样的。如果连 Optimization Toolbox 都没有那就麻烦了。不过大多数学校机房会装全套工具箱。如果你在自己电脑上装注意安装时勾选 Optimization Toolbox否则fmincon会提示未定义函数。我之前帮人排查过一个fmincon未定义的问题最后发现是安装 MATLAB 时只装了基础模块工具箱没选。要判断自己有没有工具箱执行license(test, Optimization_Toolbox)返回 1 表示可用返回 0 表示没装或没激活。这个方法比翻安装界面快得多。5. 一点个人经验最后说点写作业之外的经验。我自己做这种模型最大的感受是不要把 Realized GARCH 当成一个“高级版 GARCH”去背诵而要理解它到底解决了什么问题。传统 GARCH 只用日收益平方这一种信息Realized GARCH 把日内高频数据压缩成一个更可靠的观测值并把它显式放进模型。这个思想一旦想通后面再接触 GARCH-M、EGARCH、TGARCH 都只是在这个框架上加加减减。实际操作中还有个建议先跑普通 GARCH拿到参数后再作为 Realized GARCH 的初始值。不要一上来直接跑 Realized GARCH因为它的参数更多似然面更复杂初始值不合适很容易掉进局部最优。你可以用多组初值试一遍选对数似然最小的那个这比任何“高级调参”都管用。最后再分享一个小技巧写报告时不要只贴代码和参数一定要把 GARCH 模型的条件方差序列和 Realized GARCH 的条件方差序列画在同一张图上然后把 αβ 和 βγ 这类持续性强度的变化写清楚。老师看到你能解释“加入已实现测度后波动率持续性下降”这个现象基本就愿意给高分了。我自己当初也是靠这个点拿到了作业答辩里的不少认可。本文还有配套的精品资源点击获取
返回列表