R語言蒙特卡洛計(jì)算和快速傅立葉變換計(jì)算矩生成函數(shù)
原文鏈接：http://tecdat.cn/?p=13734
對(duì)精算科學(xué)來說，當(dāng)我們處理獨(dú)立隨機(jī)變量的總和時(shí)，特征函數(shù)很有趣，因?yàn)榭偤偷奶卣骱瘮?shù)是特征函數(shù)的乘積。?
介紹
在概率論中，讓?
?對(duì)于?
?和?
?對(duì)于?
?是一些隨機(jī)變量的累積分布函數(shù)?
，即?
。什么是矩生成函數(shù)?
，即?
?？
如何編寫?
?？
在概率教科書中，標(biāo)準(zhǔn)答案是
如果?
?是離散的
如果?
?（絕對(duì)）連續(xù)，
?
?是的密度?
。這里，?
?顯然不是離散變量。但是是連續(xù)的。需要繪制該分布函數(shù)以查看，?
， 對(duì)所有?
?
我們有一個(gè)不連續(xù)的0。因此，我們在這里必須謹(jǐn)慎一些：?
?既不是連續(xù)的也不是離散的。讓我們使用公式，
如果也可以寫?
，
這只是說總體平均值是每個(gè)子組平均值的重心。
?然后讓?
?而?
?
）。
讓我們考慮三個(gè)不同的組成部分。
?
?
（因?yàn)樗且粋€(gè)實(shí)值常量），在這里?
。
所以最后，我們計(jì)算?
。觀察一下?
?給定?
?是具有密度的（絕對(duì)）連續(xù)隨機(jī)變量。觀察所有?
，
和?
，即?
?給定?
?是指數(shù)分布。
因此，?
?是指數(shù)變量和Dirac質(zhì)量之間的混合?
。這實(shí)際上是問題的棘手部分，因?yàn)楫?dāng)我們看到上面的公式時(shí)，它并不明顯。
從現(xiàn)在開始，這是高中階段的計(jì)算，
如果?
?。如果把所有的放在一起
蒙特卡洛計(jì)算
可以使用蒙特卡洛模擬來計(jì)算該函數(shù)，
> F=function(x) ifelse(x<0,0,1-exp(-x)/3)> Finv=function(u) uniroot(function(x) F(x)-u,c(-1e-9,1e4))$root
或（以避免不連續(xù)的問題）
> Finv=function(u) ifelse(3*u>1,0,uniroot(function(x)+ F(x)-u,c(-1e-9,1e4))$root))
在這里，逆很容易獲得，因此我們可以使用
然后，我們使用
> plot(u,v,type="b",col='blue')> lines(u,Mtheo(u),col="red")
?
蒙特卡洛模擬的問題在于，僅當(dāng)它們有效時(shí)才應(yīng)使用它們。我可以計(jì)算
> M(3)[1] 5748134
有限總和始終可以通過數(shù)字計(jì)算。就算在這里?
?不存在。就像Cauhy樣本的平均值一樣，即使期望值不存在，我也總是可以計(jì)算出來
> mean(rcauchy(1000000))[1] 0.006069028
這些生成函數(shù)在存在時(shí)會(huì)很有趣。也許使用特征函數(shù)是一個(gè)更好的主意。
生成函數(shù)
首先，讓我們定義那些函數(shù)。
?
如果?
?足夠小。
現(xiàn)在，如果我們使用泰勒展開式
和
如果我們看一下該函數(shù)在0點(diǎn)的導(dǎo)數(shù)的值，那么
?可以為某些隨機(jī)矢量在更高維度上定義一個(gè)矩生成函數(shù)?
，
?
如果要導(dǎo)出給定分布的矩，則一些矩生成函數(shù)很有趣。另一個(gè)有趣的特征是，在某些情況下，此矩生成函數(shù)（在某些條件下）完全表征了隨機(jī)變量的分布。?
，
?對(duì)所有人?
， 然后?
。
?
快速傅立葉變換
回想一下歐拉公式，
因此，看到傅立葉變換就不會(huì)感到驚訝。從這個(gè)公式，我們可以寫
使用傅立葉分析中的一些結(jié)果，我們可以證明概率函數(shù)滿足
也可以寫成
如果在點(diǎn)處的分布是絕對(duì)連續(xù)的，則可以獲得類似的關(guān)系?
，
實(shí)際上，我們可以證明，
然后可以使用1951年獲得的吉爾-佩萊阿茲（Gil-Peleaz）的反演公式來獲得累積分布函數(shù)，
這意味著，在金融市場上工作的任何人都知道用于定價(jià)期權(quán)的公式（例如，參見??Carr＆Madan（1999）??）。好處是，可以使用任何數(shù)學(xué)或統(tǒng)計(jì)軟件來計(jì)算這些公式。
特征函數(shù)和精算科學(xué)
對(duì)精算科學(xué)來說，當(dāng)我們處理獨(dú)立隨機(jī)變量的總和時(shí)，特征函數(shù)很有趣，因?yàn)榭偤偷奶卣骱瘮?shù)是特征函數(shù)的乘積?？紤]計(jì)算Gamma隨機(jī)變量復(fù)合和的99.5％分位數(shù)的問題，即
?
?和?
。策略是分散損失金額，
然后，要計(jì)算的代碼?
， 我們用
?99.5％分位數(shù)
> sum(cumsum(f)<.995)
考慮以下?lián)p失金額
> print(X[1:5])[1] 75.51818 118.16428 14.57067 13.97953 43.60686
讓我們擬合一個(gè)伽瑪分布。我們可以用
shape ? ? ? ? rate1.309020256 ? 0.013090411(0.117430137) (0.001419982)
?
> alpha[1] 1.308995> beta[1] 0.01309016
無論如何，我們都有個(gè)人損失的Gamma分布參數(shù)。并假設(shè)泊松計(jì)數(shù)變量的均值為
> lambda <- 100
同樣，可以使用蒙特卡洛模擬。我們可以使用以下通用代碼：首先，我們需要函數(shù)來生成兩種感興趣的變量，
如果我們生成一百萬個(gè)變量，我們可以得到分位數(shù)的估算，
> set.seed(1)> quantile(rcpd4(1e6),.995)99.5%13651.64
另一個(gè)想法是記住Gamma分布的比例：獨(dú)立Gamma分布的總和仍然是Gamma（在參數(shù)上有附加假設(shè)，但在此我們考慮相同的Gamma分布）。因此，可以計(jì)算復(fù)合和的累積分布函數(shù)，
如果我們求解那個(gè)函數(shù)，我們得到分位數(shù)
> uniroot()$root[1] 13654.43
這與我們的蒙特卡洛計(jì)算一致?，F(xiàn)在，我們也可以在此處使用快速傅立葉變換，
> sum(cumsum(f)<.995)[1] 13654
讓我們比較獲得這三個(gè)輸出的計(jì)算時(shí)間
> system.timeuser ? ? ?system ? ? elapsed2.453 ? ? ? 0.106 ? ? ? 2.611> system.timeuser ? ? ?system ? ? elapsed0.041 ? ? ? 0.012 ? ? ? 0.361> system.timeuser ? ? ?system ? ? elapsed0.527 ? ? ? 0.020 ? ? ? 0.560
?
