你能簡單直觀地解釋一下 IRLS 方法來找到 GLM 的 MLE 嗎?
背景:
我正在嘗試關注普林斯頓對 GLM 的 MLE 估計的評論。
我了解 MLE 估計的基礎知識:
likelihood
、score
、 觀察和預期Fisher information
以及Fisher scoring
技術。而且我知道如何用 MLE 估計證明簡單的線性回歸。
問題:
我什至無法理解這種方法的第一行:(
背後的直覺是什麼工作變量定義為:
為什麼使用它們而不是估計?
它們之間的關係是什麼
response/link function
?和如果有人有一個簡單的解釋,或者可以將我引導到一個更基本的文本,我將不勝感激。
幾年前,我為我的學生寫了一篇關於這個的論文(西班牙語),所以我可以嘗試在這裡重寫這些解釋。我將通過一系列增加複雜性的示例來研究 IRLS(迭代重加權最小二乘法)。對於第一個示例,我們需要位置規模族的概念。讓 $ f_0 $ 在某種意義上是一個以零為中心的密度函數。我們可以通過定義來構造一個密度族 $$ f(x)= f(x;\mu,\sigma)= \frac{1}{\sigma} f_0\left(\frac{x-\mu}{\sigma}\right) $$ 在哪裡 $ \sigma > 0 $ 是一個尺度參數並且 $ \mu $ 是位置參數。在測量誤差模型中,通常將誤差項建模為正態分佈,我們可以使用上面構建的位置尺度族來代替該正態分佈。什麼時候 $ f_0 $ 是標準正態分佈,上面的構造給出了 $ \text{N}(\mu, \sigma) $ 家庭。
現在我們將在一些簡單的例子中使用 IRLS。首先,我們將在模型中找到 ML(最大似然)估計量 $$ Y_1,Y_2,\ldots,Y_n \hspace{1em} \text{i.i.d} $$ 隨著密度 $$ f(y)= \frac{1}{\pi} \frac{1}{1+(y-\mu)^2},\hspace{1em} y\in{\mathbb R}, $$ 柯西分佈位置族 $ \mu $ (所以這是一個位置族)。但首先是一些符號。的加權最小二乘估計量 $ \mu $ 是(誰)給的 $$ \mu^{\ast} = \frac{\sum_{i=1}^n w_i y_i} {\sum_{i=1}^n w_i}. $$ 在哪裡 $ w_i $ 是一些權重。我們將看到 ML 估計器 $ \mu $ 可以用相同的形式表示,用 $ w_i $ 殘差的一些函數 $$ \epsilon_i = y_i-\hat{\mu}. $$ 似然函數由下式給出 $$ L(y;\mu)= \left(\frac{1}{\pi}\right)^n \prod_{i=1}^n \frac{1}{1+(y_i-\mu)^2} $$ 並且對數似然函數由下式給出 $$ l(y)= -n \log(\pi) - \sum_{i=1}^n \log\left(1+(y_i-\mu)^2\right). $$ 它的導數關於 $ \mu $ 是 $$ \begin{eqnarray} \frac{\partial l(y)}{\partial \mu}&=& 0-\sum \frac{\partial}{\partial \mu} \log\left(1+(y_i-\mu)^2\right) \nonumber \ &=& -\sum \frac{2(y_i-\mu)}{1+(y_i-\mu)^2}\cdot (-1) \nonumber \ &=& \sum \frac{2 \epsilon_i}{1+\epsilon_i^2} \nonumber \end{eqnarray} $$ 在哪裡 $ \epsilon_i=y_i-\mu $ . 寫 $ f_0(\epsilon)= \frac{1}{\pi} \frac{1}{1+\epsilon^2} $ 和 $ f_0'(\epsilon)=\frac{1}{\pi} \frac{-1\cdot 2 \epsilon}{(1+\epsilon^2)^2} $ ,我們得到 $$ \frac{f_0'(\epsilon)}{f_0(\epsilon)} = \frac{\frac{-1 \cdot2\epsilon}{(1+\epsilon^2)^2}} {\frac{1}{1+\epsilon^2}} = -\frac{2\epsilon}{1+\epsilon^2}. $$ 我們發現 $$ \begin{eqnarray} \frac {\partial l(y)} {\partial \mu} & =& -\sum \frac {f_0'(\epsilon_i)} {f_0(\epsilon_i)} \nonumber \ &=& -\sum \frac {f_0'(\epsilon_i)} {f_0(\epsilon_i)} \cdot \left(-\frac{1}{\epsilon_i}\right) \cdot (-\epsilon_i) \nonumber \ &=& \sum w_i \epsilon_i \nonumber \end{eqnarray} $$ 我們使用定義的地方 $$ w_i= \frac{f_0'(\epsilon_i)} {f_0(\epsilon_i)} \cdot \left(-\frac{1}{\epsilon_i}\right) = \frac{-2 \epsilon_i} {1+\epsilon_i^2} \cdot \left(-\frac{1}{\epsilon_i}\right) = \frac{2}{1+\epsilon_i^2}. $$ 記住那個 $ \epsilon_i=y_i-\mu $ 我們得到方程 $$ \sum w_i y_i = \mu \sum w_i, $$ 這是IRLS的估計方程。注意
- 權重 $ w_i $ 總是積極的。
- 如果殘差很大,我們對相應的觀察給予較少的權重。
為了在實踐中計算 ML 估計量,我們需要一個起始值 $ \hat{\mu}^{(0)} $ ,例如,我們可以使用中位數。使用這個值我們計算殘差 $$ \epsilon_i^{(0)} = y_i - \hat{\mu}^{(0)} $$ 和權重 $$ w_i^{(0)} = \frac{2}{1+\epsilon_i^{(0)} }. $$ 的新價值 $ \hat{\mu} $ 是(誰)給的 $$ \hat{\mu}^{(1)} = \frac{\sum w_i^{(0)} y_i} {\sum w_i^{(0)} }. $$ 以這種方式繼續我們定義 $$ \epsilon_i^{(j)} = y_i- \hat{\mu}^{(j)} $$和 $$ w_i^{(j)} = \frac{2}{1+\epsilon_i^{(j)} }. $$ 通過時的估計值 $ j+1 $ 算法變為 $$ \hat{\mu}^{(j+1)} = \frac{\sum w_i^{(j)} y_i} {\sum w_i^{(j)} }. $$ 繼續直到序列 $$ \hat{\mu}^{(0)}, \hat{\mu}^{(1)}, \ldots, \hat{\mu}^{(j)}, \ldots $$ 收斂。
現在我們用更一般的位置和規模族來研究這個過程, $ f(y)= \frac{1}{\sigma} f_0(\frac{y-\mu}{\sigma}) $ , 細節較少。讓 $ Y_1,Y_2,\ldots,Y_n $ 與上面的密度無關。也定義 $ \epsilon_i=\frac{y_i-\mu}{\sigma} $ . 對數似然函數是 $$ l(y)= -\frac{n}{2}\log(\sigma^2) + \sum \log(f_0\left(\frac{y_i-\mu}{\sigma}\right)). $$ 寫作 $ \nu=\sigma^2 $ , 注意 $$ \frac{\partial \epsilon_i}{\partial \mu} = -\frac{1}{\sigma} $$ 和 $$ \frac{\partial \epsilon_i}{\partial \nu} = (y_i-\mu)\left(\frac{1}{\sqrt{\nu}}\right)' = (y_i-\mu)\cdot \frac{-1}{2 \sigma^3}. $$ 計算對數似然導數 $$ \frac{\partial l(y)}{\partial \mu} = \sum \frac{f_0'(\epsilon_i)}{f_0(\epsilon_i)}\cdot \frac{\partial \epsilon_i}{\partial \mu} = \sum\frac{f_0'(\epsilon_i)}{f_0(\epsilon_i)}\cdot\left(-\frac{1}{\sigma}\right)= -\frac{1}{\sigma}\sum\frac{f_o'(\epsilon_i)}{f_0(\epsilon_i)}\cdot \left(-\frac{1}{\epsilon_i}\right)(-\epsilon_i) = \frac{1}{\sigma}\sum w_i \epsilon_i $$ 並將其等於零給出與第一個示例相同的估計方程。然後尋找估算器 $ \sigma^2 $ : $$ \begin{eqnarray} \frac{\partial l(y)}{\partial \nu} &=& -\frac{n}{2}\frac{1}{\nu} + \sum\frac{f_0'(\epsilon_i)}{f_0(\epsilon_i)}\cdot \frac{\partial \epsilon_i}{\partial\nu} \nonumber \ &=& -\frac{n}{2}\frac{1}{\nu}+\sum\frac{f_0'(\epsilon_i)}{f_0(\epsilon_i)} \cdot \left(-\frac{(y_i-\mu)}{2\sigma^3}\right) \nonumber \ &=& -\frac{n}{2}\frac{1}{\nu} - \frac{1}{2}\frac{1}{\sigma^2} \sum\frac{f_0'(\epsilon_i)}{f_0(\epsilon_i)}\cdot \epsilon_i\nonumber \ &=& -\frac{n}{2}\frac{1}{\nu}-\frac{1}{2}\frac{1}{\nu} \sum\frac{f_0'(\epsilon_i)}{f_0(\epsilon_i)}\cdot \left(-\frac{1}{\epsilon_i}\right) (-\epsilon_i)\cdot\epsilon_i\nonumber \ &=& -\frac{n}{2}\frac{1}{\nu}+\frac{1}{2}\frac{1}{\nu}\sum w_i \epsilon_i^2 \stackrel{!}{=} 0. \nonumber \end{eqnarray} $$ 導致估計器 $$ \hat{\sigma^2} = \frac{1}{n}\sum w_i (y_i-\hat{\mu})^2. $$ 在這種情況下也可以使用上面的迭代算法。
在下文中,我們給出了一個使用 R 的數值示例,用於雙指數模型(具有已知比例)和數據
y <- c(-5,-1,0,1,5)
。對於該數據,ML 估計器的真實值為 0。初始值為mu <- 0.5
。該算法的一次通過是iterest <- function(y, mu) { w <- 1/abs(y - mu) weighted.mean(y, w) }
使用此功能,您可以嘗試“手動”進行迭代然後迭代算法可以通過
mu_0 <- 0.5 repeat {mu <- iterest(y, mu_0) if (abs(mu_0 - mu) < 0.000001) break mu_0 <- mu }
練習:如果模型是 $ t_k $ 具有尺度參數的分佈 $ \sigma $ 顯示迭代由權重給出 $$ w_i = \frac{k + 1}{k + \epsilon_i^2}. $$ 練習:如果密度是邏輯的,則顯示權重由下式給出 $$ w(\epsilon) = \frac{ 1-e^\epsilon}{1+e^\epsilon} \cdot - \frac{1}{\epsilon}. $$
暫時我會把它留在這裡,我將繼續這篇文章。