大數(shù)據(jù)分析R中泊松回歸模型實例

　　如果您知道如何以及何時使用泊松回歸，它可能是一個非常有用的工具。在大數(shù)據(jù)分析R中泊松回歸模型實例中，我們將深入研究泊松回歸，它是什么以及R程序員如何在現(xiàn)實世界中使用它。
 

　　具體來說，我們將介紹：
 

　　1)泊松回歸實際上是什么，什么時候應該使用它

　　2)泊松分布及其與正態(tài)分布的區(qū)別

　　3)使用GLM進行泊松回歸建模4)

　　5)為計數(shù)數(shù)據(jù)建模泊松回歸

　　6)使用jtools可視化來自模型的發(fā)現(xiàn)

　　7)為費率數(shù)據(jù)建模泊松回歸
 

　　什么是泊松回歸模型?
 

　　泊松回歸模型最適合用于以結果為計數(shù)的事件建模?；蚋唧w地講， 計數(shù)數(shù)據(jù)：具有非負整數(shù)值的離散數(shù)據(jù)對某事進行計數(shù)，例如在給定時間范圍內事件發(fā)生的次數(shù)或雜貨店排隊的人數(shù)。
 

　　計數(shù)數(shù)據(jù) 也可以表示為 比率數(shù)據(jù)，因為事件在時間范圍內發(fā)生的次數(shù)可以表示為原始計數(shù)(即“一天，我們吃三頓飯”)或比率(“我們以每小時0.125餐”)。
 

　　泊松回歸可以讓我們確定哪些解釋變量(X值)對給定的響應變量(Y值，計數(shù)或比率)產(chǎn)生影響，從而幫助我們分析計數(shù)數(shù)據(jù)和比率數(shù)據(jù)。例如，雜貨店可以使用泊松回歸來更好地了解和預測生產(chǎn)線中的人數(shù)。
 

　　什么是泊松分布?
 

　　泊松分布是一種以法國數(shù)學家SiméonDenis Poisson命名的統(tǒng)計理論。它假設在特定時間范圍內發(fā)生事件y的概率 ，假設 y的 發(fā)生不受y先前發(fā)生時間的影響 。可以使用以下公式以數(shù)學方式表示：

　　此處， μ (在某些教科書中，您可能會看到 λ 而不是 μ)是每單位曝光可能發(fā)生的事件的平均次數(shù) 。也稱為泊松分布參數(shù)。該曝光可以是時間，空間，人口尺寸，距離或面積，但它往往是時間，與表示 噸。如果未給出曝??光值，則假定它等于 1。
 

　　讓我們通過為不同的μ值創(chuàng)建泊松分布圖來形象化這一點。
 

　　首先，我們將創(chuàng)建6種顏色??的向量：

　　接下來，我們將為分布創(chuàng)建一個列表，該列表將具有不同的μ值：

　　然后，我們將為μ創(chuàng)建一個值向量， 并在 每個μ的分位數(shù)范圍為0-20 的值之間循環(huán) ，并將結果存儲在列表中：

　　最后，我們將使用繪制點 plot()。 plot() 是R中的基本圖形功能。在R中繪制數(shù)據(jù)的另一種常見方法是使用流行的 ggplot2 軟件包。這在AAA教育的R課程中有所介紹 。但是對于大數(shù)據(jù)分析R中泊松回歸模型實例，我們將堅持使用基本R函數(shù)。

　　請注意，我們曾經(jīng)dpois(sequence,lambda) 在Poisson分布中繪制概率密度函數(shù)(PDF)。在概率論中，概率密度函數(shù)是描述連續(xù)隨機變量(其可能值為隨機事件的連續(xù)結果的變量)具有給定值的相對可能性的函數(shù)。(在統(tǒng)計中，“隨機”變量只是其結果是隨機事件的結果的變量。)
 

　　泊松分布與正態(tài)分布有何不同?
 

　　泊松分布最常用于查找在給定時間間隔內發(fā)生事件的概率。由于我們是在討論計數(shù)，對于Poisson分布，結果必須為0或更高–事件不可能發(fā)生負次數(shù)。另一方面，正態(tài)分布 是連續(xù)變量的連續(xù)分布，它可能導致正值或負值：

　　我們可以在R中生成正態(tài)分布，如下所示：
 

　　在R中，dnorm(sequence, mean, std.dev) 用于繪制正態(tài)分布的概率密度函數(shù)(PDF)。
 

　　要了解泊松分布，請從Chi Yau的R Tutorial教科書中考慮以下問題 ：
 

　　如果平均每分鐘有12輛汽車過橋，那么在任何給定的分鐘內有十七輛或更多的汽車過橋的概率是多少?
 

　　在這里，每分鐘過橋的平均汽車數(shù)量為 μ = 12。
 

　　ppois(q, u, lower.tail = TRUE) 是R函數(shù)，可提供隨機變量小于或等于某個值的可能性。
 

　　我們必須找到擁有十七 輛或更多汽車的可能性 ，因此我們將 lower.trail = FALSE q設置為16：

　　要獲得百分比，我們只需要將此輸出乘以100。現(xiàn)在我們有了問題的答案：在任何特定分鐘內，有17輛或更多的汽車過橋的可能性為10.1%。
 

　　泊松回歸模型和GLM
 

　　廣義線性模型是其中響應變量遵循除正態(tài)分布以外的其他分布的模型。這與線性回歸模型相反，在線性回歸模型中，響應變量遵循正態(tài)分布。這是因為廣義線性模型具有分類的響應變量，例如“是”，“否”;或A組，B組，因此，范圍從-∞到+∞。因此，響應變量和預測變量之間的關系可能不是線性的。在GLM中：
 

　　ÿ 我 ?=? α ?+? β 1 X 1 我 ?+? β 2 X 2 我 ?+ ... + β p X p我 ?+? Ê 我???????????????我 ?= 1,2，...。ñ
 

　　響應變量 y i 由預測變量 和某些誤差項的線性函數(shù)建模。
 

　　泊松回歸模型是 廣義線性模型(GLM) ，用于對計數(shù)數(shù)據(jù)和列聯(lián)表建模。輸出 Y (計數(shù))是遵循泊松分布的值。它假設 可以通過一些未知參數(shù)將其建模為線性形式的期望值(均值)的對數(shù) 。
 

　　注意： 在統(tǒng)計信息中，列聯(lián)表(示例)是取決于多個變量的頻率矩陣。
 

　　要將非線性關系轉換為線性形式，可以使用鏈接函數(shù)，該 函數(shù) 是 泊松回歸的對 數(shù)。因此，泊松回歸模型也稱為對數(shù)線性模型。泊松回歸模型的一般數(shù)學形式為：
 

　　升 Ô 克(Ý)= α ?+? β 1 X 1 ?+? β 2 X 2 ?+ ... + β p X p
 

　　哪里，
 

　　1)y：是響應變量

　　2)α 和 β：是數(shù)值系數(shù)， α 為截距，有時 α 也被表示為 β 0，它是相同的

　　3)x 是預測變量/解釋變量
 

　　使用諸如最大似然估計(MLE)或最大擬似然之類的方法來計算系數(shù)。
 

　　考慮具有一個預測變量和一個響應變量的方程：
 

　　l o g(y)= α ?+? β(x)
 

　　這等效于：
 

　　Ý ?=? È (α ?+? β(X)) ?=? È α ?+? Ê β ?*? X
 

　　注意：在Poisson回歸模型中，預測變量或解釋變量可以同時包含數(shù)值或分類值。
 

　　等分散是泊松分布和泊松回歸最重要的特征 之一，這意味著分布的均值和方差相等。
 

　　方差度量數(shù)據(jù)的傳播。它是“與均值平方差的平均值”。如果所有值都相同，則方差(Var)等于0。值之間的差異越大，方差越大。平均值是數(shù)據(jù)集的平均值。平均值是值的總和除以值的數(shù)量。
 

　　假設平均值(μ)表示為 E(X)
 

　　E(X)= μ
 

　　對于泊松回歸，均值和方差如下：
 

　　v一個[R (X)= σ 2 ë(X)
 

　　其中 σ 2 是色散參數(shù)。由于 v一個[R (X)= È(X)(=方差平均值)必須持有泊松模型是完全適合， σ 2 必須等于1。
 

　　當方差大于均值時，稱為 過度分散 ，并且大于1。如果小于1，則稱為 欠分散。
 

　　使用計數(shù)數(shù)據(jù)進行泊松回歸建模
 

　　在R中，該 glm() 命令用于對廣義線性模型進行建模。這是以下內容的一般結構 glm()：

　　在大數(shù)據(jù)分析R中泊松回歸模型實例中，我們將使用這三個參數(shù)。有關更多詳細信息，請查閱 R文檔，但讓我們快速看一下它們各自指的是什么：

　　glm() 通過以下默認鏈接功能為家庭提供八個選擇：

　　讓我們開始建模吧!
 

　　我們將對與編織過程中紗線斷裂的頻率有關的泊松回歸建模。這些數(shù)據(jù)可以datasets 在R 的包中找到 ，因此我們要做的第一件事是使用安裝包install.package("datasets") 并使用 加載庫 library(datasets)：

　　該 datasets 程序包包含大量數(shù)據(jù)集，因此我們需要專門選擇紗線數(shù)據(jù)。查閱軟件包文檔，我們可以看到它稱為 warpbreaks，因此讓我們將其存儲為對象。

　　讓我們看一下數(shù)據(jù)：

　　輸出： [1] "breaks" "wool" "tension"
 

　　我們的數(shù)據(jù)中有什么?
 

　　該數(shù)據(jù)集查看每根固定長度的織機上不同類型的織機發(fā)生了多少經(jīng)紗斷裂。我們可以在此處的文檔中閱讀有關此數(shù)據(jù)集的更多詳細信息，但是這是我們將要查看的三列以及每列所指的內容：

　　在這六種類型的經(jīng)紗中，有9種織機上都有測量值，數(shù)據(jù)集中總共有54個條目。
 

　　讓我們看一下如何使用ls.str() 命令來構造數(shù)據(jù)：

　　輸出：

　　從上面可以看到數(shù)據(jù)中存在的類型和級別。 閱讀 此書可進一步了解R中的因素。
 

　　現(xiàn)在，我們將使用 data 數(shù)據(jù)框。記住，使用泊松分布模型，我們試圖找出一些預測變量如何影響響??應變量。在這里， breaks 是響應變量 wool 和 tension 是預測變量。
 

　　我們可以breaks 通過創(chuàng)建直方圖來查看因變量數(shù)據(jù)的連續(xù)性：

　　顯然，數(shù)據(jù)不是像正態(tài)分布那樣呈鐘形曲線形式。
 

　　讓我們來檢查 mean() 和 var() 的因變量：

　　輸出： [1] 28.14815

　　輸出： [1] 174.2041
 

　　方差遠大于平均值，這表明我們將在模型中出現(xiàn)過度分散現(xiàn)象。
 

　　讓我們使用glm() 命令擬合泊松模型 。

　　summary() 是用于生成各種模型擬合函數(shù)結果的結果摘要的通用函數(shù)。

　　輸出：

　　解釋泊松模型
 

　　我們剛剛獲得了很多信息，現(xiàn)在我們需要對其進行解釋。命名的第一列 Estimate 是的系數(shù)值 α (截距)， β 1 等。以下是參數(shù)估計的解釋：
 

　　1)當x = 0時，e x p(α)=對均值μ的影響

　　2)È X p(β)=在X每單位增加，預測器變量具有倍增效應 ë X p(β上Y的平均值)，即 μ

　　3)如果 β = 0，則exp(β)= 1，并且預期計數(shù)為 e x p(α)，并且Y和X不相關。

　　4)如果 β > 0，則exp(β)> 1，并且預期計數(shù)是exp(β)的乘積是X = 0時的乘積

　　5)如果 β <0，則exp(β)<1，并且預期計數(shù)是exp(β)的乘積，比X = 0時小。
 

　　如果 family = poisson 保留， glm() 則使用最大似然估計MLE計算這些參數(shù) 。
 

　　R將類別變量視為偽變量。通過為變量中的級別分配一些數(shù)字表示形式，類別變量(也稱為指標變量)被轉換為虛擬變量。一般規(guī)則是，如果 因子變量中有 k個類別，則輸出 glm() 將有 k -1個類別，其余1作為基本類別。
 

　　我們可以在上面的摘要中看到，對于羊毛，“ A”已成為基礎，未在摘要中顯示。同樣，對于張力“ L”已成為基本類別。
 

　　要查看哪些解釋變量對響應變量有影響，我們將查看 p 值。如果 p小于0.05， 則該變量會影響響應變量。在上面的總結中，我們可以看到所有的p值均小于0.05，因此， 兩個 解釋變量(羊毛和張力)均對休息產(chǎn)生重大影響。注意*** 在每個變量的末尾如何使用R輸出 。星星數(shù)表示重要性。
 

　　在開始解釋結果之前，讓我們檢查模型是過度分散還是分散不足。如果“ 殘余偏差” 大于自由度，則存在過度分散。這意味著估算是正確的，但是標準誤差(標準偏差)是錯誤的，模型無法解釋。
 

　　Null偏差顯示僅包含截距(均值)的模型對響應變量的預測程度，而包含自變量的模型對響應變量的預測程度。上面我們可以看到，添加了3個(53-50 = 3)自變量，將偏差從297.37降低到210.39。值之間的更大差異意味著不合適。
 

　　因此，要獲得更正確的標準誤差，我們可以使用 擬泊松 模型：

　　輸出：

　　比較模型：


　　現(xiàn)在我們有了兩個不同的模型，讓我們比較一下它們，看看哪個更好。首先，我們將安裝該arm 庫，因為該 庫包含我們需要的功能：

　　現(xiàn)在，我們將使用該 se.coef() 函數(shù)從每個模型中提取系數(shù)，然后 cbind() 將提取的值組合到單個數(shù)據(jù)幀中，以便進行比較。

　　輸出：

　　在上面的輸出中，我們可以看到系數(shù)相同，但標準誤差不同。
 

　　牢記這些要點，讓我們來看看羊毛的估算 。它的值是 -0.2059884，并且指數(shù) -0.2059884 是 0.8138425。

　　輸出：[1] 0.1861575
 

　　這表明從A型羊毛更改為B型羊毛會導致 斷點 減少0.8138425 倍，因為估計值-0.2059884為負。換句話說，如果我們將羊毛類型從A更改為B，則假設所有其他變量都相同，則折斷次數(shù)將減少18.6%。
 

　　從模型預測
 

　　建立模型后，我們可以使用 predict(model, data, type) 包含訓練數(shù)據(jù)以外的其他數(shù)據(jù)的新數(shù)據(jù)框來預測結果。讓我們來看一個例子。
 

　　輸出： [1] 23.68056
 

　　我們的模型預測， B型羊毛和M級張力的羊毛大約會斷裂24次。
 

　　使用jtools可視化結果
 

　　與他人共享分析時，表格通常不是吸引人們注意力的最佳方法。圖表可以幫助人們更快地掌握您的發(fā)現(xiàn)。在R中可視化數(shù)據(jù)的最流行的方法可能是 ggplot2 (在AAA教育的數(shù)據(jù)可視化課程中講授 )，我們還將使用稱為jtools的出色 R軟件包 ，其中包括專門用于匯總和可視化回歸模型的工具。讓我們 jtools 來可視化 poisson.model2。

　　jtools 提供 plot_summs() 并 plot_coefs() 可視化了模型的摘要，還允許我們與比較模型 ggplot2。

　　在上面的代碼中，該plot_summs(poisson.model2, scale = TRUE, exp = TRUE) 圖使用中的擬泊松族繪制了第二個模型 glm。
 

　　1)中的第一個參數(shù) plot_summs() 是要使用的回歸模型，它可以是一個或多個。
 

　　2)scale 有助于解決變量比例不同的問題。
 

　　3)exp 之所以設置為，TRUE是因為對于Poisson回歸，我們更可能對估計的指數(shù)值感興趣，而不是線性的。
 

　　你可以找到更多的細節(jié)jtools和 plot_summs() 這里的文件中。
 

　　我們還可以可視化預測變量之間的交互。 jtools 為不同類型的變量提供不同的功能。例如，如果所有變量都是分類變量，則可以 cat_plot() 用來更好地理解它們之間的交互。對于連續(xù)變量， interact_plot() 使用。
 

　　在 warpbreaks 數(shù)據(jù)中，我們具有分類的預測變量，因此，我們將使用cat_plot() 可視化變量 之間的交互作用，為其提供參數(shù)以指定我們要使用的模型，正在查看的預測變量以及其他預測變量與之結合以產(chǎn)生結果。

　　輸出：

　　我們可以做同樣的事情來看看tension：

　　輸出：

　　在上面，我們看到每種張力的三種不同類別(L，M和H)如何影響每種羊毛的斷裂。例如，低張力和A型羊毛的斷裂率最高。
 

　　我們還可以定義通過cat_plot() 使用 geom 參數(shù)創(chuàng)建的圖的類型 。此參數(shù)增強了對圖的解釋。我們可以這樣使用它，將其geom 作為附加參數(shù)傳遞 給 cat_plot：

　　輸出：

　　我們還可以通過添加plot.points = TRUE以下內容將觀察結果包括在圖中：

　　輸出：

　　還有許多其他設計選項，包括線條樣式，顏色等，使我們可以自定義這些可視化效果的外觀。有關詳細信息，請參閱此處的jtools文檔 。
 

　　使用速率數(shù)據(jù)進行泊松回歸建模
 

　　到目前為止，在大數(shù)據(jù)分析R中泊松回歸模型實例中，我們已經(jīng)對計數(shù)數(shù)據(jù)進行了建模，但是我們也可以對速率數(shù)據(jù)進行建模，該數(shù)據(jù)可以預測一段時間或分組內的計數(shù)數(shù)量。速率數(shù)據(jù)建模公式如下：
 

　　升 Ô 克(X / Ñ)= β 0 ?+Σ 我 β 我 X 我
 

　　這等效于：(應用對數(shù)公式)
 

　　升 Ô 克(X) - 升 Ö 克(Ñ)= β 0 ?+Σ 我 β 我 X 我
 

　　升 Ô 克(X)= 升 Ö 克(Ñ)+ β 0 ?+Σ 我 β 我 X 我
 

　　因此，可以通過將 系數(shù)為1 的log(n)項包括在內來對速率數(shù)據(jù)進行建模 。這稱為 offset。該偏移量用offset() R 建模 。
 

　　讓我們使用eba1977 從 ISwR包中調用的另一個數(shù)據(jù)集 為費率數(shù)據(jù)建立Poisson回歸模型。首先，我們將安裝該軟件包：

　　現(xiàn)在，讓我們看一下有關數(shù)據(jù)的一些詳細信息，并打印前十行以了解數(shù)據(jù)集包含的內容。

　　輸出：

　　為了對速率數(shù)據(jù)建模，我們使用 X / n ，其中 X 是要發(fā)生的事件， n 是分組。在此示例中， X =病例 (事件是癌癥病例)， n =流行病 (人群是分組)。
 

　　如上式所示，費率數(shù)據(jù)由 對數(shù)(n)解釋，在此數(shù)據(jù)中 n 為人口，因此我們將首先找到人口的對數(shù)。我們可以為案例/人群建模， 如下所示：

　　輸出：

　　現(xiàn)在，讓我們使用來對費率數(shù)據(jù)建模 offset()。

　　輸出：

　　在此數(shù)據(jù)集中，我們可以看到殘余偏差接近自由度，并且色散參數(shù)為 1.5(23.447 / 15) ，該值很小，因此該模型非常適合。
 

　　我們 fitted(model) 用來返回模型擬合的值。它使用建立模型的訓練數(shù)據(jù)返回結果。試一試吧：

　　使用該模型，我們可以使用該predict() 函數(shù)預測新數(shù)據(jù)集每1000人口的病例數(shù) ，就像我們之前對計數(shù)數(shù)據(jù)模型所做的那樣：

　　輸出： [1] 2.469818
 

　　因此， 對于40-54歲年齡段人群中的科靈市，我們可以預期每千人中大約有2或3例肺癌。
 

　　與計數(shù)數(shù)據(jù)一樣，我們也可以使用準泊松來獲得速率數(shù)據(jù)的更正確的標準誤差，但是出于大數(shù)據(jù)分析R中泊松回歸模型實例的目的，我們將不再重復該過程。
 

　　結論
 

　　泊松回歸模型在計量經(jīng)濟學和現(xiàn)實世界的預測中具有重要意義。在大數(shù)據(jù)分析R中泊松回歸模型實例中，我們學習了泊松分布，廣義線性模型和泊松回歸模型。
 

　　我們還學習了如何使用來為R中的計數(shù)數(shù)據(jù)和費率數(shù)據(jù)實現(xiàn)泊松回歸模型 glm()，以及如何將數(shù)據(jù)擬合到模型中以預測新的數(shù)據(jù)集。此外，我們研究了如何在glm()使用準泊松時獲得更準確的標準誤差， 并看到了一些可用于可視化的可能性 jtools。