數理統計12：樞軸量法、分位數、正態參數區間估計

時間 2021-02-17

標籤函數 spa code orm get class pdf 變量欄目 Java 简体版

原文原文鏈接

上篇文章中，咱們探討了區間估計的相關基本概念，也提出了Neyman置信區間，今天咱們將聚焦於如何尋找置信區間的問題上，並對最經常使用的整體：正態整體給出一些置信區間的找法。爲了方便起見，如下咱們都讓置信水平爲\(1-\alpha\)。函數

因爲本系列爲我獨自完成的，缺乏審閱，若是有任何錯誤，歡迎在評論區中指出，謝謝！spa

Part 1：樞軸量法

樞軸變量法是基於點估計量的。咱們知道，統計量是樣本的函數，這意味着統計量中不能含有未知參數，而參數的點估計量是用統計量的觀測值做爲待估參數的估計值，其分佈必定含有待估參數，樞軸量法的思想就是，經過必定的變換，讓點估計的函數的分佈不含待估參數，進而基於分佈來構造區間估計。code

舉一個簡單的例子，對於正態整體\(N(\mu,4)\)，顯然\(\bar X\sim N(\mu,4/n)\)，這裏\(\bar X\)的分佈含有未知參數\(\mu\)。構造其樞軸量，就是找到一個函數變換，使得新的隨機變量分佈不含未知參數。注意，這裏用了隨機變量這個詞而不是統計量，意味着樞軸量不是統計量，即不能由樣本觀測值計算出，這是由於雖然樞軸量的分佈不含未知參數，可是樞軸量的表現形式含有未知參數。顯然，這裏orm

\[\bar X-\mu\sim N(0,\frac{4}{n}), \]

這樣，\(\bar X-\mu\)的分佈已知，天然容易找到一個常數區間\([c,d]\)，使得這個區間有\(1-\alpha\)的機率包含\(\bar X-\mu\)的觀測值，雖然此時咱們不知道區間的端點是多少，但至少知道端點能夠是固定的數\(c,d\)。對樞軸量使用不等式變換，即\(\bar X-\mu\in[c,d]\Rightarrow \mu\in[\bar X-d,\bar X-c]\)，獲得置信水平爲\(1-\alpha\)的置信區間。這就是樞軸量法的操做步驟。get

不一樣分佈族的參數對於整體的意義是不一樣的。像正態分佈\(N(\mu,\sigma^2)\)的均值\(\mu\)，均勻分佈\(U(a,a+r)\)的起點\(a\)這種參數主要影響觀測值的大小，能夠直接經過\(X-\mu\)，\(X-a\)的手段消除，這種參數稱爲位置參數；正態分佈\(N(\mu,\sigma^2)\)的標準差\(\sigma\)，指數分佈\(E(\lambda)\)的速率\(\lambda\)這種參數主要影響觀測值的離散程度，能夠經過\(X/\sigma\)，\(\lambda X\)之類的手段消除，這種參數稱爲尺度參數。面對不一樣的參數，每每有不一樣的方式構造樞軸量，應當結合其特色選擇構造方式。qt

接下來，咱們將樞軸量法運用到一些具體實例中，體會樞軸量法的實際應用。io

Part 2：分位數

在開始樞軸量法的實際應用前，先介紹分位數這一律念，這爲咱們肯定置信區間提供了有力武器。如今，咱們先給出分位數的定義（若是不特別說明，都指整體分位數而非樣本分位數）。分位數能夠分爲下分位數和上分位數，通常稱下分位數爲分位數，但咱們更經常使用的是上分位數。class

對於某個分佈\(F\)，\(X\sim F\)，\(F\)的（下）\(\alpha\)分位數指的是這樣一個點\(x_\alpha\)，知足pdf

\[\mathbb{P}(X\le x_\alpha)=\alpha, \]

若是\(F\)有反函數\(F^{-1}(x)\)，則有\(x_\alpha=F^{-1}(\alpha)\)。變量

對於某個分佈\(F\)，\(X\sim F\)，\(F\)的上\(\alpha\)分位數指的是這樣一個點\(y_\alpha\)，知足

\[\mathbb{P}(X\ge y_\alpha)=\alpha. \]

換言之，若是分佈函數\(F\)存在反函數\(F^{-1}(x)\)，則\(F\)的上\(\alpha\)分位數是\(y_\alpha=F^{-1}(1-\alpha)\)。

也許結合圖形會更好理解。對於上面的圖形，若是這是一個密度函數，黑色部分的面積爲\(0.05\)，那麼紅色的點就是該分佈的上\(0.05\)分位數，同時也是其下\(0.95\)分位數。

顯然，對於一個徹底肯定的分佈，其各分位數都是徹底肯定的，是常數。分佈多種多樣，不可能對全部分佈都詳細地討論其分位數，可是對一些經常使用分佈，其各分位數的數值已經被製成表，這包括標準正態分佈，與正態分佈的三大衍生分佈——\(\chi^2\)分佈、\(t\)分佈、\(F\)分佈。

從此，咱們將使用\(u_\alpha\)、\(\chi^2_\alpha(n)\)、\(t_\alpha(n)\)、\(F_{\alpha}(m,n)\)，分別表明標準正態分佈、自由度爲\(n\)的\(\chi^2\)分佈、自由度爲\(n\)的\(t\)分佈、自由度爲\(m,n\)的\(F\)分佈的上\(\alpha\)分位數，給定分佈類型、分位、自由度後，這些值均可以經過查閱分位數表獲得。

Part 3：單正態整體參數區間估計

正態分佈參數的區間估計，主要是對\(\mu\)和\(\sigma^2\)給出區間估計。其中，對\(\mu\)的估計又要分紅\(\sigma^2\)已知和未知兩種狀況。不過，你們已經知道，與參數\(\mu\)關聯最緊密的點估計是\(\bar X\)，與參數\(\sigma^2\)關聯最緊密的點估計是\(S^2\)，所以樞軸量法也必定圍繞着它們做文章。

事實上，只要是均值，基本均可以轉化爲一個正態分佈；只要是方差，基本上均可以轉化爲一個卡方分佈。而後，根據\(t\)分佈和\(F\)分佈的構造方式，就能構造出服從這四種分佈的樞軸量來。

一、\(\sigma^2\)已知時\(\mu\)的區間估計

當\(\sigma^2\)已知時，

\[\bar X\sim N(\mu,\sigma^2), \]

樞軸量顯然是\(\bar X-\mu\sim N(0,\sigma^2/n)\)，已經由此肯定了一個已知分佈，可是爲了寫出區間估計的具體表達，咱們還應當將這個已知的正態分佈，轉化爲標準正態分佈。爲此，咱們選擇樞軸量爲

\[\frac{\sqrt{n}(\bar X-\mu)}{\sigma}\sim N(0,1), \]

找一個以\(1-\alpha\)機率涵蓋標準正態分佈觀測的區間，天然會找到

\[[-u_{\alpha/2},u_{\alpha/2}], \]

即

\[-u_{\alpha/2}\le \frac{\sqrt{n}(\bar X-\mu)}{\sigma}\le u_{\alpha/2}, \]

進行恆等變換，就獲得了\(\mu\)的\(1-\alpha\)置信區間爲

\[\bar X-\frac{\sigma}{\sqrt{n}}u_{\alpha/2}\le \mu\le \bar X+\frac{\sigma}{\sqrt{n}}u_{\alpha/2},\\ \left[\bar X-\frac{\sigma}{\sqrt{n}}u_{\alpha/2},\bar X+\frac{\sigma}{\sqrt{n}}u_{\alpha/2} \right]. \]

二、\(\sigma^2\)未知時\(\mu\)的區間估計

當\(\sigma^2\)未知時，咱們不能順利地將\(\bar X\)標準化，但應該容易想到，用\(S^2\)替代\(\sigma^2\)能起到相似的效果，也就是構造這樣的變量：\(\sqrt{n}(\bar X-\mu)/S\)。

神奇的是，\(\frac{\sqrt{n}(\bar X-\mu)}{S}\)能夠被改寫成\(\chi^2\)分佈的形式：

\[\frac{\sqrt{n}(\bar X-\mu)}{S}=\frac{\frac{\sqrt{n}(\bar X-\mu)}{\sigma}}{\sqrt{\frac{(n-1)S^2}{\sigma^2}/(n-1)}}\stackrel{d}=\frac{N(0,1)}{\sqrt{\chi^2(n-1)/(n-1)}}\sim t(n-1), \]

所以，天然能夠找到一個區間\(\frac{\sqrt{n}(\bar X-\mu)}{S}\in[-t_{\alpha/2}(n-1),t_{\alpha/2}(n-1)]\)，所以獲得\(\mu\)的\(1-\alpha\)置信區間爲

\[\left[\bar X-\frac{S}{\sqrt{n}}t_{\alpha/2}(n-1),\bar X+\frac{S}{\sqrt{n}}t_{\alpha/2}(n-1) \right]. \]

用R語言進行模擬，如今咱們每次生成10個服從\(N(5,9)\)的簡單隨機樣本，這一組樣本能夠生成一個對應的\(95\%\)置信區間（假設\(\sigma^2=9\)是未知的）。進行10000次試驗，生成10000個這樣的區間，觀察這個區間會涵蓋待估參數\(\mu=5\)多少次。

rm(list = ls())  # 清空工做區
dev.off()  # 清空圖形區

n <- 10
t <- qt(0.975, df=n-1)  # t的上0.025分位數
upper.bound <- c()
lower.bound <- c()

for (j in 1:10000){
  xlst <- rnorm(10, 5, 3)
  M <- mean(xlst)
  S <- sqrt(var(xlst))
  upper.bound[j] <- M+S*t/sqrt(n-1)
  lower.bound[j] <- M-S*t/sqrt(n-1)
}

plot(1:10000, rep(5, 10000), main = "Mean of Norm", ylim = c(-10, 20), cex=0.1)
points(1:10000, upper.bound, col='red', cex=0.1)
points(1:10000, lower.bound, col='blue', cex=0.1)
sum(upper.bound < 5 | lower.bound > 5)  # 輸出區間不涵蓋待估參數的次數

輸出結果爲，10000次中，有384次沒有涵蓋待估參數，圖示以下：

三、\(\mu\)未知時\(\sigma^2\)的區間估計

咱們在給出\(S^2\)的分佈信息時，實際上已經構建出了樞軸量：

\[\frac{(n-1)S^2}{\sigma^2}\sim \chi^2(n-1). \]

因此咱們能夠給出相應的\(1-\alpha\)機率區間，不過要注意的是\(\chi^2\)分佈不像正態分佈、\(t\)分佈同樣，是對稱分佈，因此機率區間爲\(\frac{(n-1)S^2}{\sigma^2}\in[\chi^2_{1-\alpha/2}(n-1),\chi^2_{\alpha/2}(n-1)]\)，獲得\(\sigma^2\)的置信區間爲

\[\left[\frac{(n-1)S^2}{\chi^2_{\alpha/2}(n-1)},\frac{(n-1)S^2}{\chi^2_{1-\alpha/2}(n-1)} \right] \]

使得機率爲\(1-\alpha\)的區間理論上有無窮多個，它們有相同的置信度，本該選擇其中精度最高（即長度最短）的，可是這樣會很麻煩，爲方便起見，選取\(\alpha/2\)分位數和\(1-\alpha/2\)分位數構造置信區間，這樣的置信區間稱爲等尾置信區間。

四、\(\mu\)已知時\(\sigma^2\)的區間估計

看起來，\(S^2\)這個樞軸量不考慮\(\mu\)是否已知，可是樞軸量法是依賴於點估計的，天然咱們首先會考慮的是充分統計量，而\(S^2\)在\(\mu\)已知時，不是\(\sigma^2\)的充分統計量。爲此，咱們要基於\(\sigma^2\)的充分統計量構建樞軸量。令

\[Q=\sum_{j=1}^n(X_j-\mu)^2, \]

它是\(\sigma^2\)的充分統計量，且咱們能夠將其轉化爲\(\chi^2\)分佈：

\[\frac{Q}{\sigma^2}=\sum_{j=1}^n\left(\frac{X_j-\mu}{\sigma} \right)^2\sim \chi^2(n), \]

這是因爲各個\(X_j\)相互獨立服從\(N(\mu,\sigma^2)\)，因而能夠構造出機率區間\(\frac{Q}{\sigma^2}\in[\chi^2_{1-\alpha/2}(n),\chi^2_{\alpha/2}(n)]\)，所以\(\sigma^2\)的置信區間是

\[\left[\frac{Q}{\chi^2_{\alpha/2}(n)},\frac{Q}{\chi^2_{1-\alpha/2}(n)} \right]. \]

Part 4：雙正態整體參數區間估計

雙整體\(X\sim N(\mu_1,\sigma_1^2)\)、\(Y\sim N(\mu_2,\sigma_2^2)\)的參數估計，主要指的是均值差和方差比，可是要分狀況討論，這是由於若是不加任何限制，它們的置信區間理論還不是很完善。

如下記\(X_1,\cdots,X_m\)是從\(X\)中抽取的樣本，相應的樣本均值和樣本方差是\(\bar X\)和\(S_x^2\)；\(Y_1,\cdots,Y_n\)是從\(Y\)中抽取的樣本，相應的樣本均值和樣本方差是\(\bar Y\)和\(S_y^2\)。

一、均值差\(\mu_2-\mu_1\)的區間估計

第一種狀況是\(\sigma_1^2,\sigma_2^2\)未知但\(m=n\)時，即樣本容量相等時。一種簡單的作法是構造\(Z_j=Y_j-X_j\)，即讓成對數據相減，這樣\(Z_1,\cdots,Z_n\)便獨立同分佈於\(N(\mu_2-\mu_1,\sigma_1^2+\sigma_2^2)\)，至關於從單整體中抽樣，並對整體均值做估計。顯然，樞軸量應該是

\[\frac{\sqrt{n}[\bar Z-(\mu_2-\mu_1)]}{S_z}\sim t(n-1). \]

因此置信區間是

\[\left[\bar Z-\frac{S_z}{\sqrt{n}}t_{\alpha/2}(n-1),\bar Z+\frac{S_z}{\sqrt{n}}t_{\alpha/2}(n-1) \right]. \]

不過要注意的是，構形成對數據時，實際上丟失了部分信息，這是由於\(Z_1,\cdots,Z_n\)並不是\(\mu_2-\mu_1\)的充分統計量。

總而言之，當\(m=n\)時，咱們將其轉化爲成對數據，本質上仍是單整體參數估計。

第二種狀況是\(\sigma_1^2,\sigma_2^2\)已知時，此時用\(\bar Y-\bar X\)來最大程度地估計\(\mu_2-\mu_1\)，則有

\[\bar Y-\bar X\sim N\left(\mu_2-\mu_1,\frac{\sigma_1^2}{m}+\frac{\sigma_2^2}{n} \right), \]

因而樞軸量是

\[\frac{\bar Y-\bar X-(\mu_2-\mu_1)}{\sqrt{\frac{\sigma_1^2}{m}+\frac{\sigma_2^2}{n}}}\sim N(0,1), \]

\(\mu_2-\mu_1\)的置信區間爲

\[\left[\bar Y-\bar X-u_{\alpha/2}\sqrt{\frac{\sigma_1^2}{m}+\frac{\sigma_2^2}{n}},\bar Y+\bar X+u_{\alpha/2}\sqrt{\frac{\sigma_1^2}{m}+\frac{\sigma_2^2}{n}} \right]. \]

第三種狀況稍微複雜一些，也是最常處理的一種狀況，即\(\sigma_1^2=\sigma_2^2=\sigma^2\)但未知時。此時仍有

\[\bar Y-\bar X\sim N\left(\mu_2-\mu_1,\sigma^2\sqrt{\frac{1}{m}+\frac{1}{n}} \right), \]

爲了解決\(\sigma^2\)未知的問題，也要用適當的樣本方差來替代，但此時如此構造：

\[Q=\sum_{j=1}^m(X_j-\bar X)^2+\sum_{j=1}^n(Y_j-\bar Y)^2, \]

能夠發現

\[\frac{Q}{\sigma^2}=\sum_{j=1}^m\left(\frac{X_j-\bar X}{\sigma} \right)^2+\sum_{j=1}^n\left(\frac{Y_j-\bar Y}{\sigma} \right)^2\stackrel{d}=\chi^2(m-1)+\chi^2(n-1)\stackrel{d}=\chi^2(n+m-2). \]

因而用\(Q/(n+m-2)\)代替\(\sigma^2\)最合適不過了，所以

\[\frac{\frac{\bar Y-\bar X-(\mu_2-\mu_1)}{\sigma\sqrt{\frac{1}{m}+\frac{1}{n}}}}{\sqrt{\frac{Q}{\sigma^2}/(n+m-2)}}=\sqrt{\frac{mn(m+n-2)}{m+n}}\frac{(\bar Y-\bar X)-(\mu_2-\mu_1)}{\sqrt{Q}}\sim t^2(m+n-2), \]

即\(\mu_2-\mu_1\)的置信區間爲

\[\left[(\bar Y-\bar X)-t_{\alpha/2}(n+m-2)\sqrt{\frac{Q(m+n)}{(m+n-2)mn}},(\bar Y-\bar X)+t_{\alpha/2}(n+m-2)\sqrt{\frac{Q(m+n)}{(m+n-2)mn}} \right]. \]

現進行模擬，每次實驗中，從\(N(3,9)\)中抽取5個樣本，從\(N(5,25)\)中抽取8個樣本，對均值差進行估計。進行10000次重複試驗。

rm(list = ls())
dev.off()

m <- 5
n <- 8
t <- qt(0.975, m+n-2)
upper.bound <- c()
lower.bound <- c()

for (j in 1:10000){
  xlst <- rnorm(m, 3, 3)
  ylst <- rnorm(n, 5, 5)
  My <- mean(ylst)
  Mx <- mean(xlst)
  Q <- (m-1)*var(xlst) + (n-1)*var(ylst)
  upper.bound[j] <- My-Mx+t*sqrt( (Q*(m+n)) / ((m+n-2)*m*n) )
  lower.bound[j] <- My-Mx-t*sqrt( (Q*(m+n)) / ((m+n-2)*m*n) )
}

plot(1:10000, rep(2,10000), xlab = "Experiment", ylab = "CI", main = "CI of the mean difference", cex=0.1, ylim = c(-13, 17))
points(1:10000, upper.bound, col = "red", cex = 0.1)
points(1:10000, lower.bound, col = "blue", cex = 0.1)
sum(upper.bound < 2 | lower.bound > 2)  # 輸出結果爲302

當咱們將樣本容量提高至\(m=200\)，\(n=100\)時，出現了929個不包含真值的區間。看起來區間估計是否包含真值，與\(m,n\)的相對大小有很大的關係。

當\(\sigma_1^2\ne \sigma_2^2\)且未知時，均值差的參數估計稱爲Behrens-Fisher問題，比較複雜，這裏不加討論。

二、方差比\(\sigma_2^2/\sigma_1^2\)的區間估計

事實上方差比比較好討論，由於方差的估計量必定是本身整體中所抽取的樣本方差，故服從某個\(\chi^2\)分佈，做比後應當能獲得\(F\)分佈，這裏簡要地討論兩種狀況。

首先是\(\mu_1,\mu_2\)都未知的狀況，此時天然\(S_x^2,S_y^2\)會被用做估計量，因爲

\[\frac{(m-1)S_x^2}{\sigma_1^2}\sim \chi^2(m-1),\quad \frac{(n-1)S_y^2}{\sigma_2^2}\sim \chi^2(n-1), \]

因此

\[\frac{\frac{(m-1)S_x^2}{\sigma_1^2}/(m-1)}{\frac{(n-1)S_y^2}{\sigma_2^2}/(n-1)}=\frac{S_x^2}{S_y^2}\frac{\sigma_2^2}{\sigma_1^2}\sim F(m-1,n-1), \]

故應有\(\frac{S_x^2\cdot\sigma_2^2}{S_y^2\cdot\sigma_1^2}\in[F_{1-\alpha/2}(m-1,n-1),F_{\alpha/2}(m-1,n-1)]\)，因此置信區間爲

\[\left[\frac{S_y^2}{S_x^2}F_{1-\alpha/2}(m-1,n-1),\frac{S_y^2}{S_x^2}F_{\alpha/2}(m-1,n-1) \right]. \]

進行模擬，每次實驗中，從\(N(3,9)\)中抽取9個樣本，從\(N(5,25)\)中抽取25個樣本，對方差比進行估計。進行10000次重複試驗。

rm(list = ls())
dev.off()

m <- 9
n <- 25
f1 <- qf(0.025, m-1, n-1)
f2 <- qf(0.975, m-1, n-1)
upper.bound <- c()
lower.bound <- c()

for (j in 1:10000){
  xlst <- rnorm(m, 3, 3)
  ylst <- rnorm(n, 5, 5)
  Sx2 <- var(xlst)
  Sy2 <- var(ylst)
  lower.bound[j] <- Sy2/Sx2 * f1
  upper.bound[j] <- Sy2/Sx2 * f2
}

plot(1:10000, rep(25/16,10000), xlab = "Experiment", ylab = "CI", main = "CI of the variance ratio", cex=0.1, ylim = c(0, 40))
points(1:10000, upper.bound, col = "red", cex = 0.1)
points(1:10000, lower.bound, col = "blue", cex = 0.1)
sum(upper.bound < 25/16 | lower.bound > 25/16)

而後是\(\mu_1,\mu_2\)均已知的狀況，這時方差的估計量會是

\[Q_x^2=\frac{1}{m}\sum_{j=1}^m(X_j-\mu_1)^2,\quad Q_y^2=\frac{1}{n}\sum_{j=1}^n(Y_j-\mu_2)^2, \]

顯然有

\[\frac{mQ_x^2}{\sigma_1^2}\sim \chi^2(m),\quad \frac{nQ_y^2}{\sigma_2^2}\sim\chi^2(n), \]

因而樞軸量爲

\[\frac{\frac{mQ_x^2}{\sigma_1^2}/m}{\frac{nQ_y^2}{\sigma_2^2}/n}=\frac{Q_x^2}{Q_y^2}\frac{\sigma_2^2}{\sigma_1^2}\sim F(m,n), \]

即\(\sigma_2^2/\sigma_1^2\)的區間估計爲

\[\left[\frac{Q_y^2}{Q_x^2}F_{1-\alpha/2}(m,n),\frac{Q_y^2}{Q_x^2}F_{\alpha/2}(m,n) \right]. \]

事實上，方差的估計，就是在\(\mu\)已知時使用對均值的離差平方和，在\(\mu\)未知時使用對樣本均值的離差平方和，注意分母與\(\chi^2\)分佈的自由度便可。由此，讀者應該也能推斷出\(\mu_1\)已知\(\mu_2\)未知，或者\(\mu_1\)未知\(\mu_2\)已知時方差比的區間估計了。

本文給出了區間估計的求法——樞軸量法的實際應用案例，並對正態整體的參數做區間估計。下篇文章將着眼於非正態整體的區間估計，因爲非正態整體不像正態整體同樣，擁有性質良好、分佈透明的點估計可使用，因此咱們可能須要尋找一種更通用的方法來求區間估計。