Regression

將分位數回歸公式化為線性規劃問題?

  • December 29, 2018

如何將分位數回歸表述為線​​性規劃問題?在查看中位數分位數問題時,我知道它是

minimize ni=1|β0+Xiβ1Yi| transforms into  minimize ni=1ei s.t. eiβ0+Xiβ1Yi ei(β0+Xiβ1Yi)

但是如何轉換其他分位數的最小化?

您使用分位數回歸估計器

ˆβ(τ):=argminθRKNi=1ρτ(yixiθ).

在哪裡 τ(0,1) 是根據需要估計的分位數和函數選擇的常數 ρτ(.) 定義為

ρτ(r)=r(τI(r<0)).

為了看清目的 ρτ(.) 首先考慮將殘差作為參數,當它們被定義為 ϵi=yixiθ . 因此,最小化問題中的總和可以重寫為

Ni=1ρτ(ϵi)=Ni=1τ|ϵi|I[ϵi0]+(1τ)|ϵi|I[ϵi<0]

使得與觀察相關的正殘差 yi 高於建議的分位數回歸線 xiθ 被賦予權重 τ 而與觀察相關的負殘差 yi 低於建議的分位數回歸線 xiθ 加權 (1τ) .

直觀地說:

τ=0.5 正負殘差以相同的權重被“懲罰”,並且相同數量的觀察值在最優的“線”上方和下方,因此該線 xiˆβ 是中位數回歸“線”。

什麼時候 τ=0.9 每個正殘差的權重是負殘差的 9 倍 1τ=0.1 因此對於“線”上方的每個觀察都是最佳的 xiˆβ 大約 9 將放置在該線下方。因此,“線”代表 0.9 分位數。(有關這方面的確切陳述,請參閱 Koenker (2005)“分位數回歸”中的 THM.2.2 和推論 2.1)

這些圖中說明了這兩種情況。左面板 τ=0.5 和右面板 τ=0.9 .

rho 函數 tau = 0.5 和 tau = 0.9

線性程序主要使用標準形式進行分析和求解

(1)  minz  cz  subject to Az=b,z0

要達到標準形式的線性程序,第一個問題是在這樣的程序中 (1) 執行最小化的所有變量 z 應該是積極的。為了實現這一點,使用鬆弛變量將殘差分解為正負部分:

ϵi=uivi

積極的部分 ui=max(0,ϵi)=|ϵi|I[ϵi0]vi=max(0,ϵi)=|ϵi|I[ϵi<0] 是負面的部分。由檢查函數分配權重的殘差總和被視為

Ni=1ρτ(ϵi)=Ni=1τui+(1τ)vi=τ1Nu+(1τ)1Nv,

在哪裡 u=(u1,,uN)v=(v1,,vN)1N 是向量 N×1 所有坐標等於 1 .

殘差必須滿足 N 約束

yixiθ=ϵi=uivi

這導致公式為線性程序

minθRK,uRN+,vRN+τ1Nu+(1τ)1Nv|yi=xiθ+uivi,i=1,,N,

如 Koenker (2005)“分位數回歸”第 10 頁方程 (1.20) 所述。

然而值得注意的是 θR 仍然不限於標準形式(1)上的線性程序中要求的為正。因此再次使用分解成正負部分

θ=θ+θ

又在哪裡 θ+=max(0,θ) 是積極的部分並且 θ=max(0,θ) 是負面的部分。這 N 約束可以寫成

y:=[y1  yN]=[x1  xN](θ+θ)+INuINv,

在哪裡 IN=diag1N .

接下來定義 b:=y 和設計矩陣 X 將自變量上的數據存儲為

X:=[x1  xN]

重寫約束:

b=X(θ+θ)+INuINv=[X,X,IN,IN][θ+ θ u v]

定義 (N×2K+2N) 矩陣

A:=[X,X,IN,IN]

並介紹 θ+θ 作為要最小化的變量,因此它們是 z 要得到

b=A[θ+ θ u v]=Az

因為 θ+θ 僅通過約束 a 影響最小化問題 0 維度的 2K×1 必須作為係數向量的一部分引入 c 可以適當地定義為

c=[0 τ1N (1τ)1N],

從而確保 $ c^\top z = \underbrace{\mathbf 0^\top(\theta^+ - \theta^-)}{=0}+\tau \mathbf 1_N^\top u + (1-\tau)\mathbf 1_N^\top v = \sum{i=1}^N \rho_\tau(\epsilon_i). $

因此 c,Ab 然後定義和程序如 (1) 完全指定。

這可能最好用一個例子來消化。要在 R 中解決這個問題,請使用 Roger Koenker 的 quantreg 包。這裡還說明瞭如何設置線性程序並使用線性程序求解器求解:

base=read.table("http://freakonometrics.free.fr/rent98_00.txt",header=TRUE)
attach(base)
library(quantreg)
library(lpSolve)
tau <- 0.3


# Problem (1) only one covariate
X <- cbind(1,base$area)
K <- ncol(X)
N <- nrow(X)

A <- cbind(X,-X,diag(N),-diag(N))
c <- c(rep(0,2*ncol(X)),tau*rep(1,N),(1-tau)*rep(1,N))
b <- base$rent_euro
const_type <- rep("=",N)

linprog <- lp("min",c,A,const_type,b)
beta <- linprog$sol[1:K] - linprog$sol[(1:K+K)]
beta
rq(rent_euro~area, tau=tau, data=base)


# Problem (2) with 2 covariates
X <- cbind(1,base$area,base$yearc)
K <- ncol(X)
N <- nrow(X)

A <- cbind(X,-X,diag(N),-diag(N))
c <- c(rep(0,2*ncol(X)),tau*rep(1,N),(1-tau)*rep(1,N))
b <- base$rent_euro
const_type <- rep("=",N)

linprog <- lp("min",c,A,const_type,b)
beta <- linprog$sol[1:K] - linprog$sol[(1:K+K)]
beta
rq(rent_euro~ area + yearc, tau=tau, data=base)

引用自:https://stats.stackexchange.com/questions/384909