R中的单变量和多变量随机波动率建模

结果如图3.14所示：建模RPeriods ahead5 10 15 20中的随机波动率-0.6-0.2 0.2 0.6图3：多步预测分布（实心、灰色和黑色）和观测值（虚线、红色）。R> obs\\u CHF<-100*logret（exrates$CHF[pred\\u ind]）R>qs<-t（apply（predy（pred\\u svl），2，分位数，c（0.05，0.5，0.95）））R>ts.plot（cbind（qs，obs\\u CHF），xlab=“Periods ahead”，lty=c（rep（1，3），2），+col=c（“gray80”，“black”，“gray80”，“red”））4.4。滚动窗口估计由ugarchroll（）在R包rugarch（Ghalanos 2020）中提出，我们引入包装函数svsample\\u roll（）、svtsample\\u roll（）、svlsample\\u roll（）和svtlsample\\u roll（），分别围绕其相应的例程svsample（）、svtsample（）、svlsample（）和svtlsample（）构建，以简化SVmodels的滚动窗口估计。在这种估计方法中，要么是固定宽度的时间窗口在时间序列中移动，要么是一系列扩展的时间窗口，具有相同的起始时间点，覆盖更大、更大的观测数据块，并且在所有时间窗口中独立估计同一模型。接下来，使用每个估计模型进行样本外预测，通常在时间窗口前一天到一周。最后，预测值集可用于评估模型。在贝叶斯统计中，评估模型预测能力的自然方法是通过其后验预测分布。其密度也称为预测密度，定义为asp（yt+1 | yo[1:t]）=ZKp（yt+1 | yo[1:t]，κ）p（κ| yo[1:t]）dκ，（12），其中κ收集所有未观察到的变量，即在SVtl最常见的情况下，κ=（u，Д，σ，ρ，ν，h，β）>，积分域K是κ的所有可能值集。Wefollow Geweke和Amisano（2010）在我们的符号中使用上标o表示观测值的向量yo[1:t]=（y，y，…，yt）>。

方程12可以看作是所有参数后验分布上预测可能性的积分，因此它解释了预测值的后验参数不确定性。方程12中的积分没有闭合形式，其维数随t增加而增加；它是可触发的。因此，我们依赖蒙特卡罗积分，并从后验预测分布进行模拟。为了评估观测值x=yot+1时的预测密度，Darjus Hosszejni，Gregor Kastner 15称为预测可能性，我们应用计算p（x | yo[1:t]）≈MMXm=1p（x | yo[1:t]，κ（m）），（13），其中κ（m）表示SV模型估计过程中的第m个后验样本。对于其他应用，后验预测分布的分位数，即预测分位数，可能会引起兴趣。我们通过从yt+1模拟的随机变量来估计q%分位数~ p（yt+1 | yo[1:t]），我们通过重复两步形式=1，…，获得，M：步骤1。从SV后部p模拟κ（m）（κ| yo[1:t]），并执行步骤2。从p（yt+1 | yo[1:t]，κ（m））模拟y（m）t+1。最后，我们取样本向量（y（1）t+1，y（2）t+1，…）的q%分位数，y（M）t+1）>作为预测密度的近似分位数。我们在stochvol中实现了预测似然和预测分位数的估计。所有四个滚动窗口例程svsample\\u roll（）、svtsample\\u roll（）、svlsample\\u roll（）和svtlsample\\u roll（）都具有相同的编程接口。他们希望输入长度为L的输入数据y[1:L]作为他们的目标。为了估计每个时间窗口中的SV模型j=1，J在移动或扩展窗口方案中，子向量y[J：（t+J-1） ，或分别为y[1：（t+j-1） ，作为数据，用于预测nahead≥ 1次stepsahead。第一个时间窗口的宽度t可以通过L、J和nahead确定。

下面的示例演示如何调用滚动窗口采样例程instochvol。R> 设置。seed（5）R>res<-svsample\\u roll（CHF\\u logret，n\\u ahead=1，forecast\\u length=30，+refit\\u window=“moving”，calculate\\u quantile=c（0.01，0.05），+calculate\\u predictive\\u likelion=TRUE）参数n\\u ahead用于设置nahead，forecast\\u length用于设置J，refit\\u window期望“moving”或“expanding”分别将滚动窗口方案设置为moving或expanding。参数calculate\\u分位数需要一个介于0和1之间的数字向量；这些数字被解释为要预测的分位数。此外，如果Calculate\\u predictive\\u Likelion设置为TRUE，则该函数将估计PredictiveLikelion。最后，输出res是长度J的列表，即每个时间窗口一个元素。它包含各自的后分位数和预测似然结果以及κ的所有后验参数图。4.5. 指定第2.2节中讨论的先验超参数，需要在估计过程开始之前指定先验分布。关于常见的模型参数u、Д和σ，所有的svsample（）、svtsample（）、svlsample（）和svtlsample（）分别通过它们的输入参数priormu、priorphi和priorsigma值（bu、pBu）、（aД、bД）和bσ进行期望。此外，所有抽样函数都接受参数priorbeta，通过提供（bβ，sβ），为回归系数设置16个建模随机波动率，其中bβ和sβ分别为常用平均值和常用标准差。对于一般的多变量非线性分布，存在specify\\u priors（）接口，我们将在本节后面详细介绍。通过将λν作为参数priornu传递，可以在svtsample（）和svtlsample（）中影响ν的优先级。

最后，svlsample（）和svtlsample（）通过输入参数priorho获取数字序列（aρ，bρ）。下面的代码段显示了先前超参数的所有默认值。R> svsample（CHF\\u logret，priormu=c（0，100），priorphi=c（5，1.5），+priorsigma=1，priorbeta=c（0，10000））R>svtsample（CHF\\u logret，priormu=c（0，100），priorphi=c（5，1.5），+priorsigma=1，priorbeta=c（0，10000），priornu=0.1）R>svlsample（CHF\\u logret，priormu=c（0，100），priorphi=c（5，1.5，+priorsigma=1，priorbeta=c（0，10000），priorho=c（4，4））R>svtlsample（CHF\\u logret，priormu=c（0，100），priorphi=c（5，1.5），+priorsigma=1，priorbeta=c（0，10000），priornu=0.1，+priorho=c（4，4））作为上述简明界面的替代，可以通过specify\\u priors（）函数创建的对象来指定更广泛的先验分布。该函数为每个模型参数都有一个输入参数：mu、phi、sigma2、nu、rho、beta，另外还有一个用于hcalled latent0\\u方差的方差。有一个在stochvol中创建分布的附带函数列表：sv\\u beta（）具有参数shape1和shape2，并且它接受phi和rho；sv\\u constant（）有参数值，它可以接受为u、phi、sigma2、nu、rho和latent0\\u方差；sv\\u normal（）有参数mean和sd，它可以用于mu和phi；sv\\u multinormal（）有参数mean和eithersd，dim或precision，beta接受该参数；sv\\u index（）具有argumentrate，可用于nu；sv\\u gamma（）有参数shape和rate，它被sigma2接受；sv\\u reverse\\u gamma（）具有参数shape和scale，可用于IGMA2；sv\\u infinity（）没有参数，它被nu接受，因此将学生的t分布转化为正态分布。此外，latent0\\u方差接受字符值“stational”。

所有四种采样方法都通过输入参数priorspec接受先前的SpecifionObject。specify\\u prior的所有输入参数都是可选的，其默认值和使用方式如下所示。R> ps<-指定\\u先验值（mu=sv\\u normal（mean=0，sd=100），+phi=sv\\u beta（shape1=5，shape2=1.5），rho=sv\\u常量（0），+sigma2=sv\\u gamma（shape=0.5，rate=0.5），nu=sv\\u无穷大（），+beta=sv\\u多正态（mean=0，sd=10000，dim=1），+latent0\\u方差=“平稳”）R>svsample（CHF\\u logret，priorspec=ps）4.6。建立马尔可夫链使用MCMC抽样方案进行贝叶斯推理时，后验分布的抽取次数、所谓老化阶段的长度、Darjus Hosszejni、Gregor Kastner 17的初始值马尔可夫链以及存储结果的各种策略都是人们普遍感兴趣的。输入参数draws和burnin解决了前两点。从后验分布中获得burnin+绘图的样本量，其中第一个burnin绘图数被丢弃。默认情况下，对于无杠杆的SVmodels，在烧入1000后绘制10000个元素，对于有杠杆的SV models，在烧入2000后绘制20000个元素，根据我们的经验，这对于大多数应用程序来说已经足够了。至于初始值，startpara和StartAtent提供了一种设置它们的方法。参数startpara应该是一个命名列表，将参数名称映射到startingvalues，starttant必须是长度为m的序列，其中包含h的起始值。默认值设置为ν、σ、ν和ρ的先验平均值，这些值对马尔可夫链只有最小影响。β的默认值是普通最小二乘估计值（X>X）-1X>y，其中X表示回归设计矩阵，y表示观测向量。

设置β后，根据贝叶斯线性回归log（yt）=u+ξt，u初始化对数方差u的水平~ N（bu，bu），（14），其中ξt~ N个(-1.27, 4.934). 方程14由方程1的第一行得出，取其固定期望值u，然后取x 7→ 两侧的对数（x）。作为对数（εt）的拉普拉斯近似，获得了齐次误差项ξ（Harvey and Shephard 1996）。最后，默认情况下，向量startatent的所有值都设置为初始值u。通常是并行启动独立的马尔可夫链，stochvol在其所有采样过程中都提供了这方面的便利。参数n\\u chains应为正整数，它设置独立链的数量。此外，参数parallel、n\\u cpu和cl可用于控制stochvol使用的并行性。为了覆盖默认的顺序执行策略，将parallel设置为“snow”，以使用所谓的“snow”集群，或设置为“multicore”，以使用“multicore”类型的计算（R Core Team 2020）。接下来，参数n\\u CPU应设置为要使用的并行处理单元的物理数量。最后，在应用“雪”的情况下，采样例程通过参数cl选择性地接受一个预先运行的“雪”簇。如前所述，stochvol中潜态h的采样算法依赖于Omori et al.（2007）和Kastner and FrühwirthSchnatter（2014）中的高斯混合近似。近似值往往非常好，因此默认设置为不更正模型错误。但是，可以通过以下svsample（）的专家参数在所有采样例程中启用此更正。最后，stochvol提供了三种在取样器执行期间和之后节约存储的方法。

将整型参数thinpara设置为ι会告诉采样器仅存储参数向量的每个ι图，为thinlatent提供一个值也会为h提供一个值。最后，通过设置keeptime=“last”，可以不存储完整向量h，而只存储其最后一个值。默认行为是在磨合阶段后存储每个绘图。5、factorstochvol软件包18建模R-6.-2 2 4 6 UD-4.-2 0 1 2 CAD2009 2010 2011 2012时间-2 0 2 4 6 8瑞士法郎-3.-1 1 2 3CZK-0.10 0.00 0.05DKK2009 2010 2011 2012时间-2 0 1 2 3英镑图4：六欧元汇率的百分比日志回报。使用factorstochvol拟合多元因子SV模型最常见的工作流程包括以下步骤：（1）准备数据，（2）确定识别结构，（3）指定先前的超参数，（4）运行采样器，（5）调查输出并将结果可视化，以及（6）预测（如果需要）。以下各节将详细描述这些步骤。5.1. 准备数据factorstochvol中的主力是采样函数fsvsample（）。它期望数据以矩阵Y=（Y，…，yn）>的形式出现，其中包含n行和m列。为了举例说明，我们使用了stochvol中的汇率数据集，其中包含3140个每日观察值，从2000年3月3日到2012年4月4日，23种货币对的汇率。为了保持分析简单且计算时间适中，我们仅按字母顺序（澳元、加元、瑞士法郎、捷克克朗、丹麦克朗、英镑）对前六个系列的最后1001天进行建模，以供进一步分析。我们不使用名义汇率，而是计算对数回报。这给我们留下了一个大小为n=1000和m=6的数据集。

数据是使用下面的代码片段准备的，并在图4中使用zoo包可视化（Zeileis和Grothendieck 2005）。R> 库（“factorstochvol”）R>库（“zoo”）R>数据（“exrates”，package=“stochvol”）R>m<-6R>n<-1000R>y<-100*logret（tail（exrates[，seq\\u len（m）]，n+1））R>y<-zoo（y，order.by=tail（exrates$date，n））R>绘图（y，main=，xlab=“Time”）5.2。决定识别结构Darjus Hosszejni，Gregor Kastner 19因子模型中的可能性对某些因子转换是不变的，例如因子及其载荷的重新排序或符号切换。除此之外，它通常是多模式的。因此，识别因子负荷远非小事。在因子SV模型中解决此问题的最常用方法是施加较低的对角线因子载荷矩阵，其中对角线以上的所有元素均设置为零（例如Aguilar和West 2000；Chib等人2006；Han 2006；Zhou、Nakajima和West 2014）。要在factorstochvol中使用此约束，可以将参数restrict=“upper”传递给主采样函数fsvsample（）。显然，这种做法带来了顺序依赖性，例如，第一个变量不允许加载到除第一个因素以外的任何其他因素上。在助手函数preorder（）中实现了一种自动排序数据的特别方法。对数据进行最大似然因子模型拟合后（使用stats包中的factanal（）和默认的varimax旋转），序列排序如下：因子1上负荷最高的变量放在第一位，因子2秒上具有最高负荷的变量（除非该变量已经放在第一位，在这种情况下，取第二高负荷的变量），等等。

对于手头的数据集，这意味着双因素模型的排序如下。R> 前序（y，因子=2）##[1]2 3 1 4 5 6根据此算法，第二个序列应放在第一位，第三个序列应放在第二位。此后，字母顺序保持不变。为了在不重新排序数据的情况下实现这一效果，可以通过restrict将大小为m×r的逻辑矩阵传递给fsvsample（），其中条目TRUE表示该元素被限制为零；FALSE表示根据数据进行估计。与preorder（）类似，函数findrestrict（）尝试自动化此过程。再次，使用静态因素分析的最大可能性估计；然而，findrestrict（）使用的算法与上述算法稍有不同：因子2、3、…、上绝对载荷最低的变量，确定r（相对于因子1）领先于第一个因子，即因子3、4、…、上绝对载荷最低的变量，r（相对于系数1和2）排在第二位，等等。下面是手头数据集的结果。R> findrestrict（y，factors=2）##[，1][，2][1，]FALSE FALSE##[2，]FALSE TRUE##[3，]FALSE FALSE##[4，]FALSE FALSE##[5，]FALSE##[6，]FALSE如果使用参数restrict=“auto”调用fsvsample（），它会自动调用findrestrict具有适当数量的因子的strict（）。使用restrict=“none”（默认）会使采样器不会对载荷矩阵施加任何约束；因此，在重新生成后验平仓的过程中，20建模随机波动率可能不稳定或受到多种局部模式的影响。然而，如果分析的主要关注点不是对因子载荷本身的推断，那么不确定因子载荷可能是首选。

当寻求协方差矩阵的推断时，尤其是这种情况，因为这只取决于通过等式9的旋转不变变换得到的∧。为了更详细地讨论这些问题，我们请读者参考Sentana和Fiorentini（2001），他们讨论了异性恋的自动识别。Kastner等人（2017）对因子SV模型在不同识别模式下的对数预测得分进行了比较；静态因素模型的相关问题，另请参见Frühwirth Schnatterand Lopes（2018）。为了继续当前示例，我们选择在因子载荷矩阵之前使用行正常伽马收缩时，不对因子载荷矩阵施加任何先验限制（参见Kastner2019）。5.3. 指定先验超参数从上述关于因子数量和识别方案的明显先验选择中，可以在factorstochvol中找到许多超参数选择。关于对数方差过程，该接口类似于svsample（）。在下面，i=1，m和j=1，r分别对特质和因子对数方差过程进行索引。这对常见的先验超参数（bβ，bβ）可以作为长度为2的序列传递给先验超参数。