中国科学技术大学：《实用统计软件》课程课件讲义（统计计算与软件）第六讲 Monte Carlo方法在统计推断中的应用

因此灭标准差的估计为 s()= a[品-时]=e- 或者也可以使用无偏估计量 1 因此，前例中的标准差估计为 = 计算程序如下 sqrt(sum((g-mean(g))2))/m #sd(g)/sqrt(m)#for unbiased estimator Code Previous Next First Last Back Forward 4

œdX¯IOOè seˆ (¯x) = 1 √ m  1 m (xi − x¯) 2 1/2 = 1 m (xi − x¯) 2 1/2 . ½ˆè屶^ÆO˛ seˆ (¯x) = 1 √ m  1 m − 1 (xi − x¯) 2 1/2 . œd, c~•IOOè seˆ (θˆ) = 1 m  Xm j=1 (θˆ(j) − θˆ) 2   . OéßSXe ↑Code sqrt(sum((g-mean(g))^2))/m #sd(g)/sqrt(m) #for unbiased estimator ↓Code Previous Next First Last Back Forward 4

1.1.2 Estimation of MSE Monte Carlo方法可以用于计算一个估计量的MSE.一个估计量的MSE定 义为MSE(©=E6-.如果从总体X中产生了m个样本x(1,·,x(m), 则e的MSE的Monte Carlo估计为 1 ME=(6)-)2, m J=1 其中U)=(x). 例2：使用Monte Carlo方法估计标准正态分布的截尾均值-1的MSE. 当样本存在异常点时，截尾的样本均值常常被用来估计总体的中心.假 设X1,·,Xn为一个随机样本，X()·,X(n)为相应的次序统计量，则一 个k水平的截尾样本均值为 1 n一 -= n-2k ∑xo i=k+1 Previous Next First Last Back Forward 5

1.1.2 Estimation of MSE Monte Carlo ê{å±^uOéòáO˛MSE. òáO˛MSE½ ¬èMSE(θˆ) = E[θˆ − θ] 2 . XJloNX•) máx (1) , · · · , x(m) , KθˆMSEMonte CarloOè MSEˆ = 1 m Xm j=1 (θˆ(j) − θ) 2 , Ÿ•θˆ(j) = θˆ(x (j) ). ~ 2: ¶^Monte Carloê{OIO©Ÿó˛äX¯ [−1]MSE. 3…~:û, ó˛ä~~^5OoN•%. b X1, · · · , XnèòáëÅ, X(1), · · · , X(n)èÉAgS⁄O˛, K ò ákY²ó˛äè X¯ [−k] = 1 n − 2k nX−k i=k+1 X(i) . Previous Next First Last Back Forward 5

本例中，目标参数为9=E了=EX-1)=0.记T=X-1,则其MSE的 一个Monte Carlo估计算法如下 1.通过如下步骤产生m个重复T),j=1,…,m: 。产生总体X的样本：x,, 。从小到大排序唱≤…≤名 ·计算T0)=高∑8 2.计算MSE MSE(T)=六∑m1(T)-P 六(T6)2. 实现程序如下 TCode n<-20 m<-1000 tmean <-numeric(m) Previous Next First Last Back Forward 6

~•, 8IÎÍèθ = EX¯ = EX¯ [−1] = 0. PT = X¯ [−1], KŸMSE òá Monte Carlo Oé{Xe 1. œLXe⁄½)má­ET (j) , j = 1, · · · , m: • )oNX: x (j) 1 , · · · , x (j) n . • lå¸S x (j) (1) ≤ · · · ≤ x (j) (n) • Oé T (j) = 1 n−2 Pn−1 i=2 x (j) (i) . 2. O éMSE MSEˆ (T) = 1 m Pm j=1(T (j) − θ) 2 = 1 m (T (j) ) 2 . ¢yßSXe ↑Code n <- 20 m <- 1000 tmean <- numeric(m) Previous Next First Last Back Forward 6

for (i in 1:m){ x <-sort(rnorm(n)) tmean[1]<-sum(x[2:(n-1)])/(n-2) mse <-mean(tmean"2) mse sqrt(sum((tmean mean(tmean))2))/m #s日 Code 截尾均值的MSE的估计为0.0504531(se÷0.007).样本均值的MSE为Var(X)/n= 1/20=0.05.另一方面，中位数本质上也是一种截尾均值，其截掉了除中间的 一个或两个点外的所有点.样本中位数的MSE估计如下 TCode n<-20 m<-1000 tmean <-numeric(m) for (i in 1:m){ x <-sort(rnorm(n)) tmean[i]<-median(x) ] Previous Next First Last Back Forward 7

for (i in 1:m) { x <- sort(rnorm(n)) tmean[i] <- sum(x[2:(n-1)]) / (n-2) } mse <- mean(tmean^2) mse sqrt(sum((tmean - mean(tmean))^2)) / m #se ↓Code ó˛äMSEOè0.0504531( ˆse .= 0.007). ˛ä MSEèV ar(X)/n = 1/20 = 0.05. ,òê°, •†Íü˛è¥ò´ó˛ä, ŸK ÿ•m òὸá: §k:. •†ÍMSEOXe ↑Code n <- 20 m <- 1000 tmean <- numeric(m) for (i in 1:m) { x <- sort(rnorm(n)) tmean[i] <- median(x) } Previous Next First Last Back Forward 7

mse <mean(tmean2) mse sqrt(sum((tmean mean(tmean))2))/m #se Code 从而样本中位数的MSE估计为0.075(s2=0.0086). 例3：比较标准正态分布与如下混合(“污染”)正态分布下的k水平戴尾均值 估计的MSE. pN(0,a2=1)+(1-p)N(0,a2=100) 我们写一个函数来对不同的k和P计算截尾均值X-的MSE.从混合正 态中产生样本时，要根据P(a=1)=p:P(a=10)=1-p来随机选 择a.注意正态随机数产生函数rnorm可以使用参数向量作为标准偏差.考 虑p=1.0,0.95,0.9以及k=0,1,·,n/2.因此程序如下 fccd加 n<-20;K<-n/2-1;m<-1000; mse <-matrix(0,n/2,6) Previous Next First Last Back Forward 8

mse <- mean(tmean^2) mse sqrt(sum((tmean - mean(tmean))^2)) / m #se ↓Code l •†ÍMSEOè0.075( ˆse = 0.0086). ~ 3: 'IO©ŸÜXe·‹(“¿/”)©ŸekY²ó˛ä OMSE. pN(0, σ2 = 1) + (1 − p)N(0, σ2 = 100) ·ÇòáºÍ5Èÿ”k⁄pOéó˛äX¯ [−k]MSE. l·‹ •) û, áä‚P(σ = 1) = p; P(σ = 10) = 1 − p5ëÅ¿ Jσ. 5øëÅÍ) ºÍrnorm屶^ÎÍï˛äèIO†.  ƒp = 1.0, 0.95, 0.9 ±9k = 0, 1, · · · , n/2. œdßSXe ↑Code n <- 20; K <- n/2 - 1; m <- 1000; mse <- matrix(0, n/2, 6) Previous Next First Last Back Forward 8