その１はこちら（極値理論の簡単な解説、Block maximaの実装を行っています）
https://qiita.com/hrkz_szk/items/43debffda9697d9dd7a9

r-block maximaモデル

Block maximaでは最大値のみを抽出しましたが、問題となるのはサンプル数。特にファイナンス系のデータの場合、数十年単位で最大値を確保できないことが多いです。そこで、r-block maximaモデルとThresholdモデルがその対処法になります。r-block maximaは、上からr個の値を抽出します。例えば、2-block maximaモデルでは、最大値とその次に大きい値を抽出します。Thresholdモデルは次回紹介したいと思いますが、これらを統合したのが点過程というモデルです。

データ

ここでは、BMWの株価の損失について、2-block maximaモデルで推定していきたいと思います。具体的に、1年ごとの最大値とその次に大きい値を抽出し、ismevのrlarg.fitを使ってGEV分布への推定を行います。
データに関しては、evirというパッケージにあるdata(bmw)を使用します。こちらは1973年1月2日から1996年7月23日におけるBMW社の株価のログリータンがvectorで入っています。

library(evir)
data("bmw")
head(bmw)

[1]  0.047704097  0.007127223  0.008883307 -0.012440569 -0.003569961  0.000000000

plot(bmw, t="l")

データの抽出

ログリータンをプラスの損失に変換し、さらにパーセント表示に変換。その後、一年間の最大値とその次に大きい値を抽出。

BmwLoss <- -1.0 * bmw * 100 
Years <- format(attr(BmwLoss, "time"), "%Y")
attr(BmwLoss, "years") <- Years
Yearu <- unique(Years)
idx <- 1:length(Yearu)
r <- 2
BmwOrder <- t(sapply(idx, function(x)
              head(sort(BmwLoss[attr(BmwLoss, "years") == Yearu[x]], decreasing = TRUE), r)))
rownames(BmwOrder) <- Yearu
colnames(BmwOrder) <- paste("r", 1:r, sep = "")

サンプルの分布

plot(Yearu, BmwOrder[, 1], col = "black", ylim = range(BmwOrder),
     ylab = "Losses BMW (percentages)", xlab = "",
     pch = 20, bg = "black")
points(Yearu, BmwOrder[, 2], col = "blue", pch = 20, bg = "blue")

推定

library(ismev)
BmwOrderFit <- rlarg.fit(BmwOrder)

$conv
[1] 0

$nllh
[1] 77.76214

$mle
[1] 4.4798274 1.7524702 0.3034724

$se
[1] 0.3411369 0.2891531 0.1587266

convのゼロは、うまく推定されたことを示しています。nllhはマイナスのログ尤度の最小値で、mleが尤度関数を用いて推定されたパラメータ、seはパラメータの標準誤差になります。$\xi$が0.3034724とプラスなので、今回の分布はFrechet分布に従うことが分かります。このため、この分布はヘビーテイルであり、損失が有限ではないということになります。

パラメータの検証

信頼区間

私の知る限りismevにはパラメータの信頼区間を自動的に計算してくれるコードを知らないので、ここでは原始的な方法で95%の信頼区間を求めてみます。

lower <- BmwOrderFit$mle - 1.96*BmwOrderFit$se
upper <- BmwOrderFit$mle + 1.96*BmwOrderFit$se
data.frame(lower=lower, upper=upper, row.names = c("mu", "sigma", "xi"))

	lower	upper
mu	3.811199061	5.1484558
sigma	1.185730115	2.3192103
xi	-0.007631649	0.6145765

muとsigmaはゼロをまたいでいませんが、最も大切なxiがゼロの可能性が。しかし、ゼロをまたいでいると行ってもとても小さな幅なので、ここでは深追いはせずに…。

モデルの検証

rlarg.diag(BmwOrderFit)

Block maximaモデルの事例と同じように、Probability PlotとQuantilte Plotとともに右上の方がラインに沿っておらず、推定したモデルでは説明しきれていないことを示しています。しかし、再現レベルはサンプルデータが95%の内側に入っているので、この点は良いのかなと思います。

続編

その３はこちら。Thresholdモデルを使った分析を行っています。
https://qiita.com/hrkz_szk/items/8b01552e0e09583d4c35

参考文献

An Introduction to Statistical Modeling of Extreme Values (Springer Series in Statistics)

An Introduction to Statistical Modeling of Extreme Values (Springer Series in Statistics) [ハードカバー]

Stuart Coles

Springer

2001-12-15

Financial Risk Modelling and Portfolio Optimization with R

Financial Risk Modelling and Portfolio Optimization with R

Bernhard Pfaff

Wiley

2016-10-03