2012-03-22 78 views
2

我想計算分佈的第n個時刻。我想利用圖書館「瞬間」在R.我已經以這種方式測試all.moments的all.moments功能:all.moments函數奇怪的結果

library(moments) 
r<-rnorm(10000) 
rr<-all.moments(r,order.max=4) 
rr 

[1] 1.000000000 0.002403360 0.962201478 -0.022694670 2.852696159 

這在我看來,這不是真的,因爲我知道,3-第4和第4時刻必須爲0。 我的錯誤在哪裏?

回答

6

第三個時刻是偏度。你是正確的:對於正態分佈,這是零。既然你只是從正態分佈抽樣,你的結果將近似爲零。

第四階矩是峯度。對於正態分佈,這是3σ^ 4。在這種情況下,σ是1,所以你的結果應該是3,它是。


要提高估算的準確性,請改善樣本量。對於1e7觀測的樣本:

> library(moments) 
> r <- rnorm(1e7) 
> all.moments(r,order.max=4) 
[1] 1.0000000000 0.0004028138 0.9995373115 0.0007276404 2.9976881271 
1

因爲這隻在預期中是正確的,而不是精確的,並且因爲高階矩有很大的差異?

(見@ Andrie的回答爲好,爲什麼第四時刻(以下V5)甚至還沒有接近零。)

> library(moments) 
> R <- t(replicate(50,all.moments(rnorm(1e4),order.max=4))) 
> summary(R) 
     V1   V2     V3    V4    
Min. :1 Min. :-0.024921 Min. :0.9714 Min. :-0.0987174 
1st Qu.:1 1st Qu.:-0.009527 1st Qu.:0.9911 1st Qu.:-0.0341950 
Median :1 Median : 0.001021 Median :0.9994 Median : 0.0067138 
Mean :1 Mean :-0.001047 Mean :1.0006 Mean :-0.0002613 
3rd Qu.:1 3rd Qu.: 0.004711 3rd Qu.:1.0147 3rd Qu.: 0.0299731 
Max. :1 Max. : 0.023356 Max. :1.0398 Max. : 0.1283456 
     V5  
Min. :2.775 
1st Qu.:2.921 
Median :3.005 
Mean :3.007 
3rd Qu.:3.092 
Max. :3.325 
0

我得到了同樣的問題,在蟒蛇的工作,我認爲計算機上浮點運算的精度有限,擁有大量的權力,也沒有幫助。我將剪切並粘貼我的Python代碼和我得到的結果。此代碼嘗試計算標準法線的前20個時刻。總之,我認爲計算分佈的高階矩是不容易的,「高階」在這裏意味着大於10左右。在另一個實驗中,我試圖通過繪製越來越多的樣本來減少我在第18次得到的差異,但是這並不實際,或者是給我的「普通」計算機。

N = 1000000 
w = np.random.normal(size=N).astype("float128") 

for i in range(20): 
    print i, mean(w**i) # simply computing the mean of the data to the power of i 

爲您提供:

0 1.0 
1 0.000342014729693 
2 1.00124397377 
3 0.000140133725668 
4 3.00334304532 
5 0.00506625342217 
6 15.0227401941 
7 0.0238395446636 
8 105.310071549 
9 -0.803915389936 
10 948.126995798 
11 -34.8374820713 
12 10370.6527554 
13 -1013.23231638 
14 132288.117911 
15 -26403.9090218 
16 1905267.02257 
17 -658590.680295 
18 30190439.4783 
19 -16101299.7354 

但正確的時刻是:1,0,1,0,3,0,15,0,105,0,945,0,10395,0 ,135135,0,2027025,0,34459425,0,654729075.

+0

它不是數值不精確的問題:它是高階矩的方差爆炸。 – 2012-11-23 00:03:47

+0

但是在這裏,我給了自己大量的N個標準正態分佈的隨機數,例如計算它們的18次方,然後計算其平均值。在這個過程中,如果我增加N,我不應該越來越接近第18個時刻的真實值嗎? – Frank 2012-11-23 00:19:25

+0

是的,但你必須增加它*很多*。 – 2012-11-23 03:02:23