R語言與抽樣技術(shù)學(xué)習(xí)筆記（Jackknife）

R語言與抽樣技術(shù)學(xué)習(xí)筆記（Randomize，Jackknife，bootstrap）

Jackknife算法

Jackknife的想法在我很早的一篇博客《R語言與點(diǎn)估計(jì)學(xué)習(xí)筆記（刀切法與最小二乘估計(jì)）》就提到過。其思想以一言蔽之就是：將樣本視為總體，在“總體”中不放回地抽取一些“樣本”來做統(tǒng)計(jì)分析。今天我們來討論Jackknife應(yīng)該怎么做以及為什么要這么做。

Jackknife的算法描述

Jackknife方法利用系統(tǒng)的劃分?jǐn)?shù)據(jù)集的辦法來推測總體樣本估計(jì)量的一些性質(zhì)。Quenouille建議用它來估計(jì)偏差，隨后John Tukey證實(shí)它用來估計(jì)估計(jì)量的方差也是極好的。
　　假設(shè)我們有隨機(jī)樣本X1,?,Xn,并從中計(jì)算統(tǒng)計(jì)量T去估計(jì)總體的參數(shù)μ。在Jackknife方法下，我們將給定數(shù)據(jù)集劃分為r組，每組數(shù)據(jù)量為k。
　　現(xiàn)在，我們移除樣本中第jth組數(shù)據(jù)，并用剩下的數(shù)據(jù)來估計(jì)參數(shù)μ,并將估計(jì)量記為T(?j)。T(?j)的均值Tˉˉˉ(?)可以用來估計(jì)參數(shù)μ，T(?j)也可以用來獲取估計(jì)量T更多的信息，但是必須要指出的是Jackknife不會提供比總體樣本更多的信息（任何抽樣技術(shù)都是不會的?。。。?；當(dāng)T是無偏估計(jì)時(shí)，T(?j)也是，T是有偏的估計(jì)量時(shí)，T(?j)也是，但是會有一點(diǎn)點(diǎn)的不同。
　　我們有時(shí)也對T與T(?j)的不同進(jìn)行加權(quán)處理，得到新的統(tǒng)計(jì)量:

這就是許多文獻(xiàn)中提到的Jackknife“偽值”，并將這些偽值的均值稱為“Jackknifed”T，記為T(J)。顯然　　通常情況下，我們?nèi)=1,r=n。在某些特定條件下，它是最優(yōu)的。
　　現(xiàn)在我們來考慮一下Jackknife對估計(jì)量偏差與方差的估計(jì)。我們這里采用偽值來考慮問題，認(rèn)為偽值的偏差就是估計(jì)量的偏差，偽值的方差就是估計(jì)量的方差。那么我們有

我們不妨取r=n,k=1，T對μ的估計(jì)的偏差肯定是，忽略高階量，那么偏差的近似估計(jì)量為：

雖然這些偽值不是獨(dú)立的，但是，我們?nèi)圆环良僭O(shè)他們是獨(dú)立的，因?yàn)檫@樣我們利用偽值估計(jì)估計(jì)量的方差變得十分簡單：

我們以正態(tài)總體為例來考慮這個(gè)問題。

    data.sim <- rnorm(30, 5, 3)  
    (mu.hat <- mean(data.sim))  

## [1] 4.339

[plain] view plain copy
 print?

    n <- length(data.sim)  
    mu.jack <- NULL  
    for (i in 1:n) {  
        mu.jack[i] <- mean(data.sim[-i])  
    }  
    (jack.estimate <- mean(mu.jack))  

## [1] 4.339

可以看到，矩估計(jì)的估計(jì)值為4.3393，Jackknife估計(jì)為4.3393,較真值5而言，估計(jì)效果還是可以的。我們還可以計(jì)算偏差：

    (bias <- (n - 1) * (mean(mu.jack) - mu.hat))  

## [1] 0

可以看到，這是一個(gè)無偏估計(jì)，不是嗎？
　　當(dāng)然，要估計(jì)我們采取的估計(jì)量的方差也是可以的：

var <- mean((mu.jack - mean(mu.jack))^2) * (n - 1)
print(var)

## [1] 0.4612

可以看到估計(jì)量的方差為0.4612與正態(tài)總體均值的標(biāo)準(zhǔn)差的理論值0.3相差不大。

Jackknife與偏差糾正

我們這里將要說明Jackknife最大的作用——糾正偏差。
　　我們都知道，我們對參數(shù)估計(jì)最常用的辦法就是矩估計(jì)與極大似然估計(jì)。然而這兩種估計(jì)不一定是無偏的。例如[0,θ]上的均勻分布U(θ)的參數(shù)θ的極大似然估計(jì)就不是無偏的，在比如正態(tài)分布N(μ,σ)中標(biāo)準(zhǔn)差的似然估計(jì)也不是無偏的。這些很容易證明，當(dāng)然，你也可以查閱王兆軍的《數(shù)理統(tǒng)計(jì)講義》，上面給出了詳細(xì)的說明。　　當(dāng)然這些偏差我們都可以通過對估計(jì)量做一些變換得到無偏估計(jì)量，但是這些常數(shù)還是不太容易通過普適的辦法得到，而大偏差估計(jì)往往被視為估計(jì)的不足，得到一種在不增加方差的基礎(chǔ)上，糾正偏差的辦法是很好的，這就是我們現(xiàn)在討論的Jackknife。
　　設(shè)T(X)是基于樣本X的參數(shù)g(θ)的估計(jì)量，且滿足，Jackknife偏差修正估計(jì)量為，Jackknife統(tǒng)計(jì)量具有如下性質(zhì)：

證明十分的簡單，注意到J(T)的表達(dá)式，利用定義即可證明。（上一小節(jié)我們也間接地證明過了這一事實(shí)）
　　我們最后來看一個(gè)例子：我們知道均勻分布U(θ)的參數(shù)θ的極大似然估計(jì)是漸進(jìn)無偏而非真正無偏的，我們假設(shè)樣本為x1,?,xn，θ的極大似然估計(jì)是x(n)（即樣本中的最大值），我們知道，所以似然估計(jì)的偏差是,我們使用Jackknife看看能不能縮減方差：

data.sim <- runif(100, 0, 7)  
    theta.hat <- max(data.sim)  
    theta.jack <- NULL  
    for (i in 1:length(data.sim)) {  
        theta.jack[i] <- max(data.sim[-i])  
    }  
    n <- length(data.sim)  
    theta.jackestimate <- n * theta.hat - (n - 1) * mean(theta.jack)  
    cat("original bias is ", 7 - theta.hat, "after jackknife the bias is ", 7 -   
        theta.jackestimate) 

## original bias is 0.002576 after jackknife the bias is -0.1084

我們從上面的運(yùn)行結(jié)果可以清楚的看到Jackknife之后，偏差確實(shí)減少了。這也是符合我們的想法的，因?yàn)楫?dāng)k變大時(shí)，Jackknife估計(jì)更接近均值的兩倍，即它更接近矩估計(jì)，而矩估計(jì)是無偏的。當(dāng)然一味追求無偏是不對的，對于一個(gè)有偏估計(jì)，其偏差趨于0時(shí)，可能導(dǎo)致方差變得很大。

Jackknife失效

若估計(jì)量不夠平滑，每次刪掉一個(gè)數(shù)據(jù)的Jackknife就會失效，估計(jì)也不再具備相合性。
　　例如：利用Jackknife方法估計(jì)從1到100中隨機(jī)抽出的10個(gè)數(shù)的中位數(shù)的標(biāo)準(zhǔn)差。

## M= 46 46 40 40 40 46 40 46 46 40

## [1] 9

我們可以看到這個(gè)估計(jì)量的標(biāo)準(zhǔn)差很不對的，因?yàn)槊看稳笔?個(gè)值，中位數(shù)也就至多會有兩個(gè)取值，無論你的樣本容量取值多大，這都是不對的。Enfro（1993）提出了delete K Jackknife算法解決這個(gè)問題。但是這個(gè)要求還是很高的，它要求n,k都比較大。那么好的估計(jì)方法應(yīng)該是什么呢？這就要用到我們接下來將要介紹的bootstrap方法。