查看全球股票指数数据的每日收盘价

原始数据由2665个对以下代表性股票指数日收盘价的观察数据组成，时间跨度为1993年4月27日至2003年7月14日:

加拿大:多伦多证券综合指数(股票代码^GSPTSE)法国:CAC 40(股票代码^FCHI)德国:DAX(股票代码^GDAXI)日本:日经225(股票代码^N225)英国:富时100(股票代码^FTSE)美国:标准普尔500(股票代码^GSPC)

负载Data_GlobalIdx1导入每日指数关闭数

下面的图表说明了每个指数的相对价格走势。每个指数的初始水平已归一化，以方便相对业绩的比较，并且没有明确考虑股息调整。

图plot(日期，ret2price(price2ret(数据)))“x”)包含(“日期”) ylabel (的索引值)标题(“相对日指数收盘价”)传说(系列,“位置”，“西北”）

为后续建模做准备，将每个指数的收盘水平转换为每日对数回报(有时称为几何回报或连续复合回报)。

returns = price2ret(数据);对数回报率%T = size(返回，1);回报率% #(即历史样本量)

由于整体建模方法的第一步涉及反复应用GARCH过滤和极值理论来表征每个股票指数回报系列的分布，因此对特定国家的细节进行检查是有帮助的。您可以将下一行代码更改为集合{1,2,3,4,5,6}中的任何整数，以检查任何索引的详细信息。

指数= 1;1 =加拿大，2 =法国，3 =德国，4 =日本，5 =英国，6 =美国图(日期(2:结束)，返回(:，索引))“x”)包含(“日期”) ylabel (“返回”)标题(“每日对数收益”）

过滤每个索引的返回值

用GPD对分布尾部建模要求观测值近似独立和同分布(i.i.d)。然而，大多数财务收益序列表现出一定程度的自相关，更重要的是表现出异方差。

例如，与所选指数相关的收益的样本自相关函数(ACF)揭示了一些温和的序列相关性。

图autocorr(returns(:，index))“收益ACF样本”）

然而，平方收益的样本ACF说明了方差的持久性程度，并表明GARCH建模可能显著地影响后续尾部估计过程中使用的数据。

图autocorr(returns(:，index).^2)“收益平方的样本ACF”）

为了产生一系列i.i.d.观察，将一阶自回归模型拟合到每个股票指数的回报的条件均值

条件方差的不对称GARCH模型

一阶自回归模型补偿自相关性，而GARCH模型补偿异方差。特别地，最后一项通过布尔指标将不对称性(杠杆)合并到方差中，如果先验模型残差为负，则取值1，否则为0(参见[3]）.

此外，每个指数的标准化残差被建模为标准化的Student's t分布，以补偿通常与股权回报相关的肥尾。这是

下面的代码段从每个股票指数的收益中提取经过过滤的残差和条件方差。

模型= arima(基于“增大化现实”技术的南,“分布”，“t”，“方差”gjr (1,1));nindexes = size(数据，2);指数% #残差= NaN(T, nIndices);预分配存储%方差= NaN(T, nIndices);fit = cell(nIndices,1);选项= optimoptions(@fmincon，“显示”，“关闭”，.．.“诊断”，“关闭”，“算法”，“sqp”，“TolCon”1 e);为i = 1: nindexes fit{i} =估计(模型，返回(:，i)，“显示”，“关闭”，“选项”、选择);(残差(:,i)、方差(:,我)]=推断(适合{我},返回(:,我));结束

对于选定的指数，比较从原始收益中过滤出来的模型残差和相应的条件标准差。下图清楚地说明了在过滤残差中存在的波动率(异方差)的变化。

图subplot(2,1,1) plot(日期(2:结束)，残差(:，索引))datetick(“x”)包含(“日期”) ylabel (“残留”)标题(过滤后的残差的) subplot(2,1,2) plot(日期(2:结束)，sqrt(方差(:，索引))datetick(“x”)包含(“日期”) ylabel (“波动”)标题(“过滤条件标准差”）

从每个返回序列中过滤出模型残差后，用相应的条件标准差对残差进行标准化。这些标准化残差代表了样本CDF尾的EVT估计所基于的潜在零均值、单位方差、i.i.d.系列。

残差=残差。/平方根(方差);

为了结束本节，请检查标准化残差和平方标准化残差的acf。

将标准化残差的acf与原始收益的相应acf进行比较，可以发现标准化残差现在大约是i.i.d，因此更容易接受后续的尾部估计。

图autocorr(残差(:，index))“标准化残差的ACF样本”)图autocorr(残差(:，指数).^2)标题(标准化残差平方的样本ACF）

估计半参数cfs

给定上一步的标准化i.i.d.残差，用高斯核估计每个指标的经验CDF。这平滑了CDF估计，消除了未平滑样本CDF的阶梯模式。尽管非参数核CDF估计非常适合于分布的内部，其中大部分数据被发现，但当应用于上尾和下尾时，它们往往表现不佳。为了更好地估计分布的尾部，将EVT应用于每个尾部的残差。

具体来说，找到上下限阈值，使每个尾保留10%的残差。然后将每个尾部的极端残差超出相关阈值的量以最大似然拟合为参数GPD。这种方法通常被称为违规情况分布或峰值超过阈值方法。

给定每个尾部的超过量，优化负对数似然函数来估计GPD的尾部指数(zeta)和尺度(beta)参数。

下面的代码段创建类型的对象paretotails，每个索引返回序列都有一个这样的对象。这些paretotails对象封装了参数GP下尾、非参数核平滑内部和参数GP上尾的估计，为每个指标构建了一个复合半参数CDF。

得到的分段分布对象允许在CDF内部进行插值，并允许在每个尾部进行外推(函数计算)。外推是非常可取的，允许在历史记录之外估计分位数，并且对于风险管理应用程序是非常宝贵的。

此外，Pareto尾对象还提供了计算CDF和逆CDF(分位数函数)的方法，以及查询分段分布各段之间边界的累积概率和分位数的方法。

nPoints = 200;CDF的每个区域中抽样点的% #tailFraction = 0.1;分配给每个尾的残差的十进制分数tails = cell(nIndices,1);帕累托尾对象的单元格数组为i = 1:nIndices尾{i} = paretotails(残差(:，i)， tailFraction, 1 - tailFraction，“内核”）;结束

在估计了复合半参数经验CDF的三个不同区域后，以图形方式拼接并显示结果。再次注意，分别显示为红色和蓝色的下尾和上尾区域适合外推，而内核光滑的内部(黑色)适合内插。

下面的代码使用拟合所基于的数据以外的数据调用感兴趣的Pareto尾部对象的CDF和逆CDF方法。具体来说，所引用的方法可以访问拟合状态，现在调用它们来选择和分析概率曲线的特定区域，作为强大的数据过滤机制。

图保存在网格在minProbability = cdf(尾巴{指数}，(min(残差(:，指数))));maxProbability = cdf(尾巴{指数}，(max(残差(:，指数))));pLowerTail = linspace(minProbability, tailFraction, nPoints);%样品下尾pUpperTail = linspace(1 - tailFraction, maxProbability, nPoints);上尾%样品pInterior = linspace(tailFraction, 1 - tailFraction, nPoints);样品内部百分比plot(icdf(tails{index}， pLowerTail)， pLowerTail，“红色”，“线宽”， 2) plot(icdf(tails{index}， pInterior)， pInterior，“黑”，“线宽”， 2) plot(icdf(tails{index}， pUpperTail)，“蓝”，“线宽”， 2) xlabel(标准化的剩余的) ylabel (“概率”)标题(“经验提供”)({传奇“帕累托下尾”“内核平滑内部”.．.“帕累托上尾”}，“位置”，“西北”）

评估GPD契合度

虽然前面的图表说明了复合CDF，但更详细地检查GPD拟合是有指导意义的。

GP分布的CDF参数化为

对于超过(y)，尾部指数参数(zeta)和规模参数(beta)。

为了直观地评估GPD拟合，绘制残差上尾超过的经验CDF以及由GPD拟合的CDF。虽然只使用了10%的标准化残差，但拟合分布紧跟超过数据，因此GPD模型似乎是一个很好的选择。

figure [P,Q] = border (tails{index});%累积概率&边界上的分位数y = sort(residuals(residuals(:，index) > Q(2)， index) - Q(2));排序超出百分比情节(y) (cdf(反面{指数},y + Q (2)) - P (2)) / P (1)) (F (x) = ecdf (y);经验CDF %持有在楼梯(x, F,“r”网格)在传奇(拟合广义帕累托CDF，“经验提供”，“位置”，“东南”）;包含(“超过数”) ylabel (“概率”)标题(“标准化残差的上尾”）

校准t Copula

给定标准化残差，现在估计标量自由度参数(DoF)和线性相关矩阵(R)的t联结使用copulafit统计和机器学习工具箱中的函数。的copulafit函数允许用户用两种不同的方法估计t联结的参数。

默认的方法是在两步过程中执行的正式最大似然方法。

具体而言，虽然可以直接使全对数似然最大化，但目标函数的最大值往往出现在多维空间中一个长而平的窄的山谷中，在高维空间中可能出现收敛困难。为了克服这些困难，copulafit最大似然估计(MLE)分两步进行。在给定一个固定的自由度值的情况下，内步使关于线性相关矩阵的对数似然最大化。这个条件最大化被放置在一个1-D的自由度最大化范围内，从而使所有参数的对数似然最大化。在外层步骤中被最大化的函数被称为自由度的轮廓对数似然函数。

相比之下，下面的代码段使用了一种替代方法，它近似于大样本量自由度参数的剖面对数似然值。虽然这种方法通常比MLE快得多，但应谨慎使用，因为对于小或中等样本量，估计和置信限可能不准确。

具体地说，近似是通过对对数似然函数与线性相关矩阵求导而得到的，假设自由度是一个固定常数。由此得到的非线性方程，然后迭代求解相关矩阵。反过来，这个迭代嵌套在另一个优化中，以估计自由度。这种方法在精神上类似于上面的轮廓对数似然方法，但不是真正的最大似然，因为对于给定的自由度，相关矩阵不会收敛到条件最大似然。然而，对于大样本量，估计值往往非常接近最大似然值(参见[1]而且[7]）.

还有一种方法通过首先计算样本秩相关矩阵(Kendall's tau或Spearman's rho)来估计线性相关矩阵，然后用稳健的正弦变换将秩相关转换为线性相关。给定线性相关矩阵的这个估计，对数似然函数就单独的自由度参数最大化。这种方法似乎也受到剖面对数似然方法的激励，但根本不是MLE方法。然而，这种方法不需要矩阵反演，因此具有在接近奇异的相关矩阵存在时数值稳定的优势(参见[8]）.

最后是Nystrom和Skoglund[6]建议保留自由度参数作为用户指定的模拟输入，从而允许用户主观上推断资产之间的尾依赖程度。特别是，他们推荐一个相对较低的自由度，介于1和2之间，因此可以更仔细地检查联合极端的可能性。这种方法对于极端相互依赖程度至关重要的压力测试是有用的。

下面的代码段首先通过上面推导的半参数经验CDF将标准化残差转换为均匀变量，然后将t copula拟合到转换后的数据。当均匀变量通过各裕度的经验CDF进行转换时，校准方法通常被称为典型最大似然(CML)。

U =零(大小(残差));为i = 1: nIndices U (:, i) = cdf(反面{我},残差(:,我));将边距转换为均匀结束[R, DoF] = copulafit(“t”U“方法”，“approximateml”）;%适合连接符

用t Copula模拟全球指数投资组合收益

给定t关联的参数，现在通过首先模拟相应的相关标准化残差来模拟联合相关的股票指数收益。

方法模拟依赖的均匀变量copularnd统计和机器学习工具箱中的函数。

然后，通过外推到GP尾部，插值到平滑的内部，通过反演每个指标的半参数边缘CDF，将均匀变量转换为标准化残差。这产生了与上述AR(1) + GJR(1,1)滤波过程所得的结果一致的模拟标准化残差。这些残差在时间上是独立的，但在任何时间点上都是相关的。模拟的标准化残差数组的每一列在孤立地观察时表示一个i.i.d单变量随机过程，而每一行共享由关联符引起的秩相关。

下面的代码段模拟了在22个交易日的一个月期限内对相关标准化指数残差进行的2000次独立随机试验。

s = RandStream.getGlobalStream();reset(s) nTrials = 2000;% #的独立随机试验地平线= 22;% VaR预测范围Z =零(地平线，nTrials, nIndices);%标准化残差数组U =共生(“t”， R, DoF，视界* nTrials);% t联结模拟为j = 1: nIndices Z (:,:, j) =重塑(icdf(反面{j}, U (:, j)),地平线,nTrials);结束

使用模拟的标准化残差作为i.i.d.输入噪声过程，通过计量经济学工具箱™重新引入原始指数回报中观察到的自相关和异方差过滤器函数。

请注意,过滤器接受来自copula的用户指定的标准化干扰，并一次为单个索引模型模拟多条路径，这意味着为每个索引连续模拟和存储所有样本路径。

为了充分利用当前信息，请指定必要的预样本模型残差、方差和返回值，以便每个模拟路径都从共同的初始状态演化而来。

Y0 =返回(end，:);%预采样回报率Z0 =残差(end，:);%预采样标准化残差V0 =方差(end，:);%预样本方差simulatedReturns =零(地平线，nTrials, nIndices);为i = 1: nIndices simulatedReturns (:,:, i) =过滤器(适合{我},Z(:,:我),.．.“Y0”, Y0(我),“Z0”Z0(我),“半”V0 (i));结束

现在，重新构造模拟的返回数组，使每个页面都表示对多元返回序列的一次尝试，而不是对单变量返回序列的多次尝试。

simulatedReturns = permute(simulatedReturns， [1 3 2]);

最后，给出每个指数的模拟收益，形成由各个指数组成的权重相等的全球指数组合(国家指数的全球指数)。由于我们处理的是每日对数收益，在风险范围内的累计收益只是每个介入期间的收益之和。还要注意，投资组合的权重在整个风险范围内都是固定的，并且模拟忽略了重新平衡投资组合所需的任何交易成本(假设每日的重新平衡过程是自融资的)。

请注意，虽然模拟的收益是对数的(连续复利)，但投资组合收益系列是通过首先将单个对数收益转换为算术收益(价格变化除以初始价格)，然后对单个算术收益加权以获得投资组合的算术收益，最后转换回投资组合的对数收益来构造的。对于每日数据和较短的VaR水平，重复转换的差异很小，但对于较长的时间周期，差异可能很大。

累计回报=零(nTrials, 1);weights = repmat(1/nIndices, nIndices, 1);%等加权组合为i = 1:nTrials cumulativeReturns(i) = sum(log(1 + (exp(simulatedReturns(:，:，i)) - 1) * weights));结束

总结结果

模拟了每个指数的收益，形成了全球投资组合，报告了一个月风险范围内的最大收益和损失，以及在各种置信水平下的VaR。此外，绘制全球投资组合的经验CDF。

VaR = 100 *分位数(累积回报率，[0.10 0.05 0.01]');disp (' ')流('最大模拟损耗:%8.4f%s\n', -100 *分钟(cumulativeReturns),“%”)流('最大模拟增益:%8.4f%s\n\n', 100 * max (cumulativeReturns),“%”)流('模拟90%% VaR: %8.4f%s\n'VaR (1),“%”)流('模拟95%% VaR: %8.4f%s\n'VaR (2),“%”)流('模拟99%% VaR: %8.4f%s\n\n'VaR (3),“%”)图h = cdfplot(累计回报);h.Color =“红色”；包含(“对数回归”) ylabel (“概率”)标题(模拟一个月全球投资组合收益CDF）

最大模拟损失:14.4218%最大模拟增益:12.4470%模拟90% VaR: -4.7476%模拟95% VaR: -6.3291%模拟99% VaR: -10.6639%