如何用r实现基于mcmc方法的garch模型参数估计

以AR(3)-GARCH(2,1)模型为例：

首先在主窗口输入

LS RR RR（-1） （-2） （-3）

得出

Variable Coefficient Std. Error t-Statistic Prob.

RR(-1) 0.007606 0.059014 0.128883 0.8975

RR(-2) 0.058005 0.058549 0.990707 0.3227

RR(-3) 0.121110 0.058985 2.053245 0.0410

然后在点estimate 在下拉选项中选择ARCH

在命令窗口中再次输入

LS RR RR（-1） （-2） （-3）

并在ARCH出填入2，GARCH处为1，得出结果

Variance backcast: ON

GARCH = C(4) + C(5)*RESID(-1)^2 + C(6)*RESID(-2)^2 + C(7)

*GARCH(-1)

Coefficient Std. Error z-Statistic Prob.

RR(-1) 0.013392 0.056863 0.235514 0.8138

RR(-2) 0.120481 0.062146 1.938671 0.0525

RR(-3) 0.095921 0.056070 1.710743 0.0871

Variance Equation

C 0.000127 3.59E-05 3.553327 0.0004

RESID(-1)^2 -0.043907 0.029463 -1.490253 0.1362

RESID(-2)^2 0.248625 0.078855 3.152960 0.0016

GARCH(-1) 0.079769 0.211942 0.376372 0.7066

R-squared 0.003674 Mean dependent var 0.001397

Adjusted R-squared -0.017908 S.D. dependent var 0.013305

S.E. of regression 0.013423 Akaike info criterion -5.819411

Sum squared resid 0.049910 Schwarz criterion -5.729472

Log likelihood 833.3564 Durbin-Watson stat 1.974819

RR是上证综合指数的周收益，用此AR（3）-GARCH（2，1）是用残差来检验超额收益的。

写在最前面：

首先需要说一下，本文的bootstrap和jackknife都算是蒙特卡罗方法(Monte Carlo method)的一种。应用广泛的的MCMC链(马尔可夫链蒙特卡洛方法Markov chain Monte Carlo)也是蒙特卡罗与马尔可夫链的结合。简单来说，蒙特卡罗方法就是从已知样本的分布中随机抽取新的样本集进行评估，然后放回，再次抽取的方法。根据具体方法的不同，抽取样本集的手段也不同。

bootstrap抽样方法将观测到的样本视为一个有限的总体，是唯一的信息来源，从中有放回的随机抽样来评估总体特征，以及对抽样总体进行推断统计。bootstrap 也分参数bootstrap和非参数bootstrap，前者的分布已完全知道。但在生信领域一般没有这种情况。所以下面讨论的是非参数bootstrap。

直接上例子：

假设现在有bootstrap包中的law数据集如下，

现在我们要计算LSAT成绩（美国法学入学考试）和GPA之间的相关系数。但因为样本量太少了，所以我们使用bootstrap重复抽样评估其标准误。

200次循环抽样后，计算得se.R标准误为0.1474629

得到如下的图：

1e6次循环抽样后，计算得se.R标准误为0.1333802

得到如下的图：

如果用bootstrap包的bootstrap函数会快一些：

bootstrap函数的用法： bootstrap(抽取样本范围,重复次数,进行bootstrap的函数,bootstrap的数据集)

偏差定义为bootstrap结果（多个数值）与原数据统计结果（单个数值）的均值：

得到bias大约为0.001817608，比较小

换一个包，boot包

这里用了三种方法计算置信区间：basic、正态和百分数。样本相关系数分布接近正态，则正态置信区间接近百分数区间。此外还有“Better Bootstrap Confivendence Interval” 更好的bootstrap置信区间，称为BCa区间，使用偏差和偏度对百分数置信区间进行矫正。设置type="bca"即可。

简单的说，bootstrap是从原有真实样本中有放回地抽取n个。jacknife就是每次都抽取n-1个样本，也就是每次只剔除一个原样本。

同样地，如果以bootstrap包中的law数据进行演示：

Jackknife计算的bias为-0.006473623。 这里jackknife的偏差公式相比于bootstrap有一个(n-1)系数，推导就不写了。

标准误se为0.1425186，与bootstrap得出的比较接近。

当统计量不太平滑的时候，Jacknife有很大误差。比如说对中位数进行统计，其变化很大。在进行Jacknife之后最好再跑一次bootstrap，看看是否相差很大。

居然还能这么嵌套着玩，针对每次bootstrap形成的数列向量计算jackknife的标准差，这样可以看出bootstrap若干次取样之间的差异。

算出来分别为0.1344824和0.08545141。后者较小，表面bootstrap取样之间的variance较小。

简单来说就是一种数据分割检验的方法，将数据分割为K份，称为"K-fold"交叉检验，每次第i个子集作为测试集来评估模型，其余的用来构建模型。Admixture使用的就是这个原理。Jackknife也属于Cross Validation的应用之一。

现在我创建一个这样的alignment：

这棵树长这样，符合遗传距离：

进行bootstrap：

phylogeny的bootstrap是对每一个节点都进行bootstrap取样并建树，比如说在9号节点，查看其bootstrap子集建的树符合系统发育关系((human2,human4,human3)(human8,human1,human6,human7,human5))的百分比（不管内部怎么样，先看这个节点）。发现Node1支持率是100（1000次都符合）。而后移到下一个节点，并且只看节点内部的分支支持率是多少。

其实原理都比较简单，计算bootstrap也会有专门的软件。

参考资料：

1）中科大张伟平教授课件

2） https://ecomorph.wordpress.com/2014/10/09/phylogenetic-trees-in-r-4/

种间竞争是指具有相似要求的物种，为了争夺空间和资源，而产生的一种直接或间接抑制对方的现象。

Lotka-Volterrra的种间竞争模型：

K是环境容纳量

N则是物种的种群数量

两者稳定共存的基本条件就是dN1/dt=0和dN2/dt=0

若dN1/dt=0，则K1-N1-αN2=0——①

若dN2/dt=0，则K2-N2-βN1=0——②

当环境全被N1占领则由②式得出N1=K2/β

当环境全被N2占领则由①式得出N2=K1/α

四种情况如图所示

---