系列文章目录
统计计算一|非线性方程的求解
统计计算二|EM算法(Expectation-Maximization Algorithm,期望最大化算法)
统计计算三|Cases for EM
统计计算四|蒙特卡罗方法(Monte Carlo Method)
统计计算五|MCMC( Markov Chain Monte Carlo)
文章目录
- 系列文章目录
- 一、Bootstrap(自助法)
- (一)基本思想
- (二)偏差的自助估计
- (三)估计量标准差的 Bootstrap 估计
- 二、基于 Jackknife 法的估计
- (一)基本思想
- (二)估计量偏差的 Jackknife 估计
- (三)估计量标准差的 Jackknife 估计
- 三、Permutation Test 置换检验
一、Bootstrap(自助法)
(一)基本思想
设 X 1 , X 2 , . . . , X n X_1, X_2, . . . , X_n X1,X2,...,Xn 为来自总体分布 F ( x ; θ ) F(x; θ) F(x;θ) 的样本,对于感兴趣参数 θ θ θ,通常可采用极大似然估计或者矩估计等估计方法得到 θ θ θ 的估计 θ ^ ( X 1 , X 2 , . . . , X n ) \hat{θ}(X_1, X_2, . . . , X_n) θ^(X1,X2,...,Xn)。但是一般情况下很难推导出 θ ^ ( X 1 , X 2 , . . . , X n ) \hat{θ}(X_1, X_2, . . . , X_n) θ^(X1,X2,...,Xn)的分布,进而得到估计量的均值和方差。
- 如果 F ( x ; θ ) F(x; θ) F(x;θ) 分布形式已知,可以通过蒙特卡洛方法模拟出与 θ ^ ( X 1 , X 2 , . . . , X n ) \hat{θ}(X_1, X_2, . . . , X_n) θ^(X1,X2,...,Xn)同分布的样本,进而根据样本的信息估计 θ ^ ( X 1 , X 2 , . . . , X n ) \hat{θ}(X_1, X_2, . . . , X_n) θ^(X1,X2,...,Xn)的分布以及分布特征。
蒙特卡洛步骤:
- 从总体分布 F ( x ; θ ) F(x;\theta) F(x;θ)中独立产生 n n n个数据 X 11 , X 12 , . . . , X 1 n X_{11},X_{12},...,X_{1n} X11,X12,...,X1n,得到的 θ \theta θ的估计,记为 θ 1 ^ = : θ ^ 1 ( X 11 , X 12 , . . , X 1 n ) \hat{\theta_1}=:\hat{\theta}_1(X_{11},X_{12},..,X_{1n}) θ1^=:θ^1(X11,X12,..,X1n)
- 重复上述步骤 m m m次,得到的估计分别为 θ ^ 1 , θ ^ 2 , . . . , θ ^ m \hat{\theta}_1,\hat{\theta}_2,...,\hat{\theta}_m θ^1,θ^2,...,θ^m
- 基于样本 θ ^ 1 , θ ^ 2 , . . . , θ ^ m \hat{\theta}_1,\hat{\theta}_2,...,\hat{\theta}_m θ^1,θ^2,...,θ^m推断 θ ^ 1 ( X 11 , X 12 , . . , X 1 n ) \hat{\theta}_1(X_{11},X_{12},..,X_{1n}) θ^1(X11,X12,..,X1n)的分布, E ( θ ^ ) E(\hat{\theta}) E(θ^)和 V a r ( θ ^ ) Var(\hat{\theta}) Var(θ^)可以分别用该样本的均值和方差估计,即 E ^ ( θ ^ ) = m − 1 ∑ i = 1 m θ ^ i \hat{E}(\hat{\theta})=m^{-1}\sum_{i=1}^m\hat{\theta}_i E^(θ^)=m−1∑i=1mθ^i, V a r ^ ( θ ^ ) = m − 1 ∑ i = 1 m ( θ ^ i − E ^ ( θ ^ ) ) 2 \hat{Var}(\hat{\theta})=m^{-1}\sum_{i=1}^m(\hat{\theta}_i-\hat{E}(\hat{\theta}))^2 Var^(θ^)=m−1∑i=1m(θ^i−E^(θ^))2
- 如果 F ( x ; θ ) F(x; θ) F(x;θ) 分布形式未知,唯一的信息只有样本 X 1 , X 2 , . . . , X n X_1, X_2, . . . , X_n X1,X2,...,Xn。不能利用蒙特卡洛方法从总体 F ( x ; θ ) F(x;\theta) F(x;θ)中产生数据,进而不能近似 θ ^ ( X 1 , X 2 , . . . , X n ) \hat{\theta}(X_1, X_2, . . . , X_n) θ^(X1,X2,...,Xn)的分布及其相关特征。Bootstrap方法是利用样本分布 X 1 , X 2 , . . . X n X_1,X_2,...X_n X1,X2,...Xn代替总体分布,从分布 X 1 , X 2 , . . . , X n X_1, X_2, . . . , X_n X1,X2,...,Xn中有放回的产生数据,进而近似 θ ^ ( X 1 , X 2 , . . . , X n ) \hat{\theta}(X_1, X_2, . . . , X_n) θ^(X1,X2,...,Xn)的分布
Bootsrap方法的步骤:
- 从样本 X 1 , X 2 , . . . , X n X_1, X_2, . . . , X_n X1,X2,...,Xn中有放回的产生数据 X 11 ∗ , X 12 ∗ , . . . , X 1 n ∗ X_{11}^*, X_{12}^*, . . . , X_{1n}^* X11∗,X12∗,...,X1n∗,得到的 θ \theta θ的估计,记为 θ ^ 1 ∗ = : θ ^ 1 ∗ ( X 11 ∗ , X 12 ∗ , . . . , X 1 n ∗ ) \hat{\theta}_1^*=:\hat{\theta}_1^*(X_{11}^*, X_{12}^*, . . . , X_{1n}^*) θ^1∗=:θ^1∗(X11∗,X12∗,...,X1n∗)
- 重复上述步骤m次, 得到的估计分别为 θ ^ 1 ∗ , θ ^ 2 ∗ , . . . , θ ^ m ∗ \hat{\theta}_1^*,\hat{\theta}_2^*,...,\hat{\theta}_m^* θ^1∗,θ^2∗,...,θ^m∗
- 利用 θ ^ 1 ∗ , θ ^ 2 ∗ , . . . , θ ^ m ∗ \hat{\theta}_1^*,\hat{\theta}_2^*,...,\hat{\theta}_m^* θ^1∗,θ^2∗,...,θ^m∗近似 θ ^ \hat{\theta} θ^的分布及其特征
用 Bootstrap 方法抽取到样本 X 1 ∗ , X 2 ∗ , . . . , X n ∗ X_{1}^*, X_{2}^*, . . . , X_{n}^* X1∗,X2∗,...,Xn∗的经验分布 F n ∗ F^*_n Fn∗是 F n F_n Fn 的逼近, F n F_n Fn 是总体分布 F F F 的逼近。这两种逼近可以表示为 F n ∗ → F n → F F^∗_n → F_n → F Fn∗→Fn→F。如果经验分布函数 F n ( x ) F_n(x) Fn(x) 没有靠近总体分布函数 F ( x ) F(x) F(x),则重复抽样下的分布也不会靠近 F ( x ) F(x) F(x)。
(二)偏差的自助估计
1、参数 θ \theta θ的估计
估计量 θ ^ \hat{\theta} θ^的偏差: B i a s ( θ ^ ) = E ( θ ^ ) − θ Bias(\hat{\theta})=E(\hat{\theta})-\theta Bias(θ^)=E(θ^)−θ;估计量 θ ^ \hat{\theta} θ^偏差的自助估计: B i a s ^ ( θ ^ ) \widehat{Bias}(\hat{\theta}) Bias (θ^)
- 对于 θ \theta θ的估计 θ ^ \hat{\theta} θ^:基于观察到的样本数据直接得到
- 对于 E ( θ ^ ) E(\hat{\theta}) E(θ^)的估计 E ^ ( θ ^ ) \hat{E}(\hat{\theta}) E^(θ^):利用自助法获得的 θ ^ 1 ∗ , θ ^ 2 ∗ , . . . , θ ^ m ∗ \hat{\theta}_1^*,\hat{\theta}_2^*,...,\hat{\theta}_m^* θ^1∗,θ^2∗,...,θ^m∗的均值来估计,即 E ^ ( θ ^ ) = m − 1 ∑ i = 1 m θ ^ i ∗ \hat{E}(\hat{\theta})=m^{-1}\sum_{i=1}^m\hat{\theta}_i^* E^(θ^)=m−1∑i=1mθ^i∗
B i a s ^ ( θ ^ ) = m − 1 ∑ i = 1 m θ ^ i ∗ − θ ^ \widehat{Bias}(\hat{\theta})=m^{-1}\sum_{i=1}^m\hat{\theta}_i^*-\hat{\theta} Bias (θ^)=m−1i=1∑mθ^i∗−θ^
如果偏差大于0,说明 θ ^ \hat{\theta} θ^平均来看过高估计了 θ \theta θ;偏差小于0则说明 θ ^ \hat{\theta} θ^平均来看过低估计了 θ \theta θ。因此,经过偏差修正的参数 θ \theta θ的估计量为:
θ ~ = θ ^ − B i a s ^ ( θ ^ ) \tilde{\theta}=\hat{\theta}-\widehat{Bias}(\hat{\theta}) θ~=θ^−Bias (θ^)
2、参数 σ 2 \sigma^2 σ2的估计
设数据 X 1 , X 2 , . . . , X n X_1,X_2,...,X_n X1,X2,...,Xn来自方差为 σ 2 \sigma^2 σ2的分布,则 σ 2 \sigma^2 σ2的估计为 σ ^ = n − 1 ∑ i = 1 n ( X i − X ˉ ) 2 \hat{\sigma}=n^{-1}\sum_{i=1}^n(X_i-\bar{X})^2 σ^=n−1∑i=1n(Xi−Xˉ)2。 σ 2 \sigma^2 σ2的偏差 B i a s ( σ ^ 2 ) = E ( σ ^ 2 ) − σ 2 = − σ 2 / n Bias(\hat{\sigma}^2)=E(\hat{\sigma}^2)-\sigma^2=-\sigma^2/n Bias(σ^2)=E(σ^2)−σ2=−σ2/n。
此时 σ ^ 2 \hat{\sigma}^2 σ^2偏差的自助估计 B i a s ^ ( σ ^ 2 ) = m − 1 ∑ i = 1 m σ ^ i 2 ∗ − σ ^ 2 \widehat{Bias}(\hat{\sigma}^2)=m^{-1}\sum_{i=1}^m\hat{\sigma}_i^{2*}-\hat{\sigma}^2 Bias (σ^2)=m−1∑i=1mσ^i2∗−σ^2是否为 − σ 2 / n -\sigma^2/n −σ2/n的无偏估计?
因此 σ ^ 2 \hat{\sigma}^2 σ^2偏差的自助估计 B i a s ^ ( σ ^ 2 ) = m − 1 ∑ i = 1 m σ ^ i 2 ∗ − σ ^ 2 \widehat{Bias}(\hat{\sigma}^2)=m^{-1}\sum_{i=1}^m\hat{\sigma}_i^{2*}-\hat{\sigma}^2 Bias (σ^2)=m−1∑i=1mσ^i2∗−σ^2不是 − σ 2 / n -\sigma^2/n −σ2/n的无偏估计。
(三)估计量标准差的 Bootstrap 估计
估计量 θ ^ \hat{\theta} θ^的标准差 S E ^ ( θ ^ ) \widehat{SE}(\hat{\theta}) SE (θ^)的自助估计就是采用 θ ^ \hat{\theta} θ^分布Bootstrap估计的标准差。
S E ^ ( θ ^ ) = ( 1 m ∑ i = 1 m ( θ ^ i ∗ − θ ^ ˉ ∗ ) 2 ) 1 / 2 \widehat{SE}(\hat{\theta})=\left(\frac{1}{m}\sum_{i=1}^m(\hat{\theta}_i^*-\bar{\hat{\theta}}^*)^2\right)^{1/2} SE (θ^)=(m1i=1∑m(θ^i∗−θ^ˉ∗)2)1/2
自助法得到的估计 V a r ^ ( θ ^ ) \widehat{Var}(\hat{\theta}) Var (θ^)和 V a r ( θ ^ ) Var(\hat{\theta}) Var(θ^)之间的关系是什么?
因此,自助法得到的估计 V a r ^ ( θ ^ ) \widehat{Var}(\hat{\theta}) Var (θ^)的条件期望为 m − 1 m n σ ^ 2 \frac{m-1}{mn}\hat{\sigma}^2 mnm−1σ^2,是方差 V a r ( θ ^ ) = n − 1 σ 2 Var(\hat{\theta})=n^{-1}\sigma^2 Var(θ^)=n−1σ2的估计。
二、基于 Jackknife 法的估计
(一)基本思想
Jackknife 估计的基本思想是,对于给定样本 X 1 , X 2 , . . . , X n X_1, X_2, . . . , X_n X1,X2,...,Xn,每次删除其中一个 (或者几个) 样本点,基于剩下的样本采用相同的估计量公式得到 θ θ θ 的估计,经过逐个删除并分别计算估计之后,便可以得到一系列估计值,基于这些估计值进而估计 θ ^ \hat{θ} θ^的分布特征。
Jackknife 方法的步骤:
- 从观测样本 X 1 , X 2 , . . . , X n X_1, X_2, . . . , X_n X1,X2,...,Xn中去掉第 i i i个数据 X i X_i Xi之后的剩余样本定义为第 i i i个Jackknife样本,记为: X ( − i ) = ( X 1 , . . . , X i − 1 , X i = 1 , . . . , X n ) X_{(-i)}=(X_1,...,X_{i-1},X_{i=1},...,X_n) X(−i)=(X1,...,Xi−1,Xi=1,...,Xn)
- 基于第 i i i个Jackknife样本 X ( − i ) , i = 1 , 2 , . . . , n X_{(-i)},i=1,2,...,n X(−i),i=1,2,...,n,得到相应的估计 θ ^ ( − i ) = θ ^ ( X ( − i ) ) , i = 1 , . . . , n \hat{\theta}_{(-i)}=\hat{\theta}(X_{(-i)}),i=1,...,n θ^(−i)=θ^(X(−i)),i=1,...,n
(二)估计量偏差的 Jackknife 估计
基于 θ ( − i ) ^ \hat{\theta_{(-i)}} θ(−i)^, θ ^ \hat{\theta} θ^的偏差 B i a s ( θ ^ ) = E ( θ ^ ) − θ Bias(\hat{\theta})=E(\hat{\theta})-\theta Bias(θ^)=E(θ^)−θ的Jackknife估计为:
B i a s ^ ( θ ^ ) = ( n − 1 ) ( 1 n ∑ i = 1 n θ ^ ( − i ) − θ ^ ) = n − 1 n ∑ i = 1 n ( θ ^ ( − i ) − θ ^ ) \widehat{Bias}(\hat{\theta})=(n-1)\left(\frac{1}{n}\sum_{i=1}^n\hat{\theta}_{(-i})^-\hat{\theta}\right)=\frac{n-1}{n}\sum_{i=1}^n(\hat{\theta}_{(-i)}-\hat{\theta}) Bias (θ^)=(n−1)(n1i=1∑nθ^(−i)−θ^)=nn−1i=1∑n(θ^(−i)−θ^)
设数据 X 1 , X 2 , . . . , X n X_1,X_2,...,X_n X1,X2,...,Xn来自方差为 σ 2 \sigma^2 σ2的分布,则 σ 2 \sigma^2 σ2的估计为 σ ^ = n − 1 ∑ i = 1 n ( X i − X ˉ ) 2 \hat{\sigma}=n^{-1}\sum_{i=1}^n(X_i-\bar{X})^2 σ^=n−1∑i=1n(Xi−Xˉ)2。 σ 2 \sigma^2 σ2的偏差 B i a s ( σ ^ 2 ) = E ( σ ^ 2 ) − σ 2 = − σ 2 / n Bias(\hat{\sigma}^2)=E(\hat{\sigma}^2)-\sigma^2=-\sigma^2/n Bias(σ^2)=E(σ^2)−σ2=−σ2/n。此时,对于每一个Jackknife估计 θ ^ ( − i ) \hat{\theta}_{(-i)} θ^(−i),基于样本量n-1的样本 X ( − i ) = ( X 1 , . . . , X i − 1 , X i = 1 , . . . , X n ) X_{(-i)}=(X_1,...,X_{i-1},X_{i=1},...,X_n) X(−i)=(X1,...,Xi−1,Xi=1,...,Xn)构造,因此:
(三)估计量标准差的 Jackknife 估计
估计量 θ ^ \hat{\theta} θ^的标准差的Jackknife估计定义为:
S E ^ J a c k ( θ ^ ) = ( n − 1 n ∑ i = 1 n ( θ ^ ( − i ) − θ ^ ˉ ( ⋅ ) ) 2 ) 1 / 2 \widehat{SE}_{Jack}(\hat{\theta})=\left(\frac{n-1}{n}\sum_{i=1}^n(\hat{\theta}_{(-i)}-\bar{\hat{\theta}}_{(·)})^2\right)^{1/2} SE Jack(θ^)=(nn−1i=1∑n(θ^(−i)−θ^ˉ(⋅))2)1/2
三、Permutation Test 置换检验
置换检验(permutation test)是统计学上一种基于反证法、重抽样原则的非参数性检验。
- H0 : F = G, 即所有样本都服从同一分布
- H1 : F ≠ \neq = G, 即样本不服从同一分布
置换检验通过对比样本置换后的检验统计量与置换前的检验统计量来决定是否拒绝零假设。p 值为假设检验中假设零假设为真时观测到的至少与实际观测样本相同的样本的概率。很小的 p 值说明在零假设下观测到的概率很小。
置换检验的步骤:
- 首先计算两样本(样本容量设为 n A n_A nA 和 n B n_B nB)之间原本的检验统计量。检验统计量可以是两样本间平均数之差 ( X ˉ A − X ˉ B ) (\bar{X}_A − \bar{X}_B) (XˉA−XˉB)、方差之差 ( S A 2 − S B 2 ) (S^2_A − S^2_B) (SA2−SB2),或 t 值 ( t ) (t) (t)、卡方检验中的卡方值 ( χ 2 ) (χ^2) (χ2) 等。
- 将两个样本打乱后再重新选出两组容量等于之前两样本的新样本(即两个样本容量同样为 nA 和 nB 的样本),并计算新的检验统计量。
- 如接受零假设 H0 : F = G,即样本源于同一分布,则随机抽样计算出的新检验统计量应不难大于最初置换前算出的两样本间检验统计量(如为双侧检验,则是其绝对值应不难大于置换前算出的两样本间检验统计量),即这个概率应大于设定的 I 型错误(假阳性)概率 α。