AIC和BIC都是用于同时衡量模型拟合度和复杂度的统计量。换句话说，我们希望利用AIC或者BIC选择出来的模型既有不错的拟合度又不至于太复杂（太多自变量）。

对于一个回归模型，如果优化目标只是最大化log-likelihood函数，即不限制变量系数或变量数量，当然是模型越大（变量越多）模型越精确。考虑到模型计算以及过度拟合，我们又倾向于一个简单的模型。所以我们需要寻找一个相对最优化的模型，来平衡模型大小与模型拟合准确度的关系。选择模型的评估标量有很多，AIC和BIC就是其中的两个。它们的相似处在于都是通过限制模型变量的数量来控制模型的大小，而不同的是惩罚变量数量的函数不同。

具体统计上来说，

$$\text{AIC}=2p - 2\log L(\hat{\theta}) $$

其中$\theta$是模型参数，p是模型参数的数量。从模型选择上说，我们要选择AIC最小的那个模型。

$$\text{BIC}=p\log(n)-2\log L(\hat{\theta}) $$

其中n是样本数量。类似地，我们优先选择BIC最小的模型。

我们可以注意到，当$1/2\times \log(n) > 1$, 也就是，样本数量大于$e^2 = 7.39$时, BIC的复杂度惩罚函数比AIC的惩罚函数大。换句话说，BIC相对于AIC 会选择一个更小的模型。

以上是通过统计角度来简单解释 AIC和BIC。在信息论角度，我们也可以用最小描述长度原则来解释这两个评估值。

p.s. BIC的概念来源于贝叶斯，有知道的童鞋可以解释一下为什么吗？

线性回归需要满足四个前提假设

1. Linearity 线性

应变量和每个自变量都是线性关系。

2. Indpendence 独立性

对于所有的观测值，它们的误差项相互之间是独立的。

3. Normality 正态性

误差项服从正态分布。

4. Equal-variance 等方差

所有的误差项具有同样方差。

这四个假设的首字母，合起来就是LINE，这样很好记。

如果这些假设不满足，最终的回归模型可能不精确乃至比较离谱。其中最重要的假设是线性，如果自变量和应变量本身并非线性关系，那么拟合出来的模型当然会有偏差。比如说，人的身高和年龄，如果是1到15岁，这个线性模型看起来似乎问题不大。那么1到60岁，身高和年龄还是线性关系吗？显然就不是了吧。

回归里出现的p值也是针对于假设检验来说的。

假设你的回归模型是$Y=aX_1+bX_2+c.$

$a$所对应的假设检验中，零假设是在$b$和$c$都是正确值得情况下$a=0$，对立假设是在$b$和$c$都正确的情况下$a\neq 0$。这一般都是采用双侧t检验。$a$所对应的p值就是这个假设检验的p值。

是对每个回归系数做T检验。

$H_0$：回归系数为0

$H_1$：回归系数不为0.

一般做T检验是需要很多样本的，然后根据样本均值和方差来求p值。在回归中，回归系数的方差是可以利用公式求出来的，而样本均值就用回归系数的值来代替，这样就可以得到p值了。

方差膨胀因子(Variance Inflation Factor，VIF)是用来检测回归模型中自变量的多重共线性(Multi-Colinearity)的。

我们可以用相关系数来比较一个自变量与另一个自变量的相关性。如果我们想判断一个自变量与其他剩下所有的自变量的多重共线性，我们就可以用VIF。

比如我们一共有5个自变量，$X_1, X_2, X_3, X_4, X_5$。我们把$X_1$当作观测值，把剩下的变量$X_2, X_3, X_4, X_5$当作自变量，进行线性回归。

$$X_1 = c_{1,0}+c_{1,2}X_2+c_{1,3}X_3+c_{1,4}X_4+c_{1,5}X_5+e_1$$

我们用$R^2_1$上表示面这个回归的$R^2$。$X_1$的方差膨胀因子

$$VIF_1=\frac{1}{1-R^2_1},$$

相应地，变量$X_k$的方差膨胀因子

$$VIF_k=\frac{1}{1-R^2_k},$$

其中$R^2_k$是把$X_k$当作观测值，把剩下的变量当作自变量的回归的$R^2$。

我们知道$R^2$越接近1，拟合越好。这里$VIF_k$等于1说明没有共线性，越接近1，$X_k$的多重共线性越小。$VIF_k$越大，$X_k$与其他变量的多重共线性越大。

关于上面的说的解答有一个问题，如果说 Rk2  越接近1 代表拟合的越好，那不是越说明其他的变量可以更好的线性替代第k个变量么？那难道不应该是第k个变量的多重共线性更严重么？这样看来岂不是VIF越大(R2越小）反而应该共线性越小么？（因为说明其他变量拟合该变量的效果不好，不能被线性替代），谁来解释一下？

对于一元简单线性回归

$$Y=\beta X + \beta_0 + \epsilon$$

我们知道

$$\beta = \frac{\text{Cov}(X,Y)}{\text{Var}(X)}=\text{Cor}(X,Y)\frac{\sigma_Y}{\sigma_X}$$

可以看出$\beta$的确是和$\text{Cor}(X,Y)$正相关的（标准差$\sigma_X$和$\sigma_Y$都是正数）。

如果我们把自变量$X$和因变量$Y$都先进行标准化处理，使得$\sigma_X=\sigma_Y=1$，那么回归系数$\beta$就等于$X$和$Y$的皮尔逊相关系数了。

确实学习的参数正比于输入和输出皮尔逊相关系数。

用列向量矩阵表示的最小二乘公式   

$$\hat{y}^T=X_{new}^T(XX^T)^{-1}Xy^T$$

$\Sigma_{XX}=XX^T$是covariance matrix,是半正定的实对称方阵，可做eigendecomposition， $XX^T=Q\Lambda Q^T$。其中$Q$是$\Sigma_{XX}$的eigen vector，把输入$X$投影到相互独立的eigen space，对角线矩阵$\Lambda$表示每个eigen vector轴上的variance，

covariance matrix逆是partial covariance matrix,  $(XX^T)^{-1}=Q\Lambda^{-1} Q^T=QD(QD)^T$，其中对角线矩阵$D=D^T, DD^T=\Lambda^{-1}$

$$\hat{y}^T=(X_{new}^TQD)(X^TQD)^Ty^T$$

令 $Z^T=X^TQD$

$$\hat{y}^T=Z_{new}^T Zy^T$$

$$\hat{y}^T=Z_{new}^T\Sigma_{Zy}$$

数据流程是：$X$经过独立化算子$Q$和标准化算子$D$得到新特征$Z$,最小二乘模型学习的是$\Sigma_{Zy}$。如果$y$经过标准化（$\sigma_y=1$），则学习的系数就是$Corr(Z,y)$。

一个理想情况是$X$本身个特征独立并输入输出都经过标准化，$X=Z$，则$w=corr(X,y)$。

----------------------------

这满足以最小二乘为代表的基于随机变量相关性的机器学习方法的流程：先把数据压缩成互不相关的信号（要求编码向量$Z$正交），再学习编码和输出的相关性（$\Sigma_Zy$），最后把新编码$Z_{new}$解码到输出空间。autuencode也符合这个流程。

自变量和因变量的相关系数和正规化后的回归系数是相等，它们的平方就等于$R^2$

y ~ x1:x2 等价于  y ~ x1x2

y ~ x1*x2 等价于 y ~ x1 + x2 + x1x2

冒号:，只表示两个变量的乘积

乘号*，表示两个变量以及它们的乘积

逻辑回归的系数和概率值并不像线性回归的系数有那么直接的线性关系。

逻辑回归模型是

$$\text{Logit}(p)=\log\left(\frac{p}{1-p}\right)=\beta_0+\beta_1x_1+\beta_2x_2+\cdots+\beta_px_p.$$

所以，Logit(p)和变量的系数是线性关系，也就是log odds和系数是成线性关系。但是当从log odds转化为odds，再转化为概率p的时候，却没有线性关系，而且它的数值变化和x的取值也有关。

举例来说，$$\text{Logit}(p)=\log\left(\frac{p}{1-p}\right)=1+\beta x$$，对于$x=2$来说，

上图可以看出，系数beta的变动和Logit(p)的变动是线性的。

但是odds $\frac{p}{1-p}$和beta就不再是线性的了。

对于最后的概率，就更难直接从beta的角度进行解释。

所以我们一般是阐释系数和log odds或者logit函数的线性关系，也就是在逻辑回归当中，x每增加一个单位，logit函数增加beta个单位。

logitic regression是一种广义线性模型，你可以理解为是几率的对数 和 输入成线性关系。

以一元logistic回归:$log( \frac{\pi (x)}{1-\pi (x)})=\alpha + \beta x$为例，

可以有两种解释：

• 在给定自变量$x = x_1$的情况下，$x$每增加一个单位，$\pi (x)$的值增加$\beta \pi (x_1)[1-\pi (x_1)]$
• 自变量$x$每增加1，其优势（即@沛县李刚 说的odds ）为原来的$e^\beta$次方倍

展开来说：

第一点：可以根据logistic的公式，将$\pi (x)$提至等式左边，令$\pi (x)$对$x$求导得到。

第二点：优势的定义是$\frac{\pi }{1-\pi}$，可以理解为一般一次失败会伴随着几次成功。推导过程如下：

    $\pi (x)\times \frac{1}{1-\pi (x)} $

$=\frac{exp(\alpha+\beta x)}{1+exp(\alpha + \beta x)}\times (1+exp(\alpha + \beta x))$

$=exp(\alpha+\beta x) $

$=e^\alpha e^{\beta x} $

$=e^\alpha (e^\beta )^x$

如果只是定性地来看的话，在特征都已经标准化（或者归一化）的情况下，系数的绝对值越大，说明这个特征越重要。

若系数为正，这个特征与目标值为1的概率正相关；若系数为负，这个特征与目标值为0的概率正相关。

在线性回归中，库克距离(Cook's Distance)描述了单个样本对整个回归模型的影响程度。库克距离越大，说明影响越大。库克距离也可以用来检测异常点。

在最理想的情况下，每个样本对模型的影响是相等的。如某个样本的库克距离非常大，我们可以视为这个样本是异常点(outlier)。通常来说，若库克距离大于1，我们就认为这个点是异常点。也有人把这个阈值设置为$4/n$。

----------如果想知道具体的计算公式，欢迎继续阅读----------

在线性回归中，

$$y = X\beta + \epsilon,$$

$y\in\mathbb{R}^n$, $X\in\mathbb{R}^{n\times p}$, $\beta\in\mathbb{R}^p$, $\epsilon\in \mathbb{R}^{p}$。 这里考虑的是$n$个观测值的样本，自变量的数量为$p$，$\epsilon$是误差向量。

我们知道这个线性回归的投影矩阵$H=X(X^TX)^{-1}X^T$，$h_i=x_i^T(X^TX)^{-1}x_i$是矩阵$H$对角线上的第$i$个元素。第$i$个样本的库克距离

$$D_i=\frac{\epsilon_i^2 h_i}{s^2 p (1-h_i)^2},$$

其中$s^2$是这个模型的均方误差MSE。

下面是几种常用的。

$y=(y_1,y_2,\cdots,y_n)$是真实值。

$\hat y = (\hat y_1, \hat y_2,\cdots, \hat y_n)$是预测值。

Mean Absolute Error(MAE)，绝对误差均值

$$MAE=\text{Mean}_{i=1}^n \left\{|y_{i}-\hat y_{i}|\right\}$$

Mean Absolute Percentage Error(MAPE)，绝对百分比误差均值

$$MAPE=\text{Mean}_{i=1}^n\left\{ \frac{|y_{i}-\hat y_{i}|}{|y_{i}|}\right\}$$

Median Absolute Error(MedianAE)，绝对误差均中位数

$$MedianAE=\text{Median}_{i=1}^n \left\{|y_{i}-\hat y_{i}|\right\}$$

Median Absolute Percentage Error(MedianAPE)，绝对百分比误差中位值

$$MedianAPE=\text{Median}_{i=1}^n\left\{ \frac{|y_{i}-\hat y_{i}|}{|y_{i}|}\right\}$$

Mean Squared Error(MSE)，平方误差均值

$$MSE=\text{Mean}_{i=1}^n\left\{ \left(y_{i}-\hat y_{i}\right)^2\right\}$$

Median Squared Error，平方误差中位数

$$MedianSE=\text{Median}_{i=1}^n \left\{\left(y_{i}-\hat y_{i}\right)^2\right\}$$

Explained Variance Score，方差解释值

$$EVS=1-\frac{\text{Var}(y-\hat y)}{\text{Var}(y)}$$

借花献佛一下。

还有修正R方(adjusted R square)。

在右侧相关问题里看见的，嘿嘿嘿。

负数的AIC是可能的，不用担心是错误。

对于离散的概率分布这是不可能的，但是对于连续的概率分布，这个是完全可以发生的。

我们知道AIC的公式是，$AIC=2p−2\log L(\hat \theta)$。

当似然函数大于1的时候，似然函数的对数就是大于0，AIC就是有可能是负数的。

修正$R_{adj}^2$的公式是

$$R_{adj}^2=1-\frac{(n-1)(1-R^2)}{n-p-1},$$

其中$n$是样本数量，$p$是模型中变量的个数。

我们知道在其他变量不变的情况下，引入新的变量，总能提高模型的$R^2$。修正$R^2$就是相当于给变量的个数加惩罚项。

换句话说，如果两个模型，样本数一样，$R^2$一样，那么从修正$R^2$的角度看，使用变量个数少的那个模型更优。使用修正$R^2$也算一种奥卡姆剃刀的实例。

adjusted $R^2$是考虑到了自由度下的$R^2$.

$$R^2=1-\frac{SS_{\text{Res}}}{SS_{\text{Total}}}$$

考虑到残差的平方和$SS_{\text{Res}}$的自由度为$n-p-1$，总体平方和$SS_{\text{Total}}$的自由度为$n-1$，那么我们修正后的$R^2_{\text{Adj}}$的公式为

$$R^2_{\text{Adj}}=1-\frac{SS_{\text{Res}}/(n-p-1)}{SS_{\text{Total}}/(n-1)}=1-\frac{SS_{\text{Res}}}{SS_{\text{Total}}}\frac{(n-1)}{(n-p-1)}=1-(1-R^2)\frac{n-p-1}{n-1}$$

公式中$n$是样本的个数，$p$是变量的个数

R方可以用来评价模型的拟合程度。当我们在评价拟合程度的同时，也考虑到模型的复杂程度，那么就是修正R方。

原理上与AIC和BIC类似。

这里的确是在做F检验。

零假设是：所有变量的回归系数都是0，也就是这个模型是个常数模型。

你的图里的p值是很小的，所以可以认为拒绝零假设。

HuberLoss是一个用于回归问题的带参数的损失函数。

$$L_{\delta}(y,f(x))=\begin{cases}\frac{1}{2}(y-f(x))^2,\text{ if }|y-f(x)|\leq \delta \\  \delta|y-f(x)| -\frac{1}{2}\delta^2, \text{ if }|y-f(x)|> \delta  \end{cases}$$

$\delta$是HuberLoss的参数，$y$是真实值，$f(x)$是模型的预测值。

当预测偏差小于$\delta$时，它采用平方误差，当预测偏差大于$\delta$的时候，采用的线性误差。相比于最小二乘的线性回归，HuberLoss降低了对outlier的惩罚程度，所以HuberLoss是一种常用的robust regression的损失函数。

最后贴一张维基上的图。

绿色是参数$\delta=1$的Huber Loss，蓝色是平方误差。

huber loss是最小二乘法和最小绝对误差的结合体。

当数据中有不少离群值或者不确定数据中是否有离群值时，可以用huber loss。

如果预测值的分布是长尾的，看起来不大正态，不大对称的话，可以考虑用Huber。Huber的鲁棒性更强一些。

Quantile Regression可以看作广义的LAD Regression.

假设$\text{y}=(y_1,y_2,\cdots, y_n)$为真实值，$\hat{\text{y}}=(\hat{y_1},\hat{y_2},\cdots,\hat{y_n})$为预测值，那么参数为$q$的Quantile Regression的目标函数为

$$L_q(\text{y},\hat{\text{y}})=\frac{1}{n}(1-q)\sum_{i:\hat{y_i}\geq y_i}(\hat{y_i} - y_i)+\frac{1}{n}q\sum_{j:\hat{y_j}< y_j}(y_j-\hat{y_j} )$$

当$q=0.5$时，Quantile Regression就退化为LAD Regression，也就是所谓的中位数回归。

$$L_{0.5}(\text{y},\hat{\text{y}})=\frac{1}{2n}\sum_{i=1}^n|\hat{y_i}-y_i|$$.

当$q>0.5$时，目标函数$L_q$对预测值偏小的结果惩罚更大。

当$q<0.5$时，目标函数$L_q$对预测值偏大的结果惩罚更大。

泊松回归是GLM中的一种，link function是log，但是并不等价于“对y取对数，然后做线性回归”。

对log(y)做线性回归，是表示log(y)服从正态分布。泊松回归是需要y服从泊松分布的，也就是

$$P(y_i|x_i, \theta)=\frac{1}{y_i!}\text{exp}(y_i\theta^Tx_i-e^{\theta^Tx_i})$$

$y_i$的均值是

$$\hat y_i = e^{\theta ^T x_i}$$

$\theta$是回归系数。求解$\theta$是根据对数最大似然估计，

$$l(y|x,\theta)=\sum_{i=1}^n\left(y_i\theta^T x_i-\text{exp}(\theta^Tx_i) - \log(y_i!)\right)$$

Poisson regression实际上就是把y取log，然后做线性回归。

log(y) = ax + b

训练出a和b，然后再做预测

log(y') = ax' + b

得到预测值y' = exp(ax' + b).

是的。它是凸的。这就是为什么我们可以用SGD来解。

线性回归是个更广的概念，只要能够写成$y = ax + b$的形式都是线性回归模型。

和很多机器学习算法一样，损失函数才是模型的核心。

当我们希望预测值和真实值的平方和尽量小时，也就是损失函数为

$$\sum_{i=1}^n(y_i-\hat y_i)^2,$$

这个就是最小二乘线性回归。

当我们希望预测值和真实值的绝对值的和尽量小时，也就是损失函数为

$$\sum_{i=1}^n|y_i-\hat y_i|,$$

这个就是LAD线性回归。

我们还有其他的不同的更复杂的损失函数，造就了不同的线性回归模型，比如分位数回归，Huber回归，等等。

由于最常用的、最简单的线性回归模型就是最小二乘法，所以通常说到线性回归的时候其实就是指的最小二乘线性回归。

最小二乘法(Ordinary Least Squares, OLS)只是一种线性回归模型，但是是最最常用的，所以当我们说线性回归的时候通常实际上都是指的最小二乘线性回归。

当然有其他线性回归模型，比如最小绝对偏差(Least Absolute Deviation, LAD)回归模型。

这里的$L_1$和$L_2$是指的向量范数的种类。

向量$x=(x_1, x_2,\cdots, x_n)$的$L_1$范数表示为

$$\|x\|_1=\sum_{i=1}^n|x_i|$$

比如说$x=(3, -4)$，那么$\|x\|_1=3+4=7$

向量$x=(x_1, x_2,\cdots, x_n)$的$L_2$范数表示为

$$\|x\|_2=\sqrt{\sum_{i=1}^nx_i^2}$$

比如说$x=(3, -4)$，那么$\|x\|_2=\sqrt{3^2+4^2}=5$

这些范数在Lasso和Ridge中的具体体现，可以参考Ridge，Lasso目标函数的表达式

L0 norm是一组向量里不为0的元素的个数

比如一个向量 v = (5, 2, 0, 3, 4)

那么L0(v) = 4.

从正则化的角度来说，L0的作用是控制模型引入特征的数量。相比于L1（Lasso）和L2（Ridge），L0更加严格地控制特征的数量。

严格来说，L0 norm不是norm，因为它不满足范数的定义。

$L_0$从不满足范数定义上的数乘：

$$\|\lambda x\|=|\lambda|\|x\|$$

假设$x=(1,0,0,1)$，$L_0$范数是2。如果乘以2，得到$(2,0,0,2)$，$L_0$范数依然是2，而不是4。

所以$L_0$范数不是严格意义上的范数。

怎么说呢，CS的人搞机器学习的时候，喜欢用损失函数（想要最小化）；统计的人搞机器学习的时候，一般都是通过式子推导到likelihood（想要最大化）。

所以常常用likelihood取负数或者log likelihood取负数来表示损失函数。

对于这里的poisson regression我们也是这样处理的。它的损失函数就是$-\log (L(\theta))$。具体而言，就是

$$\text{Cost}(\theta)=-\log (L(\theta))=-\sum_{i=1}^n\left(y_i\theta^T x_i -e^{\theta^T} x_i\right)$$

其中$n$是样本个数，$y_i$真实值，$\theta$是$p$元向量，$p$是特征的个数。

在多元回归当中，容忍度是描述多重共线性的。

容忍度就是方差膨胀因子VIF的倒数，根据VIF的定义，化简一下就是一减R方。

不需要满足。因为正态性的假设不是针对y，是指残差$Y_{real}-Y_{pred}$应该满足正态性。

你看到的这个“epoch”不是针对的回归模型，而是针对的随机梯度下降算法的。

算法要收敛到最优，RMSE当然会不当下降。

推荐阅读下自己动手用python写梯度下降

这是把rmse作为loss/cost function，并且learning rate足够小，才能使training rmse逐渐下降。但是testing rmse在达到overfitting后会上升。还有当loss function中加入一些prior/regularization项时，training rmse也会到达拐点后上升。

根据https://stats.stackexchange.com/questions/187824/how-to-model-non-negative-zero-inflated-continuous-data

如y是整数，可以用Zero-inflated Poisson models。

如果是实数：

1. Tobit regression

2. 先做0/非0的分类，然后做非0的回归。

3.Tweedie distributions

简单的方法就是第2个，其余两个在R有package。

我觉得可以先做一个逻辑回归（如果大于0，则label为1；等于0，label为0）。

然后再对所有大于0的数值，做一个线性回归。

最后两者相乘。

对y+1取log，然后再回归呢

线性回归$y=Xw+\epsilon$其中$X$是输入，$y$是输出，$w$是未知参数,$\epsilon$是噪声。所谓的解析解就是最小二乘法($w=argmin_w|Xw-y|^2$)意义下的对$X$的pseudoinverse $w=(X^TX)^{-1}X^Ty$.这有个前提是$X$和$y$服从联合正态分布。如果$X$和$y$不是联合正态分布，那么这个解析解就不是精确的。实际中X和y的分布有可能不知道或者太复杂，不能化简一步得到解析解，所以还是要靠SGD之类的循环算法去逼近最优解。

加一个blog https://wiseodd.github.io/techblog/2017/01/05/bayesian-regression/

还可以去看Murphy, Kevin P. Machine learning: a probabilistic perspective.

--------------------

1.如果X中含有非独立的变量，$rank(\Sigma_{XX})=rank(X^TX)<d$,d是X维度，则$X^TX$不可逆。此时最小二乘无解。

2.如果X中的变量非常相关，$X^TX$的condition number很大，$(X^TX)^{-1}$中误差会被放大，而且会非常大。GD不用求$X^TX$的逆，所以结果比较可信。

我觉得主要有两个原因。

第一个原因：速度快。这个当然是最主要的原因之一了。Zealing的回答里已经提到了，精确解是需要求矩阵乘积和矩阵逆的，这个计算量是比较恐怖的。SGD比这个快多了。

第二个原因：SGD是线上算法（online）。有新数据进入的时候，不需要重新计算，而精确的矩阵求解是做不到线上算法的。

之前好像在某本书里看到“当训练集特征维度数较大或者样本数过多时，使用GD应该可以提高训练的速度”，不知是不是这样。

复杂模型拟合的高阶连续函数应该都是非线性的吧

线性模型

$$y=X\beta + \epsilon$$

$\beta$最小二乘估计的解析解是

$$\hat \beta = (X^TX)^{-1}X^Ty$$

相关解答：最小二乘线性回归的公式推导

http://sofasofa.io/forum_main_post.php?postid=1002094

https://blog.csdn.net/robin_xu_shuai/article/details/75127386

R方可以是负数。如果拟合结果比平均值好，就是正数，比平均值差就是负数。

$$R^2=1-\frac{\sum_{i=1}^ne_i^2}{\sum_{i=1}^n(y_i-\bar y)^2}$$

如果我们用均值表示我们的预测值，那么残差$e_i=y_i-\bar y$，所以此时$R^2=0$。如果拟合度比均值还差，也就是$e_i^2 > (y_i-\bar y)^2$，那么$R^2$自然就是负数了。

当然可以是负的，如果你用平均值取拟合，那么$R^2=0$；如果拟合效果比平均值还差，那么$R^2<0$。

可以。记得在回归时如果不添加常数项有可能出现这种情况。

如果是用的线性回归模型，理论上不会出现$R^2<0$，最差不过接近于0。

但是使用一些其他方法做回归预测，可能会导致$R^2<0$，比如一个非常过拟合的xgboost之类的；或者当算法的优化函数不是RMSE，比如log后的回归，或者huber回归，也是有可能导致$R^2<0$的。

参考Covariance_matrix，线性回归的假设是输入$X$和输出$Y$是联合正态分布。$X$是列向量，每一列为一个数据点。

$ \mu_{X,Y} = \begin{pmatrix} \mu_{X} \\ \mu_{Y} \end{pmatrix} $

$ \Sigma_{X,Y} = \begin{pmatrix} \Sigma_{XX} & \Sigma_{XY} \\ \Sigma_{YX} & \Sigma_{YY} \end{pmatrix} $

线性回归是求conditional mean

$\mu_{Y|X}=\mu_Y+\Sigma_{YX}\Sigma_{XX}^{-1}(X-\mu_X)$

$=\Sigma_{YX}\Sigma_{XX}^{-1}X   + (\mu_Y-\Sigma_{YX}\Sigma_{XX}^{-1}\mu_X)$

$=wX+b$

其中$w=\Sigma_{YX}\Sigma_{XX}^{-1}$,

$b=\mu_Y-w\mu_X$

如果要截距$b=0$，一个充分条件是$\mu_Y=0$，$\mu_X=0$。一般$X$要作normalization，可以保证$\mu_X=0$，还需要让$\mu_Y=0$。几何上意义是通过对原始数据$X,Y$的平移变换(Translation)，让拟合的直线过坐标的原点。

如果$X,Y$没有减去均值，而且强行令$b=0$ ，则偏差为$\mu_Y-w\mu_X$

我觉得你不应该人为的限制截距为0。

首先，如果常数项为0，那么它本质上已经不是一个正确的线性回归了。

其次，如果它本身的特征决定了常数项是0，应该不是由人为所决定，而是通过数据学习得到的。

看了下Comparison of SVM and LS-SVM for Regression，有几点想法：

1.V-SVMR中，只把$\epsilon$带之外的点作为support vector。也就是说$\epsilon$选太大，会扔掉太多training data，所以不准确。而LS_SVMR用到了所有training data。

2. V-SVMR要在所有data中找SV（support vector），是quadratic programming；而LS_SVMR只解一个线性方程组。所以在训练时，LS_SVMR要花更少计算量。

3.在预测时，V_SVMR的只用到SV，而LS_SVMR用到所有数据，所以V_SVMR计算量更小。论文里提到用变通的方法，在LS_SVMR中只留下$|\alpha _k|$大的点，从而减小LS_SVMR在预测时的计算量。

4.理论上并不能证明LS_SVMR比V_SVMR准确，比如在这篇论文的实验，LS_SVMR有更小的MSE。但也在其他论文中见过V_SVMR有更小的MSE，可能和调参有关系。只能说LS_SVMR的训练需要更小的计算量。

------------------------------------------------

5.

. 

LS_SVMR的loss function里第二项就是SSE，而V_SVMR没有。当$\gamma$很大时，loss近似正比于MSE。

所以如果用MSE去衡量V_SMR和LS_SVMR，那么LS_SVMR的MSE很可能会低一点。

《机器学习技法》的第六课讨论过这个问题。

贴个课件的图：

意思大概是，LS-SVMR使用了全部的数据作为support vector，一来计算很慢，二来容易受到噪声、异常值影响。可能是在实际使用中，和准确度相比，这两个因素也很重要吧。

可以通过画图来判断，不过需要先拟合好模型。

比如上面的图看起来就是等方差，下面的图看起来就是异方差。

横轴一般就是某个独立变量。

线性回归的残差

$$e=y-X\hat{\beta}$$

残差的平方和为

$$e^Te=(y-X\hat{\beta})^T(y-X\hat{\beta})=y^Ty-2\hat{\beta}^TX^Ty+\hat{\beta}^TX^TX\hat{\beta}$$

这个是二次型，所以肯定是凸的，对$\beta$求导，可得

$$\frac{d e^Te}{d\hat\beta}=-2X^Ty+2X^TX\hat\beta$$

令导数为$0$，那么

$$X^Ty = X^TX\hat\beta$$

所以

$$\hat\beta=(X^TX)^{-1}X^Ty$$

线性回归的推导也可以从最大似然估计(MLE)入手。

$$y=X\beta+\epsilon$$

根据线性回归的假设，我们知道$\epsilon$是多元正态分布，并且可以认为协方差矩阵为$\sigma_0^2I$，$I$是单位矩阵。

所以变量$y_i$的条件分布是均值为$X_i\beta$并且方差为$\sigma_0^2$的正态分布，所以概率密度函数为

$$F_Y(y_i|\beta, X_i)=\frac{1}{\sqrt{2\pi\sigma_0^2}}e^{-\frac{(y_i-X_i\beta)^2}{2\sigma_0^2}}$$

所以似然函数就是

$$L(\beta, \sigma; X, y)=\prod_{i=1}^N\frac{1}{\sqrt{2\pi\sigma^2}}e^{-\frac{(y_i-X_i\beta)^2}{2\sigma^2}}$$

为了方便计算，两边同时取对数，

$$\log L(\beta, \sigma; X, y)=-N\log\sqrt{2\pi\sigma^2}+\sum_{i=1}^N(-\frac{(y_i-X_i\beta)^2}{2\sigma^2})$$

对于给定的$\sigma$，为了最大化似然函数，我们就要最大化$-(y-X\beta)^2$，也就是最小化$(y-X\beta)^2$。

所以从MLE的角度可以推出我们的损失函数就是差的平方和（最小二乘）。

接下来的推导就是@sasa写的那样了。

$\beta_{LASSO}=\underset{\beta}{\mathrm{argmin}}||y-X\beta||^2+\lambda|\beta|_1$

当$\lambda=0$时，没有限制（shrinkage），此时无偏；

当$\lambda$增大时，bias增大，variance减小；

当$\lambda$减小时，bias减小，variance增大。

看到有个说法（我没验证过），当添加先验知识从而减小自由度（degree of freedom）时，估计就是有偏的。LASSO的自由度是非零参数个数。

一个简单的例子，当y和X都是高斯分布，第一个项（似然函数项）得到$\beta$的估计是无偏，再加一个有偏的第二项（先验知识项，会迫使所有$\beta=0$，显然是有偏的），最后结果就是有偏的。

bias是从estimator的角度来说的。

其实LASSO和Ridge，甚至说所有有penalty的回归都是biased，不然会违背高斯-马尔可夫定理。

LOWESS (locally weighted scatterplot smoothing) 本质上是一种回归，类似于分段回归。比较常见的应用是，当有一群散点图时，我们需要一条曲线来表示整个散点图的整体趋势。这时LOWESS就可以做到这一点。

与传统的回归模型不同，LOWESS是在数据的局部进行回归。基本思想是对于某个数据点$x_0$，找到它的k近邻，然后利用它的k近邻的数值，在$x_0$的邻域上做线性回归（或者多项式回归），此时回归模型的损失函数是加权的损失函数（与$x_0$距离越远，权重越小）。对多个数据点做本地回归，把多个回归模型的结果连接起来就得到类似上图的曲线。

你所谓的分割value也是根据metric来的。

可以参考一下这个对于数值特征，决策树是如何决定分割阈值的？

用左子节点的均值作为左子节点的预测值，用右子节点的均值作为右子节点的预测值，对比真实值，可以得到MSE或者MAE

遍历所有特征的所有分割点，找到MSE或者MAE最小的那个分割点

最小二乘是求残差$e=Xw-y$的最大似然，也就是求最小-log似然。为简化问题，假设$e$是样本独立的nx1高斯分布,n是数据个数，$e$的nxn协方差矩阵$\Sigma_{ee}$是对角线矩阵，每个数据点在loss function中的权重是其残差的方差的倒数$1/var(e_i)$。换句话说，$e_i$方差越大，越不可信，其权重越小。而$e_i$方差的估计由实际问题的统计模型决定，我觉得是最小二乘框架中最重要，最体现对实际问题理解程度的地方。

$$\Sigma_{ee}=\begin{bmatrix} var(e_1) & & \\ & \ddots & \\ & & var(e_n) \end{bmatrix}$$

$$\Sigma_{ee}^{-1}=\begin{bmatrix} 1/var(e_1) & & \\ & \ddots & \\ & & 1/var(e_n) \end{bmatrix}$$

数据点加权的-logloss function是

$$-\log L(w)=(Xw-y)^T\Sigma_{ee}^{-1}(Xw-y)$$

因为$\Sigma_{ee}^{-1}$是对角线矩阵，有

$$-\log L(w)=\Sigma_{ee}^{-1}(Xw-y)^T(Xw-y)$$

其中$w$是待求参数，$X$是输入数据，$y$是输出，$\Sigma_{ee}$是输出的covariance matrix。每个数据对应的-logloss的权重是$1/var(e_i)$。

没有权重的线性回归

$$y_i=X_i\beta + \beta_0+\sigma_i$$

$\sigma_i$是拟合的残差。根据线性回归的基本假设，残差是相等的，所以为了确保这一点，我们给每个样本增加权重

$$w_i=\frac{1}{\sigma_i^2}$$

LASSO自由度是非零参数个数，或者应该说非零参数个数是LASSO自由度的无偏估计。论文在这里。摘要里有结论。

-----------

自由度是是经过算法处理后数据中独立变量的个数。LASSO的非零参数是还活着，还能变化的参数，其余参数都被限制为零。这些非零参数个数是自由度的一个估计。就像鱼缸养很多鱼，最后活下来鱼的数量就是鱼缸养鱼容量的一个估计。论文的贡献是证明这个估计是无偏的。

可以用线性回归的方式，比如你有$m$个模型，它们的结果分别是$Y_1,Y_2,\cdots,Y_m$，真实的结果为$Y$，你的线性组合就是

$$Y=w_0+w_1Y_1+w_2Y_2+\cdots+w_mY_m$$

回归系数$w_i$就是第$i$个模型的权重

有时候可以用一个简单的随机森林把预测结果汇总起来做个最终预测

这个可以根据岭回归的损失函数推导的，

$$Loss(\beta)=(Y-X\beta)^T(Y-X\beta)+\lambda \beta^T\beta$$

$\lambda$是正则项的系数

这个是凸问题，所以导数为0时取得解

$$\frac{\partial Loss}{\partial \beta}=2X^TX\beta-2X^TY+2\lambda I\beta=0$$

所以

$$(X^TX+\lambda I)\beta=X^TY$$

解析解为

$$\beta=(X^TX+\lambda I)^{-1}X^TY$$

Ridge是有解析解的，假如正则项的惩罚系数为$\lambda$，那么解析为

$$\hat{\beta}=(X^TX+\lambda I)^{−1}X^Ty$$

里面的$I$是单位矩阵

LASSO没有一步到位的解析解。

$$l(\beta)=(X\beta-y)^T(X\beta-y)+\lambda|\beta|_1$$

$$\frac{\partial l}{\partial \beta}=X^T(X\beta-y)+\lambda sign(\beta)=0$$

$$\beta=(X^TX+\lambda sign(.))^{-1}X^Ty$$

因为右面和当前估计值$sign(\beta)$有关，只能用迭代的方法求出。不过求矩阵逆的计算量很大，还不如用gradient descent类方法。

---------------

不好意思，看成LASSO的了。Ridge的解析解应该是strong.man写的。

如果有两个变量$X$和$Y$，并且它们的尺度都是-1到1的，那么$X$和$Y$的皮尔逊相关系数就是$X$和$Y$做一元线性回归的回归系数（解释），从这点我们可以看明显的看出，皮尔逊相关系数就是描述的线性关系。

另外一个反例是，如果$Y=X^2$，并且$X$的范围是-1到1，尽管它们是有关联的，但是它们的皮尔逊系数为0，因为它们并没有任何线性关系。

"If a relationship between two variables is not linear, the rate of increase or decrease can change as one variable changes, causing a "curved pattern" in the data. This curved trend might be better modeled by a nonlinear function, such as a quadratic or cubic function, or be transformed to make it linear. Plot 4 shows a strong relationship between two variables. However, because the relationship is not linear, the Pearson correlation coefficient is only +0.244."

source: https://support.minitab.com/en-us/minitab/18/help-and-how-to/statistics/basic-statistics/supporting-topics/basics/linear-nonlinear-and-monotonic-relationships/

ARCH是在AR模型上的提升。AR就是AutoRegression，自回归。自回归中残差是独立同分布的，也就是等方差的。

现实情况下残差未必如此，残差在某些条件下未必是等方差的，所以就有了AutoRegressive Conditional Heteroskedasticity model，ARCH自回归条件异方差模型。

想要进一步了解的话，可以看看这个教程

你的问题和下面这个问题应该是一样的

回归中自变量和因变量的相关系数和回归系数（斜率）有什么关系？

你要把x和y都正规化一下，然后它们就应该是一样的了

如何是线性回归的话应该是一样的。

比如有有两组数a组和b组，线性回归的自变量只有一个，是一个indicator，数据是否来自a组，应变量就是数据的值，这种线性回归下回归系数的p值应该是和t-test下的p值一样的。

如果每个自变量之间是独立的，那么x2的p值应该是不变的。

可是往往自变量之间是有关联的，所以会造成x1的变化导致x2的p值也跟着变。

不一定一样，按照实际情况来说，很可能不一样。因为逐步选择特征的方法是贪婪算法，贪婪算法一般都不会达到全局最优。向前和向后应该会分别到达各自的局部最优。

效果不一样。

我感觉forward stepwise好一点。

随机森林回归和决策树回归的本质是一样的，而你说到的ID3,C4.5,GINI等分类原则和回归无关的。

随机森林回归的原则一般就是MSE或者MAE。

随机森林的最终结果就是所有树的回归结果的均值。

相关的讨论：决策回归树的分割节点

很显然是凸函数。@cabbage 推导了如何得到目标函数。

对于线性回归，目标函数很简单，就是平方损失

假设是单变量的线性回归$y=ax+b+\epsilon$，那么损失函数

$$L=\frac{1}{n}\sum_{i=1}^n(ax_i+b-y_i)^2$$

要证明这个损失函数是关于$a,b$的凸函数，我们就只需要求二阶偏导

$$\frac{\partial L}{\partial a}=\frac{1}{n}\sum_{i=1}^n2x_i(ax_i+b-y_i)$$

$$\frac{\partial^2 L}{\partial a^2}=\frac{1}{n}\sum_{i=1}^n 2x_i^2\geq 0$$

所以$L$关于$a$的二阶导是非负的

$$\frac{\partial L}{\partial b}=\frac{1}{n}\sum_{i=1}^n2(ax_i+b-y_i)$$

$$\frac{\partial^2 L}{\partial b^2}=\frac{1}{n}\sum_{i=1}^n 2 \gt 0$$

所以$L$关于$b$的二阶导也是非负的

所以损失函数$L$是关于$a$和$b$的凸函数。

对于多变量线性回归也是一样的推导方法，可以证明是凸函数。

我们用$\epsilon^{(i)}$代表误差，则预测函数可以写为

$$y^{(i)} = \theta^Tx^{(i)}+\epsilon^{(i)}$$

$$\epsilon^{(i)} = y^{(i)} -\theta^Tx^{(i)}$$

其中，我们假设误差是随机分布的，均值为0，服从高斯分布$N(0,\sigma)$，因为根据中心极限定理，服从高斯分布也是对误差项分布的合理猜想。

所以

$$P(y^{(i)}|x^{(i)}; θ) = \frac{1}{\sqrt{2\pi}\sigma}\text{exp}(- \frac{(y^{(i)}-\theta^Tx^{(i)})^2}{2\sigma^2})$$

$P(y^{(i)}|x^{(i)}; θ)$表示：在$\theta$为给定的参数的情况下，概率$y^{(i)} $以$x^{(i)} $为随机变量的概率分布，注意$\theta$不是随机变量。

由于$\epsilon^{(i)}$是独立的同分布`IID：independentlyidentically distribution`，所以以$\theta$为变量的似然函数为：

$$L(θ)=L(θ;X,Y)=p(Y|X;θ) = \prod_{i=1}^{m}\frac{1}{\sqrt{2\pi}\sigma}\text{exp}(- \frac{(y^{(i)}-\theta^Tx^{(i)})^2}{2\sigma^2})$$

对 $L(θ) $取对数有：

$$l(\theta)=\log L(\theta) = \log\prod_{i=1}^{m}\frac{1}{\sqrt{2\pi}\sigma}\text{exp}(- \frac{(y^{(i)}-\theta^Tx^{(i)})^2}{2\sigma^2})$$

$$= m\log\frac{1}{\sqrt{2\pi}\sigma} - \frac1{2\sigma^2}\sum_{i=1}^{m}(y^{(i)}-\theta^Tx^{(i)})^2$$

最大化$l(\theta)$即是最小化$\frac1{2\sigma^2}\sum_{i=1}^{m}(y^{(i)}-\theta^Tx^{(i)})^2$，这样就是`loss function`

基本上就是和分布式随机梯度下降一样。

看这个问题随机梯度下降(SGD)可以被并行计算吗？

我感觉你的预测值就是平均值，没有去进行收敛的。

问题可能是

1）你的输入值没有进行标准化处理

2）你的学习率设置不当，你可以减小试试

RMSE是均方误差根，RMSLE是均方对数误差根。L是取对数的意思。

$$RMSLE=\sqrt{\frac{\sum_{i=1}^n (\log(y_i+1)-\log(\hat y_i +1))^2}{n}}$$

有时候会针对一些长尾分布的y，我们会用RMSLE。RMSE是对称的，RMSLE对高估的惩罚比低估更大。可以根据自己的需求酌情选择。

我理解的门槛回归(threshold regression)就是分段回归，piecewise regression。

找一些临界点把原数据分成几个部分，然后每个部分单独拟合出一个回归模型，比如下图，就分成了两段

这两个模型的系数相差这么小，如果非要判断其显著不同，只能设定非常高的显著水平了。。。

暂且用假设检验的方法测试一下：

两个模型都是线性模型，但是并不知道x1, x2的分布情况，那我假设x1，x2均为非正态分布的随机数，那两个模型应该不符合正态分布，同时按照楼主的要求，既然显著不同，那假设其方差也不齐，应该选择秩和检验，零假设选择为两个模型相同，即来自同一分布。

import numpy as np
from scipy.stats import wilcoxon

x1 = np.random.randint(1, 100, 10)
# array([60, 59, 22, 18, 45, 37, 17, 89, 37, 89])
x2 = np.random.randint(1, 100, 10)
# array([24, 47, 95, 48, 73, 48, 62, 44, 34, 88])
x = 0.95*x1 + 0.85*x2 + 1.01
y = 0.96*x1 + 0.84*x2 + 0.99
wilcoxon(x, y)
# WilcoxonResult(statistic=20.0, pvalue=0.4445867389117455)

我测试了几次，pvalue都远远大于0.1，通常显著水平都是设置在0.05，有时也会根据需要取0.1或0.01，但是要设置成远大于0.1的情况，并不常见。除非强行设置极大的显著水平，否则我觉得无法判断这两个模型显著不同。

误差项不符合独立性的情况通常会发生在时间序列模型里，这种情况你的模型通常也是需要提高精度的。

误差项不符合正态分布的话，你可以通过对因变量进行非线性变换，或者使用不同link function的广义线性回归。

如果是误差项的方差不同的话，可以考虑加权。参考如何计算加权最小二乘法的样本权重？

将fit = glm(wage.logit ~ bs(age.logit, df = 3), family = binomial(link = "logit")) 改成fit = lrm(wage.logit ~ bs(age.logit, df = 3), family = binomial(link = "logit"))就可以了

比如训练集是$\{(X_1,y_1),(X_2,y_2),(X_3,y_3),\ldots,(X_n,y_n)\}$

训练得到一个线性回归模型$\hat y = X\beta_1$。按照boosting的思想，我们要用真实值减去模型的预测值，得到第二轮训练的真实标签，所以第二轮的训练集是$\{(X_1,\epsilon_1),(X_2,\epsilon_2),(X_3,\epsilon_3),\ldots,(X_n,\epsilon_n)\}$。这时你可以再训练一个线性回归得到$\hat \epsilon = X\beta_2 $。

如果此时停止迭代，那么你得到的最终模型就是应该是$X(\beta_1+\beta_2)$，本质上还是线性回归模型。这个新模型是可以最小化平方误差的。而$X\beta_1 $也是最小化平方误差的，所以两者是等价的。所以boosting没有意义。实际上$\beta_2=0$。

boosting是bagging的进阶版，实际上连bagging对线性回归都没有效果，何况是boosting呢？

线性回归和bagging的讨论可以看这个线性回归的bagging

第一个问题可以参考一下Lasso和岭回归的正则项包含截距（常数项）吗？

第二个问题：所有的最小二乘回归方程都会经过$(\bar x, \bar y)$这一个点，如果你对所有的自变量以及y进行中心化，也就是说自变量的均值$\bar x=0$，因变量的均值$\bar y=0$。一个线（或者超平面）要经过$(0,0)$这个点，那么截距一定是0。

快速排除的话，有几个方法：

1）可以考虑单因子模型，单因子预测效果最差的特征可以排除。

2）看特征的缺失值占比，缺失值占比过高的，可以直接剔除

3）看特征方差，特征是常数或者接近为常数的，可以直接剔除

4）做LASSO，被LASSO剔除的变量可以不考虑

5）重复的特征，比如两个特征几乎一样，也可以直接剔除

如果复杂一点、精确一点的话，可以做逐步特征选择（向前或者向后），AIC、BIC，feature importance，VIF方法等等，做交叉验证来看哪些特征需要剔除掉。