一文讀懂最大似然估計

作者：阿尼·辛格

翻譯： 陳之炎

校對：丁楠雅

本文約4200字，建議閱讀10+分鐘。

本文將研究MLE是如何工作的，以及它如何用於確定具有任何分布的模型的係數。

簡介

解釋模型如何工作是數據科學中最為基本最為關鍵的問題之一。當你建立了一個模型之後，它給了你預期的結果，但是它背後的過程是什麼呢？作為一個數據科學家，你需要對這個經常被問到的問題做出解答。

例如，假設您建立了一個預測公司股價的模型。您注意到夜深人靜的時候，股票價格上漲得很快。背後可能有多種原因，找出可能性最大的原因便是最大似然估計的意義所在。這一概念常被用於經濟學、MRIs、衛星成像等領域。

來源：YouTube

在這篇文章中，我們將研究最大似然估計(以下簡稱MLE)是如何工作的，以及它如何用於確定具有任何分布的模型的係數。理解MLE將涉及到概率和數學，但我將嘗試通過例子使它更通俗易懂。

註：如前所述，本文假設您已經了解概率論的基本知識。您可以通過閱讀這篇文章來澄清一些基本概念：

《每個數據科學專業人員都應該知道的概率分布常識》：

https://www.analyticsvidhya.com/blog/2017/09/6-probability-distributions-data-science/

目錄

為什麼要使用最大似然估計(MLE)？

通過一個實例了解MLE

進一步了解技術細節

分布參數

似然

對數似然

最大似然估計

利用MLE確定模型係數

R語言的MLE實現

為什麼要使用最大似然估計(MLE)？

假設我們想預測活動門票的銷售情況。數據的直方圖和密度如下。

你將如何為這個變數建模？該變數不是正態分布的，而且是不對稱的，因此不符合線性回歸的假設。一種常用的方法是對變數進行對數、平方根（sqrt）、倒數等轉換，使轉換後的變數服從正態分布，並進行線性回歸建模。

讓我們試試這些轉換，看看結果如何：

對數轉換：

平方根轉換：

倒數轉換：

所有這些都不接近正態分布，那麼我們應該如何對這些數據進行建模，才能不違背模型的基本假設？如何利用正態分布以外的其他分布來建模這些數據呢？如果我們使用了不同的分布，又將如何來估計係數？

這便是最大似然估計(MLE)的主要優勢。

舉一個例子來加深對MLE的理解

在研究統計和概率時，你肯定遇到過諸如x＞100的概率，因為x服從正態分布，平均值為50，標準差為10。在這些問題中，我們已經知道分布（在這種情況下是正態分布）及其參數（均值和標準差），但在實際生活問題中，這些參數是未知的，並且必須從數據中估計出來。MLE可以幫助我們確定給定數據的分布參數。

讓我們用一個例子來加深理解：假設我們用數據來表示班級中學生的體重(以kg為單位)。數據如下圖所示(還提供了用於生成數據圖的R代碼)：

圖 1

ggplot(x, aes(x = x)) + geom_dotplot()

這似乎遵循正態分布。但是我們如何得到這個分布的均值和標準差呢？一種方法是直接計算給定數據的平均值和標準差，分別為49.8公斤和11.37公斤。這些值能很好地表示給定的數據，但還不能最好地描述總體情況。

我們可以使用MLE來獲得更穩健的參數估計。因此，MLE可以定義為從樣本數據中估計總體參數(如均值和方差、泊松率(Lambda)等)的方法，從而使獲得觀測數據的概率(可能性)最大化。

為了加深對MLE的理解，嘗試猜測下列哪一項會使觀察上述數據的概率最大化？

1. 均值=100，標準差=10

2. 均值=50 ，標準差=10

顯然，如果均值為100，我們就不太可能觀察到上述數據分布圖形。

進一步了解技術細節

知道MLE能做什麼之後，我們就可以深入了解什麼是真正的似然估計，以及如何對它最大化。首先，我們從快速回顧分布參數開始。

分布參數

首先，來了解一下分布參數。維基百科對這個詞的定義如下：「它是一個概率分布的量化指數」，可以將它視為樣本總數的數值特徵或一個統計模型。通過下面的圖表來理解它：

圖 2

鐘形曲線的寬度和高度的兩個參數決定均值和方差。這就是正態分布的分布參數。同樣，泊松分布由一個參數lambda控制，即事件在時間或空間間隔內發生的次數。

圖 3

大多數分布都有一個或兩個參數，但有些分布可以有多達4個參數，比如4參數β分布。

似然

從圖2和圖3中我們可以看到，給定一組分布參數，一些數據值比其他數據的概率更大。從圖1中，我們已經看到，當平均值是50而不是100時，給定的數據更有可能發生。然而，在現實中，我們已經觀察到了這些數據。因此，我們面臨著一個逆向問題：給定觀測數據和一個感興趣的模型，我們需要在所有概率密度中找到一個最有可能產生數據的概率密度函數/概率質量函數(f(x_θ)。

為解決這一逆向問題，我們通過逆轉f(x=θ)中數據向量x和(分布)參數向量θ來定義似然函數，即：

L(θ;x) = f(x| θ)

在MLE中，可以假定我們有一個似然函數L(θ；x)，其中θ是分布參數向量，x是觀測集。我們感興趣的是尋找具有給定觀測值(x值)的最大可能性的θ值。

對數似然

如果假設觀測集(Xi)是獨立的同分布隨機變數，概率分布為f0(其中f0=正態分布，例如圖1)，那麼手頭的數學問題就變得簡單了。似然函數可以簡化為：

為了求這個函數的極大值/極小值，我們可以取此函數w.r.tθ的導數，並將其設為0(因為零斜率表示最大值或極小值)。因為我們這裡有乘積，所以需要應用鏈式規則，這對乘積來說是相當麻煩的。一個聰明的訣竅是取似然函數的對數，並使其最大化。這將乘積轉換為加法，並且由於對數是一個嚴格遞增的函數，因此不會影響θ的結果值。所以我們有：

最大化似然

為找到對數似然函數LL（Th；x）的極大值，我們可以：

取LL(θ；x)函數w.r.tθ的一階導數，並將其等價於0；

取LL(θ；x)函數w.r.tθ的二階導數，並確認其為負值。

在許多情況下，微積分對最大化似然估計沒有直接幫助，但最大值仍然可以很容易地識別出來。在尋找最大對數似然值的參數值時，沒有任何東西比一階導數等於零具有更為 「優先」或特殊的位置。當需要估計一些參數時，它僅僅是一個方便的工具而已。

在通常情況下，能對函數求最大值的argmax的方法都可能適合於尋找log似然函數的最大值。這是一個無約束的非線性優化問題。我們尋求一種優化演算法以下列方式工作：

從任意起始點可靠地收斂到局部最小化

速度儘可能快

使用優化技術來最大化似然是非常常見的，可以有多種方法來實現(比如：牛頓法、Fisher評分法、各種基於共軛梯度的方法、最陡下降法、Nder-Mead型(單純形)方法、BFGS法和多種其他技術)。

結果表明，當模型假設為高斯時，MLE估計等價於一般的最小二乘法。

你可以參考下面這篇文章來證明它：

link：

http://people.math.gatech.edu/~ecroot/3225/maximum_likelihood.pdf

用MLE確定模型係數

現在讓我們看看如何使用MLE來確定預測模型的係數。

假設我們有一個n個觀測量y1，y2，…，yn的樣本，它們可以被看作是獨立的泊松隨機變數：Yi ～ P(μi)。此外，假設我們希望讓均值i(方差也同樣！)依賴於變數xi組成的向量。我們可以構成如下簡單的線性模型：

中θ是模型係數的向量。這個模型的缺點是右邊的線性預測器可以假定為任何實際值，而表示期望計數的左側泊松均值必須是非負的。這個問題的一個簡單的解決辦法是用線性模型來模擬平均值的對數。因此，我們考慮了一個具有對數鏈接對數的廣義線性模型，它可以寫成如下所示：

我們的目的是利用MLE找到θ。

現在，泊松分布如下：

利用在上一節中學到的對數似然概念來求θ。取上述方程的對數，忽略包含log(y！)的常數，我們得到的對數似然函數是：

其中，μi依賴協變數xi和θ係數的向量。我們可以用變元μi = exp(xi』θ)代替，求解方程，得到最大似然值θ。得到了θ向量之後，我們就可以通過乘以xi和θ向量來預測平均值的期望值。

用R語言實現MLE

在本節中，我們將採用一個真實的數據集，運用前面學到的概念，來解決一個問題。您可以從此鏈接下載數據集：

https://s3-ap-south-1.amazonaws.com/av-blog-media/wp-content/uploads/2018/07/ Train_Tickets.csv

數據集中的示例如下：

售票日期 計票

25-08-2012 00:00 8

25-08-2012 01:00 2

25-08-2012 02:00 6

25-08-2012 03:00 2

25-08-2012 04:00 2

25-08-2012 05:00 2

它有從2012年8月25日到2014年9月25日每小時售出的門票數量(約18K記錄)。我們的目的是預測每小時售出的門票數量。這是本文第一節中討論的同一個數據集。

這個問題可以用回歸、時間序列等技術來解決。在這裡，我們將使用我們已經學習到的統計建模技術，用R語言實現。

首先，分析一下數據。在統計建模中，我們更關心的是目標變數是如何分布的。讓我們看一看計數的分布：

hist(Y$Count, breaks = 50,probability = T ,main = "Histogram of Count Variable")

lines(density(Y$Count), col="red", lwd=2)

這可以看作是泊松分布，或者我們甚至可以嘗試擬合指數分布。

由於手頭的變數是票的計數，泊松分布是一個更合適的模型。指數分布通常用於模擬事件之間的時間間隔。

讓我們來計算一下這兩年售出的票的數量：

看上去隨著時間的推移，門票的銷售有了很大的增長。為了將問題簡單化，讓我們僅將時間作為一個因子來建模，其中時間定義為2012年8月25日以來幾個星期。我們可以把它寫成：

其中，μ（售出的票數）是泊松分布的平均值，而θ0和θ1是我們需要估計的係數。

組合方程1和2，我們得到如下的對數似然函數：

我們可以使用Rstats 4包中的mle() 函數來估計係數θ0和θ1。它需要下列主要參數：

需要最小化的負似然函數：這與我們剛剛導出的函數相同，但前面有一個負號[最大對數似然度等於最小化負對數似然]。

係數向量的起點：這是係數的初始預測值。結果可以根據這些值而變化，因為函數可能達到局部最小值。因此，通過運行不同起點的函數來驗證結果是很好的辦法。

BFGS是默認的對似然函數進行優化的方法。

在我們的示例中，負對數似然函數的編碼如下：

nll

x

y

mu = exp(theta0 + x*theta1)

-sum(y*(log(mu)) - mu)

}

我將數據分為訓練集和測試集，以便客觀地評價模型的性能。idx是測試集中行的指數。

set.seed(200)

idx

接下來，調用mle函數來獲得參數：

est

summary(est)

Maximum likelihood estimation

Call:

stats4::mle(minuslogl = nll, start = list(theta0 = 2, theta1 = 0))

Coefficients:

Estimate Std. Error

我們得到了係數的估計值，利用RMSE作為獲得測試集結果的評估度量：

pred.ts

rmse(pred.ts, Y$Count[idx])

86.95227

現在，讓我們看看我們的模型和標準線性模型(正態分布的誤差)的對比，本模型是用對數計數建模。

lm.fit

Coefficients:

Estimate Std. Error t value Pr(>|t|)

(Intercept) 1.9112992 0.0110972 172.2

age 0.0414107 0.0001768 234.3

pred.lm

rmse(exp(pred.lm), Y$Count[idx])

93.77393

可以看出來，標準線性模型的RMSE比我們的泊松分布模型要高。讓我們比較這兩個模型在樣本上的殘差圖，看看這些模型在不同的區域中表現如何：

與常規線性回歸相比，泊松回歸的誤差更接近於零。

在R語言中，使用stats包中的glm 函數建模更加簡單。它支持泊松，伽瑪，二項分布，Quasi，逆高斯，擬二項分布，擬泊松分布等等。對於上面所示的示例，可以使用以下命令直接獲得係數：

glm(Count ~ age, family = "poisson", data = Y)

Coefficients:

Estimate Std. Error z value Pr(>|z|)

(Intercept) 2.669 2.218e-03 1203

age 0.03278 2.612e-05 1255

在Python中也可以使用pymc.glm()函數，並設置為pm.glm.Familes.Poisson()系列。

章節附註

對上述例子的一種思考是，參數空間中是否存在比標準線性模型估計更好的係數。正態分布是預設分布，也是最廣泛使用的分布形式，但如果採用其它更為正確的分布，則可以得到更好的結果。最大似然估計是一種可以用於估計分布參數而不考慮所使用的分布的技術。因此，下次當您手頭有建模問題時，首先看看數據的分布情況，看看有沒有比正態分布更有意義的分布！

詳細的代碼和數據在我的Gizub存儲庫中可以找到。

有關使用年齡變數的數據讀取、格式化和建模的示例，請參閱「ModelingSingleVariables.R」文件。此外，我還使用了多個變數進行建模，它保存於「ModelingMultipleVariables.R」文件中。

原文標題：

An Introductory Guide to Maximum Likelihood Estimation (with a case study in R)

https://www.analyticsvidhya.com/blog/2018/07/introductory-guide-maximum-likelihood-estimation-case-study-r/

譯者簡介

陳之炎，北京交通大學通信與控制工程專業畢業，獲得工學碩士學位，歷任長城計算機軟體與系統公司工程師，大唐微電子公司工程師，現任北京吾譯超群科技有限公司技術支持。目前從事智能化翻譯教學系統的運營和維護，在人工智慧深度學習和自然語言處理（NLP）方面積累有一定的經驗。業餘時間喜愛翻譯創作，翻譯作品主要有：IEC-ISO 7816、伊拉克石油工程項目、新財稅主義宣言等等，其中中譯英作品「新財稅主義宣言」在GLOBAL TIMES正式發表。能夠利用業餘時間加入到THU 數據派平台的翻譯志願者小組，希望能和大家一起交流分享，共同進步

喜歡這篇文章嗎？立刻分享出去讓更多人知道吧！