
开头做金融时间序列的风险建模绕不开两个让人头疼的问题第一波动率到底该怎么估计和预测第二组合里多个资产之间的相关性到底该怎么刻画。标题里把Copulas、GARCH、EWMA、EqWMA这些模型放在一起再用CVaR、极值理论、蒙特卡洛模拟做市场风险管理这套组合拳我在Matlab上完整踩过一遍今天把整个实现思路和代码心得拆开讲。这篇内容适配两类人一类是刚接触量化风控的研究生想搭一套可运行的组合风险度量框架另一类是在金融机构做风险归因的从业者想在VaR基础上把尾部风险、资产相依结构这层东西补上。项目本身解决的核心问题很简单——单资产的风险度量已经烂大街了但当你手里握着十几条资产收益率序列想回答“组合在极端行情下到底能亏多少”的时候Copula加上蒙特卡洛模拟这条路是目前最实用的选择。1. 方案选型为什么是CopulaGARCH蒙特卡洛1.1 从“单资产风险”到“组合风险”的跳跃很多刚入门的朋友都写过单资产的GARCH模型用Matlab的estimate函数拟合一下然后在条件方差基础上算VaR一步到位。但真实场景里你管理的从来不是一只股票而是一篮子资产。这里立刻冒出一个问题组合风险不等于每个资产风险的简单加总。直接从历史收益率序列里算组合的相关系数矩阵看起来简单实际坑很多——金融收益率数据普遍存在波动率聚集效应相关系数在不同的市场状态下会漂移尤其在市场大跌时资产之间的相关性会急剧上升这就是常说的“相关性的脆弱性”。如果沿用恒定相关系数下尾风险会被严重低估。Copula模型的价值就在于把“边缘分布”和“相依结构”拆开处理。我可以先给每一条资产序列配一个GARCH模型把波动率特征、厚尾特征拟合得漂漂亮亮再用Copula函数把多条序列之间的尾部依赖关系接起来。这正是Copula能站在现代风险管理舞台中央的原因。1.2 Matlab做这件事的优势Matlab在这个项目里有两个不可替代的方便之处。第一Econometrics Toolbox提供了garch、egarch等现成的模型类配合estimate、infer、simulate一套流程省去了自己写极大似然估计迭代求解的功夫。第二Statistics and Machine Learning Toolbox里的copulafit和copularnd函数让Copula参数拟合和模拟变得极其顺手我自己早期用Python写这套流程时Copula拟合往往要借助第三方库Matlab这边则是开箱即用。当然别把Matlab当黑盒。用arma2ma转换、用filter滤波这些底层操作还是要自己掌控的。我觉得最好的定位是把Matlab当作一个高效的实验平台业务逻辑和模型判断还是得自己在代码里一条一条写清楚。1.3 项目核心流程概览整个项目的运行路径非常清晰数据从Excel读入先做对数收益率变换接着对每条资产序贯建立GARCH等波动率模型提取标准化残差然后把这些残差序列通过copulafit拟合为Copula参数接着开始蒙特卡洛模拟生成Copula场景下的联合分布随机数再反变换为标准化残差乘以波动率预测值得到未来收益路径最后统计组合损失分布并计算CVaR。2. 波动率模型GARCH、EWMA、EqWMA三选一2.1 三种模型背后的数学逻辑先聊GARCH(1,1)公式长这样$$ \sigma_t^2 \omega \alpha \varepsilon_{t-1}^2 \beta \sigma_{t-1}^2 $$它的核心思想是“今天的波动率由昨天的收益率冲击和昨天的波动率共同决定”。$\alpha$度量新信息的影响$\beta$度量记忆性两者相加小于1时模型平稳。金融资产日收益率序列里这种“大波动之后还是大波动”的聚集效应用GARCH拟合效果极好。EWMA模型则是RiskMetrics的老牌方法形式上可以写成$$ \sigma_t^2 \lambda \sigma_{t-1}^2 (1-\lambda) r_{t-1}^2 $$这里的$\lambda$在RiskMetrics框架里通常取0.94它本质上是一个指数衰减权重模型比GARCH少一个$\omega$项因此波动率永远不会均值回复。你要是用EWMA做长期预测会发现预测波动率逐渐压向一个极低的值这在实际风控中是个麻烦事。EqWMA就是等权重移动平均窗口内每个观测的权重一样$$ \sigma_t^2 \frac{1}{m} \sum_{i1}^{m} r_{t-i1}^2 $$这个模型很朴素胜在稳健不会有参数估计误差。它没有记忆衰减的概念窗口外的数据完全不参与窗口内的数据一视同仁适合样本量不多但需要快速给出基线估计的场景。2.2 Matlab里拟合GARCH的实操细节Matlab里拟合GARCH(1,1)非常直接% 假设 logRet 为对数收益率列向量 MdL garch(1,1); [EstMdL,EstParamCov,logL] estimate(MdL, logRet); % 提取条件波动率序列 [sig, logL] infer(EstMdL, logRet);这里有个小细节estimate函数默认采用极大似然估计但金融收益率的分布往往是厚尾的直接假设正态分布会让参数估计偏掉。可以指定Distribution为t让模型估计时同时优化自由度参数MdL garch(1,1); MdL.Distribution t;我实测下来工作日数据学生t分布的GARCH(1,1)组合样本外VaR的表现比纯正态假设好不少。另外提一嘴infer返回的条件波动率序列长度和原序列相同但前几期可能因为缺少滞后数据返回NaN做后续处理时记得先去掉这些NaN行。2.3 EWMA与EqWMA的Matlab实现EWMA没有工具箱函数手写也就三五行lambda 0.94; var_t(1) var(logRet); % 初始方差 for t 2:length(logRet) var_t(t) lambda * var_t(t-1) (1-lambda) * logRet(t-1)^2; end sigma_EWMA sqrt(var_t);EqWMA更简单用movmean函数就搞定了window 30; % 30日等权重窗口根据实际需求调 eq_sigma movstd(logRet, [window-1 0]) ; % 或者用 movvar 再开根号注意movstd窗口参数给的是[window-1 0]表示只取当前时刻及往回window-1个样本这样不会引入未来数据。任何风控模型里“不用未来数据”这条原则在代码里就必须落实到每一个窗口函数上。2.4 模型选择和滚动预测对三种模型做样本外对比我习惯采用滚动窗口预测的方式用过去250个交易日估计模型预测下一个交易日的波动率依次滚动得到一长串预测波动率序列。然后计算预测误差常用的指标是均方根误差RMSE或者用QMLE提出的异方差调整后的误差% 假设 realizedVol 是由日收益率平方开根号得到的代理值 MRMSE sqrt(mean((pred_sigma.^2 - realizedVol.^2).^2));实操中我通常把三种模型的预测结果叠在一张图上看谁更贴近后续真实波动率的走势。GARCH因为能捕捉均值回复中期预测更稳EWMA对近期突变反应快但长期预测会偏离EqWMA则是三者中最平滑的适合做基准线。3. Copula模型关键参数与拟合流程3.1 常用Copula族与适用场景Copula函数可以理解成一个连接函数把多个边缘分布“粘”起来从而得到一个联合分布。Matlab支持高斯Copula、t Copula和几类阿基米德CopulaClayton、Gumbel、Frank。它们的区别在于刻画相依结构的“侧重点”不同高斯Copula对称相依尾部渐近独立适合相关性较弱、没有明显极值共动特征的资产组合t Copula对称且具有尾部相依性同时捕捉上下尾的同时极端值对金融资产适配度相当高Clayton Copula上尾低、下尾高适合债市或股票市场中“暴跌时齐跌”的特征Gumbel Copula上尾高、下尾低更多用于捕捉牛市共涨的场景Frank Copula对称但尾部趋近独立适合相关性温和且无明显极值倾向的数据。我自己的做法是先对标准化残差做经验CDF变换之后把高斯Copula、t Copula、Clayton、Gumbel都试一遍看AIC和BIC哪个最小。实战中t Copula往往常胜——真实金融数据的尾部相依确实存在而t Copula在拟合尾部方面比高斯好一个量级。3.2 copulafit的两种参数估计路径在Matlab中拟合Copula一般走两步搞定。第一步把每条序列的标准化残差用经验累积分布函数转成均匀分布U(0,1)上的变量。第二步用copulafit估计Copula参数u zeros(size(resid)); for i 1:size(resid,2) u(:,i) ksdensity(resid(:,i), resid(:,i), function,cdf); end % 也可以使用 ecdf注意处理尾部边界 rho copulafit(Gaussian, u); % 或 [rho, nu] copulafit(t, u);这里有个容易踩的坑ksdensity做经验CDF变换后极值会被“压扁”到接近0或1而这种极小/极大的值在某些Copula似然函数里会让数值优化失效。常见解法是给u做一层截断比如限定在(0.001, 0.999)之间边界外就不管了。另一种是用ecdf结合pseudo-observation的思路直接通过排序位置转换成均匀分布。3.3 自由度的作用与拟合优度t Copula比高斯Copula多一个自由度参数$\nu$。$\nu$越小尾部越厚极端事件同时发生的概率越大当$\nu$趋向无穷时t Copula退化为高斯Copula。copulafit(t, u)返回的就是相关矩阵和自由度标量。判断模型好不好不要只看似然值我会做一次Bootstrap重抽样的稳定性检验。简单讲就是把原始数据有放回地抽样N次每次重新拟合Copula参数看参数的变异程度。如果自由度和相关矩阵波动很大说明拟合不稳定模型的样本外可靠性自然存疑。3.4 为什么“先GARCH、再Copula”能避免伪相关直接拿原始收益率做Copula会因为每条资产波动率结构不同而得到虚假的相关性。举个例子两只股票在同一时期的波动率都升高但它们本身可能毫无业务交集仅仅是整个市场的波动率环境抬升了计算出来的Pearson相关系数就会虚高。先通过GARCH对每条序列做“降噪”把标准化残差提取出来此时得到的是纯粹地去除了异方差效应的冲击信号再用Copula刻画这些冲击信号之间的相依关系相关性度量才真正反映了资产之间的本质关联。4. 极值理论与CVaR的衔接4.1 CVaR到底在算什么CVaR的全称是条件风险价值也叫期望损失Expected Shortfall。它在统计上定义为损失分布中超过VaR阈值的条件期望$$ CVaR_\alpha E[L | L VaR_\alpha] $$说人话就是“如果最坏的那5%场景真的发生了平均会亏多少”。这比VaR信息量足因为VaR只告诉你1%分位数的亏损数字却没告诉你一旦跌破这个数损失会有多惨。CVaR把尾部的平均损失答清楚监管和内部风控当中应用也更顺滑。在蒙特卡洛模拟框架下CVaR的计算非常简单——有了N次模拟的组合亏损序列把亏损从高到低排序取前5%或按置信水平调整的样本求平均即可。4.2 EVT在边缘分布建模中的角色EVT解决的是“尾部到底有多厚”的问题。拟合边缘分布时我虽然可以用t分布来拟合标准化残差但t分布对整个分布的描述还是“伞形”的尾部的精确形状未必刻画到位。这时EVT出场。常用做法是POT超阈值法选定一个高阈值$u$超过$u$的超出量用广义帕累托分布GPD拟合$$ F(y) 1 - \left(1 \frac{\xi y}{\beta}\right)^{-\frac{1}{\xi}} $$在Matlab中可以使用gpfit和gpinv做GPD的参数估计与分位数计算。阈值选择是关键——阈值太高则超出样本太少估计方差大阈值太低则又混入非极端的部分模型偏误大。可以选择以90%到95%的分位数作为初始阈值的起点然后画一下“超出量均值随阈值变化的曲线”找到线性区域的起点。实际把EVT嵌入Copula框架可以在GARCH的标准化残差边缘分布上用EVT建模再用Copula连接尾部也可以在Copula模拟完成之后对模拟得到的组合损失序列用EVT做CVaR外推。两种路径我都试过第二种操作起来更容易效果也不差。4.3 风险因子映射与维度选择“风险因子”这个词在标题里出现不是偶然。实际做组合风控的时候资产数量可能很多但真正驱动组合波动的因子往往只有几个——利率因子、股票市场因子、商品因子、汇率因子等。把几十条资产收益率降维成几个风险因子然后用这些因子的联合分布做蒙特卡洛模拟模拟速度更快模型稳定性也更好。Matlab的pca函数可以快速做主成分分析通过前几个主成分解释大部分方差。这一步处理得好Copula模型的维度从30维降到5维计算量和过拟合风险都会显著下降。5. 蒙特卡洛模拟的Matlab实现5.1 模拟路线图蒙特卡洛模拟在整个项目里是“算账”的最后一步。整体模拟线路如下从拟合好的Copula中生成N条相依均匀变量U(0,1)用每条资产的边缘分布的分位数函数把U变换成标准化残差结合GARCH模型预测的未来条件标准差得到未来一日的模拟收益率按组合权重加权得到未来一日的组合收益重复N轮得到组合损失的完整分布从分布中提取VaR和CVaR。5.2 关键代码与参数选择用t Copula模拟的代码如下Nsim 10000; u_sim copularnd(t, rho, nu, Nsim); % Nsim x d 的相依均匀变量 % 假设边缘分布是t分布根据GARCH标准化残差拟合出的t自由度 simResid zeros(size(u_sim)); for i 1:d simResid(:,i) tinv(u_sim(:,i), df_resid(i)); end % 预测未来一天的条件波动率以GARCH为例 [forecastSigma, ~] forecast(EstMdL, 1, Y0, logRet); % 模拟收益率矩阵 simRet simResid .* forecastSigma; % 组合模拟收益 w ones(d,1)/d; portfolio_ret simRet * w; % 组合损失取负收益为正 loss -portfolio_ret; var_95 quantile(loss, 0.95); cvar_95 mean(loss(loss var_95));这里要特别留意forecast函数的使用。forecast(EstMdL, 1, Y0, logRet)表示基于当前全部历史数据预测下一期的方差。如果你用的是EWMA或者EqWMA那你需要自己把最后一条条件方差存下来然后手动递推下一步。5.3 模拟次数与收敛判断10万次模拟在Matlab里是几秒到几十秒的事但很多研究者只做1万次省时间却可能不收敛。判断标准很简单把模拟批次增加到20万重新算CVaR如果和10万次的差不超过1%左右说明结果已经稳定如果还在跳就得加次数或者检查是不是随机数种子的问题。rng的种子设置也很重要。合规的蒙特卡洛模拟应该能复现否则第二天打开代码算出不同的数字同事会怀疑你的模型是不是有bug。项目里建议统一在一开始执行rng(2025)固定随机种子。5.4 方差缩减的几个实用技巧蒙特卡洛模拟最让人烦心的就是“坏运气”——恰好抽到一个特别极端的样本导致CVaR失真。实测中两个办法比较有效第一是拉丁超立方抽样保证样本在联合分布空间里分布更均匀Matlab里可以用lhsdesign不过要转成有相依结构的样本需要额外操作第二是对偶变量法把每一次随机数取反生成两条路径让收益率分布两侧对称对降低方差效果极好。金融场景下对偶变量法实现起来几乎零成本强烈建议试试。6. 实例演算用沪深300成分股做一个简化组合6.1 数据准备与预处理考虑到合规和数据可得性我用沪深300里选5只不同行业成分股作为例子演示流程。数据区间选最近三年的日收益率Label对齐到同一个交易日历。读入方式T readtable(daily_returns.xlsx); ret T{:, 2:end}; % 每列是一支股票的日收益率 % 剔除包含NaN的行 ret(any(isnan(ret),2),:) [];对数收益率的计算建议直接用对数差值logRet diff(log(price_table{:, 2:end}),1,1);如果拿到的数据本身就是收益率序列那就直接用。实证中差别不大但理论上对数收益率的时间加总性质更好。6.2 序贯拟合与残差提取对5条资产序列我写一层循环挨个处理d size(logRet, 2); sigMat zeros(size(logRet)); residMat zeros(size(logRet)); dfVec zeros(1, d); for i 1:d MdL garch(GARCHLags,1,ARCHLags,1,Distribution,t); EstMdL estimate(MdL, logRet(:,i), Display,off); [sig, ~] infer(EstMdL, logRet(:,i)); resid logRet(:,i) ./ sig; sigMat(:,i) sig; residMat(:,i) resid; % 记录t分布自由度 dfVec(i) EstMdL.Distribution.DoF; end这一步得到三样东西条件波动率序列、标准化残差序列、每只股票的t分布自由度。标准化残差应当近似独立同分布可以用autocorr函数做个简单白噪声检验如果残差还有显著自相关说明GARCH阶数不够得考虑GARCH(2,1)或者加入ARMA均值方程。6.3 拟合Copula与模拟标准化残差矩阵准备好后转成均匀变量然后拟合t Copulafor i 1:d [u(:,i), ~] ecdf(residMat(:,i)); u(u1, i) 1 - 1e-6; % 防止边界1导致优化失败 u(u0, i) 1e-6; % 防止边界0 end [rho, nu] copulafit(t, u);这里插一句ecdf返回的最后一项必然是1必须手动缩一下否则copulafit会弹出“数据必须在开区间(0,1)”的报错。模拟之后组合权重设为等权算出N次模拟下的组合收益分布排序之后计算95% VaR和CVaR。一套流程跑完结果可能显示95% VaR是-2.1%CVaR是-3.4%意思就是最坏5%的情景里平均会亏3.4%。如果只看单资产VaR再线性加总这个数字往往小得多多出来的那部分差就是因为资产尾部相依性带来的真实风险。6.4 结果比较仅用历史模拟 vs Copula蒙特卡洛要验证Copula方法加进来是否有效果可以把Copula蒙特卡洛模拟得到的CVaR与历史模拟法得到的CVaR做个对比。历史模拟就是从历史已实现收益序列中直接取分位数缺点在于历史数据中可能根本没出现过极端联合亏损所以会低估风险。Copula方法的额外价值在于即使样本期内没出现过的极端组合情景也能通过相依结构“造”出来提前把尾部暴露查清楚。7. 常见问题与排查技巧实录7.1 copulafit报错“非有限参数”或“NaN”遇到这个情况八成是均匀变换后出现了0或1。无论如何换一种边缘分布拟合方式或者在调用copulafit之前打印一下min(u(:))和max(u(:))一眼定位。7.2 estimate GARCH模型时梯度过大或估计时间过长先把Display,off关上默认为off然后检查是不是初始值设得不好。可以用简单矩估计给出初始值MdL garch(GARCHLags,1,ARCHLags,1); MdL.Constant var(logRet)*0.1; MdL.GARCH{1} 0.7; MdL.ARCH{1} 0.2;实操中这样给初值能大幅缩短迭代时间也能避免陷入局部极值。7.3 t Copula拟合时自由度发散如果copulafit返回的$\nu$数值极大比如超过200说明数据实际上更接近高斯Copula直接用高斯Copula替代即可不丢精度还省一个参数。7.4 蒙特卡洛模拟结果和直觉差太远模拟次数不够最常见。另一个可能是Copula的相关矩阵用了非正定的估计值导致copularnd生成的随机数莫名其妙。用nearestSPD或者correct_rho (rho rho)/2修正一下对称性再用eig检查最小特征值是否大于0。7.5 标准化残差分层看是否还有波动率聚集模拟之前务必检查标准化残差是否白噪声。写一段小代码对残差平方做Ljung-Box检验h lbqtest(resid.^2, Lags, 10);如果h1说明还有自相关需要回到GARCH模型阶数选择上重做。7.6 计算CVaR时置信水平的选择95%置信水平是最常用的但巴塞尔协议框架下现在更强调97.5%或99%尾部度量。置信水平越高尾部样本越少CVaR估计的方差就越大。应对措施就是增加模拟次数至少5万次或者结合EVT对尾部概率分布做平滑外推。8. 代码框架的整体封装思路8.1 用一个主函数串联所有环节项目写到后面建议不要把所有代码堆在一个脚本里。用四个函数模块把流程拆开fitVolModels.m输入收益率矩阵输出三种波动率模型参数、标准化残差矩阵与预测波动率fitCopulaModel.m输入标准化残差输出Copula参数与拟合优度统计量runMonteCarlo.m输入Copula参数、边缘分布参数、组合权重输出组合损失分布与VaR/CVaRplotRiskReport.m把三种波动率模型下的CVaR结果和损失分布画在一起。这样做的最大好处是换了数据只需要改第一段读数据的代码后面不用动。8.2 关于Matlab版本与工具箱我用的是MATLAB R2023agarch、copulafit、copularnd在R2018a之后的版本都稳定可用。如果你没有Econometrics Toolbox可以考虑用最简实现方案GARCH模型自己写一个简单的极大似然估计迭代Copula自己写一个高斯Copula的解析表达式。但如果条件允许直接用工具箱更省事数值性能也可靠。8.3 扩展方向动态Copula与时变参数目前用的Copula参数是静态的一次性估计得到固定相关矩阵。实际市场中资产之间的相关结构是随时间漂移的。可以做滚动窗口的Copula参数估计每个时间点都用过去500个交易日重新拟合一次Copula参数画一条相关性的动态曲线。实测中你可以明显看到市场波动率抬升的区间相关矩阵的均值会同步上升。这条动态曲线对于压力测试非常有价值能直观说明“危机时相关性上升”这一风险传导现象。9. 每次必查的几个细节最后的最后分享几个我自己每次跑这套流程都会固定检查的点。数据齐不平齐是一个重要检查项——多资产组合里最常见的坑是不同资产的交易日历不一致比如A股和港股、美股混在一起节假日差异直接导致收益率矩阵里出现大量NaN很多新手在这里吃哑巴亏。处理办法是不要用fillmissing填出假数据直接用intersect把共同交易日提取出来宁少勿假。另一个是风险因子与组合映射逻辑。如果你的组合不是简单的等权组合而是有明确的持仓市值那必须把持仓权重向量准确传进蒙特卡洛环节权重加起来要等于1且不能有负权重除非允许做空。我曾经在回测脚本里漏改了一次权重导致两根指数收益序列的CVaR差了30%排查了半天才发现是权重矩阵粘贴时错了一行。再有一个是模型对比的一致性。GARCH、EWMA、EqWMA三个模型做预测对比时样本外区间必须完全一致预测起点必须一样否则结果没有可比性。我在项目里统一设定前60%的数据用来训练后40%做样本外测试滚动窗口长度固定为250个交易日这样三张图放在一起才讲得清谁优谁劣。遇到大样本的数据函数前面加tic和toc观测一下各步耗时也很有必要。GARCH拟合一般很快最耗时的往往是大维度的copulafit如果资产维度超过20建议先用PCA降维或者用Algorithm,trust-region调整优化算法来提速。翻过这些坑之后整套项目跑下来最大的感触是金融模型的精度提升很多时候不靠换一个更复杂的模型而靠把每个环节的假设梳理干净、把数据处理的细节抠到位。Copula也好EVT也好它们都是工具真正决定风控报告是否靠谱的还是你对数据和流程的把控。