Processing math: 1%
跳到主要內容

Maximum Likelihood Estimation As Optimization Problem


在統計學裡, Maximum Likelihood Estimation (MLE) 是常見對於母體參數的估計,本文詳細介紹它數學的詳細大致長相,以及它是一個無約束最佳化問題(unconstrained optimization problem) ,所以可以用最佳化演算法去逼近 MLE
[Notation]
母體未知參數為 \theta = (\theta_1,\theta_2,...\theta_m) \in \Theta \subset \mathbb{R}^{m}
隨機向量為 X = (X_1,X_2,....X_p)
給定向量為 x = (x_1,x_2,...,x_p)
母體分布 (joint probability density function) f = f(x,\theta)
(1)\forall \theta \in \Theta, \int_{x\in \mathbb{R}^{p}} f(x,\theta) = 1 \qquad (2) \forall (x , \theta) \in \mathbb{R}^{p}\times \Theta \quad f(x,\theta) \in [0,1]

[樣本的聯合分配]
考慮獨立同分配 identical independent distributed (iid) 抽樣,n 個樣本為 \mathcal{D}=\{(X^k)\}^{n}_{k=1}d = \{(x^k)\}^{n}_{k=1},其中 X^{k} = (X^{k}_1,X^{k}_2,...X^{k}_p) \{(X^k)\}^{n}_{k=1} \sim  \prod^{n}_{k=1} f(x^{k},\theta)  =  \underbrace{f(x^{1},\theta) \cdot f(x^{2},\theta) ... \cdot  f(x^{k},\theta)}_{\text{number of variables } np + m} =(*)  

[Likelihood Function]
當樣本已知 \mathcal{D} = d 給定,視 \theta 為變數,則統計學術語稱(*)為 Likelihood function L(\theta|d) ,則我們可以找極大化這個函數,找到最大的\theta,即  \theta^{*}(d):=\underset{\theta \in \Theta}{argmax}L(\theta|d)
而當視為 \mathcal{D} (樣本尚未給定)時,我們稱之為Maximum Likelihood Estimation(MLE)  \theta_{MLE}(\mathcal{D}):=\underset{\theta \in \Theta}{argmax}L(\theta|\mathcal{D}) = \underset{\theta \in \Theta}{argmax}\prod^{n}_{k=1} f(X^{k},\theta)
其中 \theta_{MLE}(\mathcal{D}) = \theta_{MLE}(X^1,X^2,.....X^{n}) 是隨機變數

[計算 MLE 的值]
X^{k} = x^k (樣本觀察到的時候),我們要計算統計量相當於解一個最佳化問題即
\theta_{MLE}(\mathcal{D}=d) =  \underset{\theta \in \Theta}{argmax}\prod^{n}_{k=1} f(x^{k},\theta)
我們可以取對數 log 讓連乘變連加(p.s 因為 log 是嚴格遞增函數)
\underset{\theta \in \Theta}{argmax}\prod^{n}_{k=1} f(x^{k},\theta) = \underset{\theta \in \Theta}{argmax}\sum^{n}_{k=1} log f(x^{k},\theta) 
我們分別對 \theta_i 微分 =0 (first derivative test),寫成聯立方程組即
\forall i =1,2,....m  \quad \frac{\partial}{\partial \theta_i}\sum^{n}_{k=1} log f(x^{k},\theta) = \underbrace{\sum^{n}_{k=1}\frac{\frac{\partial}{\partial \theta_i}f(x^{k},\theta)}{f(x^{k},\theta)}}_{g_i(\theta)} = 0  
相當於 m 個 unknowns , m 個 constraints 的非線性方程組(nonlinear system)
\left\{\begin{array}{c} g_1(\theta_1,....\theta_m) = 0 \\ g_2(\theta_1,....\theta_m) = 0 \\ g_3(\theta_1,....\theta_m) = 0\\ g_4(\theta_1,....\theta_m) = 0\\ ..... \\ g_m(\theta_1,....\theta_m) = 0\\ \end{array}\right.
可以使用牛頓法(Newton Method)計算近似的 MLE 值 !!

+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
[Classical Example :  p = 1 , m =2 ]
我們用以上的例子證明 一維的常態分布(Normal Distribution) (\mu,\sigma^2) 的 MLE  是 (\bar{X},\frac{n-1}{n}S^2)  ,其中 \theta = (\mu , \sigma^2)\theta^{max}(\mathcal{D})= (\bar{X},\frac{n-1}{n}S^2)
其中 \bar{X} =\frac{\sum^{n}_{k=1}X^k}{n} (sample mean) ,
其中 S^2 = \frac{1}{n-1}\sum (X^k- \bar{X})^2 (unbiased sample variance)

考慮 Normal Distribution f(x^k,\mu,\sigma^2) = \frac{1}{\sqrt{2\pi\sigma^2}} \exp\left[{-\frac{(x^k-\mu)^2}{2\sigma^2}}\right]
計算 \frac{\partial}{\partial \theta_i}f(x^{k},\theta)
\left\{\begin{array}{l} \frac{\partial}{\partial \mu}f(x^k,\mu,\sigma^2)  = \frac{\partial}{\partial \mu}\left( \frac{1}{\sqrt{2\pi\sigma^2}} \exp\left[{-\frac{(x^k-\mu)^2}{2\sigma^2}}\right] \right) \\ = \frac{1}{\sqrt{2\pi \sigma^2}} e^{-\frac{(x-\mu)^2}{2\sigma^2}} \cdot \left(\frac{x^k-\mu}{\sigma^2}\right) = \left(\frac{x^k-\mu}{\sigma^2}\right) f(x^k,\mu,\sigma^2) \\  \frac{\partial}{\partial \sigma^2}f(x^k,\mu,\sigma^2)  = \frac{\partial}{\partial \sigma^2} \left( \frac{1}{\sqrt{2\pi\sigma^2}} \exp\left[{-\frac{(x^k-\mu)^2}{2\sigma^2}}\right]\right) \\ = \frac{1}{\sqrt{2\pi \sigma^2}} e^{-\frac{(x^k-\mu)^2}{2\sigma^2}}\left[ \frac{(x^k-\mu)^2}{\sigma^3}-\frac{1}{\sigma} \right] =\left[ \frac{(x^k-\mu)^2}{\sigma^3}-\frac{1}{\sigma} \right]   f(x^k,\mu,\sigma^2)  \\  \end{array}\right.

再來計算 First Derivative Test  \forall i \qquad \sum^{n}_{k=1}\frac{\frac{\partial}{\partial \theta_i}f(x^{k},\theta)}{f(x^{k},\theta)} = 0
會得到兩個式子 :
\begin{array}{l} (1)  \sum^{n}_{k=1}\left(\frac{x^k-\mu}{\sigma^2}\right)  = 0  \overset{\sigma > 0}{\Longrightarrow} \mu^{max}(d) = \frac{1}{n}{\sum_{k=1} x^k} \\ \Longrightarrow   \bar{X}(\mathcal{D}) := \frac{1}{n}\sum_{k=1} X^k\text{ 為 } \mu \text{ 的 MLE} \\  \\ (2) \sum^{n}_{k=1}\left[ \frac{(x^k-\mu)^2}{\sigma^3}-\frac{1}{\sigma} \right]  = 0   \Longrightarrow  {(\sigma^{max})}^2 = \sum^{n}_{k=1}(x^k-\mu)^2     \\   \underset{by (1)}{\Longrightarrow}  \frac{n-1}{n}S^2(\mathcal{D}):=   \frac{1}{n}\sum^{n}_{k=1}(X^k-\bar{X})^2 \text{ 為 } \sigma^2 \text{ 的 MLE}    \\  \end{array}

注意 : MLE of \sigma^2 就是 \frac{1}{n}(X^k-\bar{X})^2 但它是 biased estimator !!
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++





[以上純為學術經驗交流知識分享,如有錯誤或建議可留言~~] 
by Plus & Minus 2017.08


留言

這個網誌中的熱門文章

Nash Equilibrium & Best Responce Function (BRF) In Continuous Strategies

經濟學重要的賽局理論( Game Theory )領域,用數學描述人與人之間的理性互動,最重要的就是尋找奈許均衡( Nash equilibrium ), 本篇介紹其數學規劃與非線性方程組!!  假設有 p 名玩家(player i),i=1,2,3,4,5,....p , 正在玩一場遊戲(Game)~~,完全不合作,各自獨立作決策 每個人有決策向量 x_i \in \Omega_i \subseteq R^{n_i} (有n_i個決策變數)  定義長向量: \underbrace{x =  (x_1,x_2,x_3,....x_p)}_{\# \text{ of } \sum^{p}_{i=1}n_i \text{ variables }} \in  \prod^{p}_{i=1} \Omega_i = \Omega 對於每個 player i ,長向量可以寫成 x = (x_i , x_{-i})x_{-i} 代表其他人(不是 player i) 能做的決策向量。 所有人各自作決策後,每個人都會個自的存在報酬效用函數 f_i (x)  \in \mathbb{R}  (報酬函數皆為公開已知資訊) 假設每位玩家是理性人(會極大化自己效用) 即 \forall i = 1,2,3,4....p \qquad  \underset{x_i \in \Omega_i}{\text{max }}f_i(x)   [註: 如果為合作可視為多目標規劃問題( multiobjective ),即 x_1,x_2,...x_p 可以由領導人一起決定] [註: 如果為合作而且把效用加總,即目標式變成 \sum_{i=1}^{p} f_i(x) ,可能對集體效益有更大的幫助,但是如何分配效益給 ( player i )會是個議題,可以查關鍵字 fair optimization ] 我們可以定義每個 player i 的 Best Response Function (BRF) or Best Reponce Set S_i(x_{-i}) \subset \Omega_i $$  S_i(x...

Lattice & Multinomial Theorem

本文介紹格子點(Lattice) 幾何意義與多項式定理(Mutinomial Theorem) 的關係,並可協助我們理解計算一些機率問題。 [符號定義] 非負整數 / 非負實數:  \mathbb{Z}_{\geq 0} := \{0,1,2,3,4,......\}  \subseteq [0,\infty) =: \mathbb{R}_{\geq 0} 離散機率向量:  p_{I} := (p_{i})_{i \in I} \text{ s.t } \sum_{i\in I}p_i =1 ,|I|<\infty  發生事件 i \in I 的累積次數向量: k_{I} := (k_i)_{i \in I} \in \mathbb{Z}^{|I|}_{\geq 0} \mathbb{Z}^{|I|}_{\geq 0} 就是 |I| 維格子點 !! [格子點情境] 出發點定義為 k^{start}_{I}:= \overbrace{(0,0...,0)}^{|I|},今發生一次 p_{I} 分布隨機互斥事件,等價於"點的移動"(state transition),數學定義如下:   \text{Event } i  \text{ happens }  \Longleftrightarrow  \overbrace{(\color{red}{k_i},k_{-i})}^{k^{old}_{I}}  \underset{\text{with probability }p_{i}}{\longrightarrow}   \overbrace{(\color{red}{k_i+1},k_{-i})}^{ k^{new}_{I}}    PS1: 其中  k_{-i} := (k_{i'})_{i' \in I-\{i\}} PS2: 不管怎麼走都在第一象限,也就是只能往右,往上,往高.... 當發生 n 次獨立同分布 p_{I} (iid) 的事件後,所有可能點位置在以下的集合上 $$  S_{n}(\col...

Linear Regression By Using Linear Programming

當拿到一筆資料準備玩統計,往往會想要做線性迴歸( Linear Regression ),找出一個模型( mathematical model )來解釋變數間的關係,一般都是使用平方距離,但是如果我們採用絕對值距離呢?? 而剛好在工業工程( Industrial Engineering ),作業研究( Operation Research ) 領域,發展成熟的線性規劃( Linear Programming ) 恰好可以來解決,是一個跨領域的應用 !! 已經存在有許多商業或open source 軟體,如: Gurobi , Cplex , Xpress , Mosek , SCIP  可以輕易求解大型的線性規劃問題。而不僅如此也可以利用整數規劃( Integer Programming )來做特徵選擇 ( Feature Selection ),甚至可以偵測離群值( Detect Outlier ) !! 本文只介紹最小絕對值和,關於 Feature Selection , Detect Outlier 可以參考 Mixed-Integer Linear Programming Robust Regression with Feature Selection , Oleksii Omelchenko , 2010 的論文。 [Data Fitting Problem] 給定n筆實數型訓練資料 (training data) \{(x^{k},y^{k})\}^{n}_{k=1} = \mathcal{D} , x^{k} =(x^{k}_1,x^{k}_2, ... , x^{k}_{p})\in \mathbb{R}^{p} , y^{k} \in \mathbb{R} , 我們目標是想要找到一個函數 f_{\mathcal{D}} : \mathbb{R}^p \rightarrow \mathbb{R} 使得  \forall x \in \mathbb{R}^{p} , f_{\mathcal{D}}(x) \approx y , 精確來說: $$ \text{Find } f_{\mathcal{D}} \text{ such that } f_{\mathcal{D}}(x)\...