[問題] 信賴區間和覆蓋率的問題?

看板Statistics作者 (5566)時間7年前 (2017/02/21 02:10), 編輯推噓1(1049)
留言50則, 2人參與, 最新討論串1/1
假設我們要估計參數θ,找到一個估計量θ^, 並且用bootstrap的方式找出stand error 接著建立95%信賴區間C.I.,重複產生100組樣本, 並且算出100個θ^的值,建立100組C.I.,發現 100組裡卻只有80多組有包含θ(因為是模擬所以θ已知), 請問會造成這種現像的可能原因有哪些? 我自己想來想去只想到估計量θ^可能是有偏的估計量, 不知還有沒有其他因素,請教版上大大,謝謝 -- ※ 發信站: 批踢踢實業坊(ptt.cc), 來自: 36.224.17.7 ※ 文章網址: https://www.ptt.cc/bbs/Statistics/M.1487614226.A.7DB.html

02/21 11:20, , 1F
可否請你把抽樣和建立CI的過程說明清楚一點?我認真看
02/21 11:20, 1F

02/21 11:20, , 2F
了幾次還是不能確定你怎麼做的。
02/21 11:20, 2F

02/21 13:52, , 3F
首先先用模擬的方式產生一組我有興趣研究的樣本,用這組
02/21 13:52, 3F

02/21 13:54, , 4F
樣本求出θ^,因為θ^很複雜,所以沒辦法用手算求出Var,只
02/21 13:54, 4F

02/21 13:57, , 5F
好用bootstrap的方式求stand error代替,假設bootstrap
02/21 13:57, 5F

02/21 14:00, , 6F
1000次,再每bootstrap一次下也可求出θ^,假設叫θb^,共
02/21 14:00, 6F

02/21 14:02, , 7F
有1000個θb^,這1000個θb^的s.d.,即bootstrap下的
02/21 14:02, 7F

02/21 14:12, , 8F
stand error,接這建立C.I.,我建立C.I.的方式有四種,
02/21 14:12, 8F

02/21 14:17, , 9F
以下用R的語法說明較快,先說明一些代號,th0=θ^,
02/21 14:17, 9F

02/21 14:20, , 10F
th=那1000組的θb^,alpha<-c(0.025,0.975),
02/21 14:20, 10F

02/21 14:23, , 11F
1.normalC.I.:th0 + qnorm(alpha) * sd(th)
02/21 14:23, 11F

02/21 14:25, , 12F
2.basic C.I.:2*th0-quantile(th,rev(alpha),type=1)
02/21 14:25, 12F

02/21 14:27, , 13F
3.percentile C.I.:quantile(th, alpha, type=6)
02/21 14:27, 13F

02/21 14:30, , 14F
剛剛忘了說我bootstrap的方式是採抽後放回的抽樣方式
02/21 14:30, 14F

02/21 14:48, , 15F
第4種 BCa C.I.最麻煩,zalpha <- qnorm(alpha)
02/21 14:48, 15F

02/21 14:49, , 16F
z0<- qnorm(sum(th0 < th) / length(th))
02/21 14:49, 16F

02/21 14:49, , 17F
L<- mean(th.jack) - th.jack
02/21 14:49, 17F

02/21 14:49, , 18F
a <- sum(L^3)/(6 * sum(L^2)^1.5)
02/21 14:49, 18F

02/21 14:50, , 19F
adj.alpha <- pnorm(z0 + (z0+zalpha)/(1-a*(z0+zalpha)
02/21 14:50, 19F

02/21 14:50, , 20F
))
02/21 14:50, 20F

02/21 14:50, , 21F
quantile(th, adj.alpha, type=6)--->BCa C.I.
02/21 14:50, 21F

02/21 14:53, , 22F
這裡面的th.jack就是把我產生的那組資料
02/21 14:53, 22F

02/21 14:58, , 23F
假設叫T(n筆data),逐一的剔掉一筆data後,在剩餘的data下
02/21 14:58, 23F

02/21 15:00, , 24F
所算出的估計量,假設叫θj^,所以從第一筆開始剔除到最後
02/21 15:00, 24F

02/21 15:02, , 25F
一筆,會得到n筆的θj^,這n筆θj^就叫作th.jack
02/21 15:02, 25F

02/21 15:04, , 26F
更正:是在剩餘的data下所算出的"估計值",不是"估計量"
02/21 15:04, 26F

02/21 15:05, , 27F
---------------------------------------------------
02/21 15:05, 27F

02/21 15:06, , 28F
以上是產生一組data下所算出的4種C.I.,那我重複產生100
02/21 15:06, 28F

02/21 15:08, , 29F
組data下可算出100組4種C.I.,在分別去計算4種C.I.有沒有
02/21 15:08, 29F

02/21 15:10, , 30F
包含我所設定的參數θ的比率,即為文中的覆蓋率
02/21 15:10, 30F

02/21 15:13, , 31F
我比較疑惑的是我的4種C.I.都是建立在95%信心水準下
02/21 15:13, 31F

02/21 15:15, , 32F
照理說做100組應該要有95組左右的C.I都會包到設定的參數
02/21 15:15, 32F

02/21 15:16, , 33F
怎麼會掉到80多組?所以我想是不是估計量是有偏的,所以
02/21 15:16, 33F

02/21 15:17, , 34F
覆蓋率才會有落差?
02/21 15:17, 34F

02/21 15:28, , 35F
不好意思~我說了一大堆~我也不知道有沒有講到重點,實際
02/21 15:28, 35F

02/21 15:30, , 36F
上我的θ是一個向量,不過就只是每個分量的元素都做和我
02/21 15:30, 36F

02/21 15:31, , 37F
上面所說的一樣的事而已
02/21 15:31, 37F

02/21 15:32, , 38F
附上原始的codehttp://pastebin.com/jgUBxtvj
02/21 15:32, 38F

02/21 17:03, , 39F
我猜你boot.CI有問題:哪要那麼多額外的參數?
02/21 17:03, 39F

02/21 17:04, , 40F
一般來說只要一組n個數字,頂多再要求alpha=0.95就好了
02/21 17:04, 40F

02/21 17:05, , 41F
不管你要哪種CI,就只要這二個參數啊。為什麼要額外輸入
02/21 17:05, 41F

02/21 17:05, , 42F
th0 th stat 這些東西?
02/21 17:05, 42F

02/21 17:41, , 43F
推文裡前3種C.I.只要用到th0和th,th就是您說的一組數字
02/21 17:41, 43F

02/21 17:45, , 44F
x就是輸入產生的data,而stat是算th.jack會用到的,簡單來
02/21 17:45, 44F

02/21 17:46, , 45F
我想的求boot 95% CI就只要一個參數,就是資料本身,還
02/21 17:46, 45F

02/21 17:46, , 46F
x和stat是算BCa C.I.會用到的東西所以我要輸入
02/21 17:46, 46F

02/21 17:46, , 47F
是沒搞懂額外參數的作用是什麼…
02/21 17:46, 47F

02/21 17:47, , 48F
oh謝謝。
02/21 17:47, 48F

02/21 17:48, , 49F
boot.CI回傳的會是4種C.I.的值,我把4種C.I.的算法都寫
02/21 17:48, 49F

02/21 17:49, , 50F
在一boot.CI這個函數了
02/21 17:49, 50F
文章代碼(AID): #1Ogp4IVR (Statistics)