本文想在R软件中更好地了解分位数回归优化。在查看分位数回归之前，让我们从样本中计算中位数或分位数。

由Kaizong Ye，Liao Bao撰写

考虑一个样本。要计算中位数，请求解

可下载资源

完整程序、数据和文档（word）

可以使用线性编程技术解决。更确切地说，这个问题等同于

为了说明，考虑对数正态分布的样本，

分位数回归提出的原因，就是因为不希望仅仅是研究y的期望，而是希望能探索y的完整分布状况，或者说可能在某些情况下我们更希望了解y的某个分位数。

首先随机变量X的分布函数为：

$F(x)=P(X \leq x)$

则它的 $\tau$ 分位数可定义为： $Q_{\tau}(X)=arginf \left\{ x \in \mathbb{R}; F(x)\geq\tau \right\} (0<\tau<1)$

若将分布函数F（x）的逆定义为：

$F _{X}^ {-1}(\tau)=inf[y\in \mathbb{R}; F(y)\geq\tau]$

则

$Q_{\tau}(X)=F_{X}^{-1}(\tau)$

设有随机向量（X，Y），其中Y在给定X=x的情况下的条件累积分布函数为 $F_{Y|X=x}(y|x)$ ，则将该条件随机变量Y|X=x的 $\tau$ 条件分位数定义为： $Q_{\tau}(Y|X = x)=arginf \left\{ y\in \mathbb{R}:F(y|x)\geq\tau \right\} (0<\tau<1)$

假定我们有样本序列 $\left\{ (X_{i},Y_{i} ),(i=1,...,n)\right\}$ 满足下列回归模型，即

$Y=m(X)+\epsilon \quad X \in \mathbb{R}$

假定误差项 $\epsilon_{i}\left\{ i=1,...,n \right\}$ 为独立同分布的序列，且分布情况未知，这Y的 $\tau$ 条件分位数 $m_{\tau}(x)=argmin_{\theta \in \mathbb{R}}\mathbb{E}\left\{ \rho_{\tau}(Y-\theta)|X=x\right\}$

式中， $\rho_{\tau}(u)=u[\tau I(u \geq 0)-(1-\tau)I(u<0)]$ 为分位回归领域的损失函数，有些地方也称之为检验函数。 $I(·)$ 为示性函数。不包含示性函数的损失函数为

$\rho_{\tau}(x)=\left\{ \begin{aligned} \tau x ,u\geq0 \\ (\tau-1)x,u<0 \\ \end{aligned} \right.$

损失函数同样可以用下面一个公式进行表示：

$\rho_{\tau}(u)=u(\tau-I(u<0))$

则我们希望 $E\rho_{\tau}(x-\hat X)$ 尽可能的小，最好等于0，则： $E\rho_{\tau}(x-\hat X)=(\tau-1)\int_{-\infty}^{\hat X}(x-\hat X)dF(x)+\tau\int_{\hat X}^{\infty}(x-\hat X)dF(x)$

最优解为

$0=(1-\tau )\int_{-\infty}^{\hat{X}}dF(x)-\tau\int_{\hat{X}}^{\infty}dF(x)=F(\hat{ X})-\tau$

用样本经验分布函数代替F（x）则

$F(x)=\frac{\sum_{i=1}^{n}{I(X_{i}\leq x)}}{n}$

$\int \rho_{\tau}(x-\hat X)dF_{n}(x)=n^{-1}\sum_{i=1}^{n}\rho_{\tau}(x_{i}-\hat X)$


n = 123 
set.seed(132)
y = rlnorm(n)
median(y)
[1] 1.01523

对于优化问题，使用具有3n个约束和2n + 1参数的矩阵形式，


 
r = lp("min", c(rep(1,2*n),0),
 
tail(r$solution,1) 
[1] 1.01523

分位数

当然，我们可以将之前的代码改编为分位数

课程

R语言数据分析挖掘必知必会

从数据获取和清理开始，有目的的进行探索性分析与可视化。让数据从生涩的资料，摇身成为有温度的故事。

立即参加 ➜


tau = .3
quantile(x,tau)
      30% 
0.674124

线性程序

R代码

 
r = lp("min", c(rep(tau,n),rep(1-tau,n),0),
 
 
[1] 0.674124

分位数回归（简单）

考虑一个数据集，该数据集是一个主要城市的单位租金与面积，建筑年龄等的函数。

分位数回归的线性程序

与ai，bi≥0和

在这里使用

require(lpSolve) 
 
r = lp("min",
       c(rep(tau,n , rep(1-tau,n),0,0 , rbind(A1, A2 ,
       c(rep( =", 2*n , rep("=", n) , c(rep(0,2*n), y 
tail(r$solution,2)
[1] 147.845234   3.273453

我们可以使用R函数来拟合该模型

 
rq(ren~are , tau=tau 
Coefficients:
(Intercept)        are 
 147.845234   3.273453

我们可以使用不同的概率水平来获得图

plot( area, rent,xlab=expression
tau = .9
r = lp("min",
       c(re au,n), rep(1-tau  rbind(A1 2),
       c(rep , 2*n), rep("=", n)), c( ,2*n) y))

最受欢迎的见解

1.R语言多元Logistic逻辑回归应用案例

2.面板平滑转移回归(PSTR)分析案例实现

3.matlab中的偏最小二乘回归（PLSR）和主成分回归（PCR）

4.R语言泊松Poisson回归模型分析案例

5.R语言回归中的Hosmer-Lemeshow拟合优度检验

6.r语言中对LASSO回归，Ridge岭回归和Elastic Net模型实现

7.在R语言中实现Logistic逻辑回归

8.python用线性回归预测股票价格

9.R语言如何在生存分析与Cox回归中计算IDI，NRI指标

多元分位数回归

现在，我们尝试使用两个协变量呢，例如，让我们看看是否可以将单位的租金解释为面积的（线性）函数和建筑年龄。

 
r = lp("min",
       c(rep(ta n), rep(1- au,n),0,0, , rbin 1, A2),
        (r p("&  ,  n), rep("=  n)),  (rep(0 *n), y)) 
tail(r$sol ,3)
[1] 0.000  3.224  0.073

 
Coefficients:
 (Intercept)         are         year 
-5322.503252     3.428135     2.637234

结果是完全不同的。可以用IRLS –迭代加权最小二乘确认后者

for(s in 1:500){
 
  reg = lm(rent ~area+year ,
weigts= tau*(eps t;0 1-tau) eps&lt ))/ s(e ))
 
}
reg$coefficients
 (Intercept)         area        year
-5485.433043     3.932134     2.842943

我们可以使后者拟合多元回归，

 
lp("min",c,A consttype,b)
beta = r$sol[1:K  -  r$sol (1:K+K) 
beta
[1] -5542.633252     3.958135     2.857234

与之比较

 
rq(rent~ area + year, tau=tau 
 
Coefficients:
 (Intercept)         area        yearc 
-5542.633252     3.958135     2.857234
 
Degrees of freedom: 4571 total; 4568 residual

可下载资源

关于作者

Kaizong Ye是拓端研究室（TRL）的研究员。

本文借鉴了作者最近为《R语言数据分析挖掘必知必会 》课堂做的准备。

非常感谢您阅读本文，如需帮助请联系我们！