具有實際概率的邏輯回歸∈(𝑎,𝑏)∈(a,b)in(a,b)在哪裡0<𝑎<𝑏<10<a<b<10<a<b`<1
使用邏輯回歸對概率建模時 $ ^1 $ ，擬合概率的範圍是 $ (0,1) $ . 邏輯函數 $ ^2 $ 漸近線 $ 0 $ 和 $ 1 $ ，所以這是一個很好的匹配。然而，在某些應用中，實際範圍 $ ^3 $ 的概率可以是 $ (a,b) $ 和 $ 0<a<b<1 $ ，導致尾部出現（潛在的大量）不匹配。
例如，考慮一個所有人都工作並獲得報酬的人群。付款取決於技能、努力和其他我們衡量的東西。每個人還參與抽獎，其貨幣結果正態分佈，期望值為零。個人的收入是工作收入和抽獎結果的總和。對於任何給定的常數 $ c $ , 條件概率 $ p $ 個人收入大於 $ c $ , 以技巧、努力和其他東西為條件，而不是彩票的結果，將滿足 $ 0<a<p<b<1 $ . （我相信那裡有更好的例子，但這是首先想到的。）
問題：
- 如何調整邏輯回歸模型 $ a $ 和 $ b $ 給出？（或者不需要調整？）
- 如何調整邏輯回歸模型 $ a $ 和 $ b $ 沒有給出，但我們知道 $ 0<a<b<1 $ ?
$ ^1 $ 或者就此而言，一個概率模型
$ ^2 $ 或者就此而言，標準的普通 CDF
$ ^2 $ 感興趣人群中個體的範圍
什麼時候 $ a $ 和 $ b $ 沒有給出，只需使用通常的邏輯模型（或任何合適的），因為（如果它使用合適的鏈接函數）保證擬合概率不小於下限 $ 0 $ 並且上限不大於 $ 1. $ 這些界限給出了區間估計 $ a $ 和 $ b. $
有趣的問題是什麼時候 $ a $ 和 $ b $ 是已知的。 您正在娛樂的那種模型似乎如下。您想到了一個單參數分佈系列 $ \mathcal{F} = {F_\theta} $ 在哪裡 $ \theta $ 對應於一些“概率”參數。例如， $ F_\theta $ 可能是伯努利 $ (\theta) $ 響應時的分佈 $ Y $ 是二進制的。
對於與解釋變量向量相關的觀察 $ x, $ 響應模型 $ Y_x $ 然後採取形式
$$ Y_x \sim F_{\theta(x)};\quad \theta(x) = g(x\beta) $$
對於一些“反向鏈接函數” $ g $ 我們必須指定：它是模型的一部分。例如，在邏輯回歸中， $ g $ 經常被認為是由定義的邏輯函數
$$ g(x) = \frac{1}{1 + \exp(-x)}. $$
無論細節如何，製作時 $ n $ 獨立觀察 $ y_i $ （每個都與一個向量相關聯 $ x_i $ ) 假設符合這個模型，他們的可能性是
$$ L(\beta) = \prod_{i=1}^n \Pr(Y_{x_i} = y_i\mid \theta(x_i) = g(x_i\beta)) $$
你可以像往常一樣繼續最大化它。（豎筆僅表示它後面的參數值決定使用哪個概率函數：它不是條件概率。）
讓 $ \hat\beta $ 是相關的參數估計。的預測條件分佈 $ Y_i $ 因此是
$$ Y_i \sim F_{\hat\theta(x_i)};\quad \hat\theta(x_i) = g(x_i\hat\beta). $$
當圖像 $ g $ 包含在區間內 $ [a,b], $ 那麼顯然每個 $ \hat\theta(x) $ 無論如何，也在那個區間 $ x $ 也許。（也就是說，這個結論既適用於 $ x $ 在數據集中並用於外推到其他 $ x. $ )
一個有吸引力的選擇 $ g $ 簡單地重新調整通常的邏輯函數，
$$ g(x;a,b) = \frac{g(x) - a}{b-a}. $$
**將此視為出發點：**像往常一樣，探索性分析和擬合優度測試將幫助您確定這是否適合 $ g. $
供以後使用，**請注意 $ g $ 和 $ g(;a,b) $ 有一個比看起來更複雜的關係，**因為最終它們被用來確定 $ \hat\beta $ 通過他們的論點 $ x\beta. $ 因此，這種關係的特徵在於函數 $ x\to y $ 取決於
$$ g(x) = g(y;a,b) = \frac{g(y) - a}{b-a}, $$
有溶液（如果 $ g $ 是可逆的，通常是這樣）
$$ y = g^{-1}((b-a)g(x) + a). $$
除非 $ g $ 最初是線性的，這通常是非線性的。
**為了解決此線程中其他地方表達的問題，**讓我們比較使用獲得的解決方案 $ g $ 和 $ g(;a,b). $ 考慮最簡單的情況 $ n=1 $ 觀察和需要估計參數向量的標量解釋變量 $ \beta=(\beta_1). $ 認為 $ \mathcal{F} $ 是二項式的家庭 $ (10,\theta) $ 分佈，讓 $ x_1 = (1), $ 並想像 $ Y_i = 9 $ 被觀察到。寫作 $ \theta $ 為了 $ \theta(x_1), $ 可能性是
$$ L(\beta) = \binom{10}{9}\theta^9(1-\theta)^1;\quad \theta = g((1)(\beta_1)) = g(\beta_1). $$
$ L $ 最大化時 $ g(\beta_1) = \theta = 9/10, $ 用獨特的解決方案$$ \hat\beta = g^{-1}(9/10) = \log(9/10 / (1/10)) = \log(9) \approx 2.20. $$
現在讓我們假設 $ a=0 $ 和 $ b=1/2: $ 也就是說，我們假設 $ \theta \le 1/2 $ 不管什麼價值 $ x $ 可能有。隨著縮放版本 $ g $ 我們和以前一樣計算，只是代入 $ g(;a,b) $ 為了 $ g: $
$$ L(\beta;0,1/2) = \binom{10}{9}\theta^9(1-\theta)^1;\quad \theta = g((1)(\beta_1);0,1/2) = g(\beta_1;0,1/2). $$
這不再最大化 $ \theta=9/10, $ 因為這是不可能的 $ g(\theta;0,1/2) $ 超過 $ 1/2, $ 按設計。 $ L(\beta;0,1/2) $ 被最大化的任何 $ \beta $ 那會讓 $ \theta $ 盡可能接近 $ 9/10; $ 這發生在 $ \beta $ 任意變大。那麼，使用受限反向鏈接函數的估計是
$$ \hat\beta = \infty. $$
顯然兩者都沒有 $ \hat\theta $ 或者 $ \hat\beta $ 是原始（無限制）估計的任何簡單函數；特別是，它們與任何重新縮放無關。
**這個簡單的例子暴露了整個程序的一個危險：**當我們假設什麼 $ a $ 和 $ b $ （以及關於模型的所有其他內容）與數據不一致，我們可能會得到模型參數的奇怪估計 $ \beta. $ 這就是我們付出的代價。
但是，如果我們的假設是正確的，或者至少是合理的呢？ 讓我們用 $ b=0.95 $ 代替 $ b=1/2. $ 這次， $ \hat\theta=9/10 $ 確實最大化似然性，由此估計 $ \beta $ 滿足
$$ \frac{9}{10} = g(\hat\beta;0,0.95) = \frac{g(\hat\beta) - 0}{0.95 - 0}, $$
所以
$$ g(\hat\beta) = 0.95 \times \frac{9}{10} = 0.855, $$
蘊涵
$$ \hat\beta = \log(0.855 / (1 - 0.855)) \approx 1.77. $$
在這種情況下， $ \hat\theta $ 不變但_ $ \hat\beta $ 發生了複雜的變化（ $ 1.77 $ 不是重新縮放的版本 $ 2.20 $ ).
在這些例子中， $ \hat\theta $ 當原始估計不在區間內時必須更改 $ [a,b]. $ 在更複雜的例子中，它可能必須改變才能改變其他觀察值的估計值 $ x. $ 這是效果之一 $ [a,b] $ 限制。另一個影響是，即使限制改變了估計的概率， $ \hat\theta, $ 原始反向鏈接之間的非線性關係 $ g $ 和受限鏈接 $ g(;a,b) $ 在參數估計中引起非線性（並且可能是複雜的）變化 $ \hat\beta. $
**為了說明，**我根據這個模型創建了數據 $ \beta=(4,-7) $ 和限制 $ a=1/10 $ 和 $ b=1/2 $ 為了 $ n $ 解釋值的等距值 $ x $ 之間 $ 0 $ 和 $ 1 $ 包含，然後使用普通邏輯回歸（無約束）擬合它們一次，並使用縮放反向鏈接方法再次使用已知約束。
以下是結果 $ n=12 $ 二項式 $ (8, \theta(x)) $ 觀察結果（實際上反映了 $ 12\times 8 = 96 $ 獨立的二進制結果）：
這已經提供了洞察力：模型（左圖）預測接近上限的概率 $ b=1/2 $ 對於小 $ x. $ 隨機變化導致某些觀測值的頻率大於 $ 1/2. $ 在沒有任何限制的情況下，邏輯回歸（中圖）傾向於預測更高的概率。類似的現象發生在大 $ x. $
受限模型極大地改變了估計的斜率 $ -3.45 $ 到 $ -21.7 $ 為了使預測保持在 $ [a,b]. $ 發生這種情況的部分原因是它是一個小數據集。
*直觀地說，*更大的數據集應該產生更接近底層（真實）數據生成過程的結果。人們也可能期望不受限制的模型運行良好。可以？為了檢查，我創建了一個大一千倍的數據集： $ n=1200 $ 二項式的觀察 $ (80,\theta(x)) $ 回复。
當然，正確的模型（右面板）現在非常適合。然而，觀察到的頻率的隨機變化仍然導致普通邏輯模型超出限制。
顯然，當假定值 $ a $ 和 $ b $ 是（接近）正確的並且鏈接函數的形狀大致正確，最大似然效果很好 - 但它絕對不會產生與邏輯回歸相同的結果。
為了提供完整的文檔，這裡是
R生成第一個圖的代碼。更改12為1200和8生成80第二個圖形。# # Binomial negative log likelihood. # logistic.ab <- function(x, a=0, b=1) { a + (b - a) / (1 + exp(-x)) } predict.ab <- function(beta, x, invlink=logistic.ab) { invlink(cbind(1, x) %*% beta) } Lambda <- function(beta, n, k, x, invlink=logistic.ab, tol=1e-9) { p <- predict.ab(beta, x, invlink) p <- (1-2*tol) * p + tol # Prevents numerical problems - sum((k * log(p) + (n-k) * log(1-p))) } # # Simulate data. # N <- 12 # Number of binomial observations x <- seq(0, 1, length.out=N) # Explanatory values n <- rep(8, length(x)) # Binomial counts per observation beta <- c(4, -7) # True parameter a <- 1/10 # Lower limit b <- 1/2 # Upper limit set.seed(17) p <- predict.ab(beta, x, function(x) logistic.ab(x, a, b)) X <- data.frame(x = x, p = p, n = n, k = rbinom(length(x), n, p)) # # Create a data frame for plotting predicted and true values. # Y <- with(X, data.frame(x = seq(min(x), max(x), length.out=101))) Y$p <-with(Y, predict.ab(beta, x, function(x) logistic.ab(x, a, b))) # # Plot the data. # par(mfrow=c(1,3)) col <- hsv(0,0,max(0, min(1, 1 - 200/N))) with(X, plot(x, k / n, ylim=0:1, col=col, main="Data with True Curve")) with(Y, lines(x, p)) abline(h = c(a,b), lty=3) # # Reference fit: ordinary logistic regression. # fit <- glm(cbind(k, n-k) ~ x, data=X, family=binomial(link = "logit"), control=list(epsilon=1e-12)) # # Fit two models: ordinary logistic and constrained. # for (ab in list(c(a=0, b=1), c(a=a, b=b))) { # # MLE. # g <- function(x) logistic.ab(x, ab[1], ab[2]) beta.hat <- c(0, 1) fit.logistic <- with(X, nlm(Lambda, beta.hat, n=n, k=k, x=x, invlink=g, iterlim=1e3, steptol=1e-9, gradtol=1e-12)) if (fit.logistic$code > 3) stop("Check the fit.") beta.hat <- fit.logistic$estimate # Check: print(rbind(Reference=coefficients(fit), Constrained=beta.hat)) # Plot: Y$p.hat <- with(Y, predict.ab(beta.hat, x, invlink=g)) with(X, plot(x, k / n, ylim=0:1,, col=col, main=paste0("Fit with a=", signif(ab[1], 2), " and b=", signif(ab[2], 2)))) with(Y, lines(x, p.hat, col = "Red", lwd=2)) with(Y, lines(x, p)) abline(h = c(a, b), lty=3) } par(mfrow=c(1,1))
