线性回归

基础假设

There is a linear relationship between the independent variables(X) and the dependent variables (y)
Independence: Independence assumes that there is no relationship or correlation between the errors (residuals) of different observations.
Normality: The residuals of the linear regression model are assumed to be normally distributed.
Homoscedasticity: Homoscedasticity assumes that the variability of the errors (residuals) is constant across all levels of the independent variables.
No Multi-collinearity between features

1. 最小二乘法

最小二乘法(ordinary least squares)对线性回归的假设是：误差（residual）符合正态分布。

最小二乘法的理念是承认观测误差的存在，误差是围绕真值上下波动的，从几何角度看，当平方误差最小时即为真值。具体应用时，算法的目标是要求线性方程中未知权重w和b，当经验风险最小化时的权重认为是最优的。为了简化，我们经常在推导中忽略b，因为b可以认为存在一个常数列x=1对应的权重w的一个分量，不必单独另求。

以下向量定义为列向量，注意不同步骤需要清楚定义谁被看作是谁的函数

$y=w^Tx +b$

这里，有一个转置T，是为了求两个向量的点积。

$w=\mathop {argmin} _{w}Loss(w)$

Loss function:

$L(w)=\frac{1}{n}\sum_{i=1}^{n}{(y-\hat{y})^2}$

从标准量中计算其梯度(只有一个样本，多样本加连加符号):

$\frac{\partial{L(w)}}{\partial{w_j}} = {(y-\hat{y})}{x_j}$

梯度是让损失函数下降的方向。先看单个参数，就决定了下一步的大小(学习率还会拦一道)和方向(也就是正负)。对多个特征需要的多参数，主要是决定了不同参数大小的对比，也就成了多向量空间中的方向。对于多个样本来说，可以认为每个样本都会让参数"优化"移动一次。不同参数的移动大小，可以说决定了这个参数的重要程度。因为本身需要norm，让这个移动更纯粹一些。

其矩阵形式如下:

$L= (y-XW)^{T}(y-XW)$

$L=y^Ty-y^TXW-W^TX^Ty+y^Ty+W^TX^TXW$

继续化简需要用到矩阵求导的几个规律：

$\frac{\partial{L}}{\partial{w}}=2X^TXw-2X^Ty$

损失函数为凸函数，可知道全局最优解的取值为导数等于0处的极值
解析解的存在需要满足满秩矩阵。加L2 norm之后有了先验约束，则不需要
梯度下降的优化方法：梯度决定了参数的变化方向。模型中每一个参数都看做是损失函数的自变量。结合逻辑回归、深度学习的前向反向传播综合理解

解析解-> 导数为0:

$W=(X^TX)^{-1}X^Ty$

2. 最大似然法 Maximum likelihood estimation

最小二乘与最大似然的等价，同属于频率派(Frequentists)；最大后验则属于贝叶斯派(Bayesians)。

什么是似然（likelihood）？在似然中，所有数据对（x,y）集合，被看成是一个随机变量X，但X已经被观察到了，已经固定。以最常见的抛硬币模型来说，硬币出现正面或反面就是一个随机变量，正常硬币出现正面的概率为1/2。而似然是给定一个硬币，我们抛了10次之后出现了1次正面，随机变量出现1/10正面这个事件已经被观察到了并固定以后，这个硬币是正常硬币的可能性就被称为似然。

$w,b= argmax(logP(x|w))$

在线性回归中，贝叶斯公式可以写成:

$P(w|X)=\frac{P(X|W).P(w)}{P(X)}$

P(X|w)是似然函数，
P(w)被称为先验(prior)，是发生在我们观察样本之前，例如硬币抛出之前对硬币的假设。
P(w|X)就称为后验概率(posterior)

最大似然：

$w_{MLE}=\mathop{argmax}_w P(x|w)$

最大似然需要假设分布，对于一般线性回归的假设即是: 误差是一个高斯分布。(注意区分后验概率中的w本身是一个高斯分布)

关于变量x的高斯分布公式：

$P(x)=\frac{1}{\sqrt{2\pi}\sigma}exp(-\frac{(x-\mu)^2}{2\sigma^2})$

误差满足高斯分布，而误差是真值和预测的差值，所以：

$P(\hat{y}-xw) \sim N(0，\sigma^2)$

由于误差是均值为0，方差$\sigma^2$的高斯分布，则y的预测值是满足均值为wx，方差为$\sigma^2$的高斯分布。在似然中的x其实代表的是一个样本，即（x,y)数据对，而在特征已知的情况下，似然函数可以进一步表示为：

$P(\hat{y}|x,w)$

因此，似然也可以表达成一个高斯分布

$P(y|x,w) \sim N(w^Tx,\sigma^2)$

可以转化为：

$P(y|x,w)=\frac{1}{\sqrt{2\pi}\sigma}exp(-\frac{(y-w^Tx)^2}{2\sigma^2})$

3. 最大后验法与正则 Maximum a posteriori estimation

贝叶斯认为变量w也是一个随机变量，随机变量自带一个概率分布。因此在观察数据之前，我们可以提前根据经验预估一个w的分布，这个预估就是先验，并成为一个对参数的约束，即为最大后验法。

$w_{MAP}=\mathop{argmax}_w P(x|w)P(w)$

L1与L2理论涉及到朗格朗日优化。

$w=argmax(P(w|X))=argmax(P(X|w).P(w))$

从损失函数的角度，加了一个针对w大小的惩罚。

L1正则

假设w的先验满足拉普拉斯分布
在损失函数中引入了绝对值，绝对值的导数是不连续的。可采用坐标下降法优化，即每次更新一个w的分量

L2正则

假设w的先验服从零均值正态分布

解析解变为

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

4. 贝叶斯线性回归与高斯过程回归

贝叶斯线性回归要预测出一个关于y的分布来。稍微拓展一下贝叶斯线性回归到高斯过程回归(Gaussian Process Regression)

贝叶斯线性回归最大的不同在于预测的不是关于w的点估计，而是w的分布，即w的期望和协方差。

5. 问答

pros:
- explainable method
- interpretable results by its output coefficients
- faster to train than other machine learning models
cons:
- assumes linearity between inputs and outputs
- sensitive to outliers
- can under fit with small, high-dimensional data
为什么需要对数值类型特征进行归一化
- 使用梯度下降优化的模型，归一化容易更快通过梯度下降找到最优解，包括线性回归、逻辑回归、支持向量机、神经网络。
如何判断是否该用线性回归模型
- EDA画图
- 检查预测值与真实值的residual残差是否为高斯分布
Lasso regression如何导致参数稀疏
- 图像切点
- 朗格朗日优化
异常点如何处理
- 检测并去除
- Huber Loss is a loss function that combines MSE and MAE. It is less sensitive to outliers than a squared error loss and has the advantage of not needing to remove outlier data from the dataset
共线性问题 Multi-collinearity
- 自变量之间由于存在高度相关关系而使模型的权重参数估计失真或难以估计准确的一种特性
- 共线性不影响模型的预测而是影响对模型的解释
- 可以尝试l2正则；可以使用PCA，将特征转为独立的变量
- feature selection, ridge regression, PCA

6. 代码

梯度下降

对多元函数的参数求∂偏导数，把求得的各个参数的偏导数以向量的形式写出来，就是梯度
从几何意义上讲，梯度就是函数变化增加最快的方向。减梯度就是朝着减小的方向

def update_weights(m, b, x, y, learning_rate):
    m_deriv = 0
    b_deriv = 0
    N = len(x)
    for i in range(N):
        m_deriv += -2 * x[i] * (y[i] - (m * x[i] + b))
        b_deriv += -2 * (y[i] - (m * x[i] + b))
    
    m -= (m_deriv / float(N)) * learning_rate
    b -= (b_deriv / float(N)) * learning_rate    
    return m, b

牛顿法

在现有极小点估计值的附近对 f(x) 做二阶泰勒展开，进而找到极小点的下一个估计值

线性回归

# 2018.01.26 我学习机器学习的第一天，learning_rate too large cause the gradient explosion
import numpy as np
import matplotlib.pyplot as plt
np.random.seed(1)

class LinearRegression(object):
    def __init__(self):
        self.w = None
        self.b = None

    def fit(self, x, y):
        batch_size, feature_size = x.shape
        self.w = np.random.random((feature_size, 1))
        self.b = np.random.random(1)
        self._gradient_descent(x, y)

    def eval(self):
        pass

    def predict(self, x):
        y_hat = np.dot(x, self.w) + self.b
        return y_hat

    def plot(self, x, y, y_pred):
        plt.plot(x, y, 'r-')
        plt.plot(x, y_pred, 'b.')
        plt.show()

    def _gradient_descent(self, x, y, learning_rate=0.0001, epoch=100):
        for i in range(epoch):
            y_hat = self.predict(x)  # => batch_size*1
            loss = np.dot((y_hat-y).T,(y_hat-y))
            w_gradient = np.dot(x.T, x).dot(self.w) - np.dot(x.T, (y-self.b))
            b_gradient = self.b-np.mean(y-np.dot(x,self.w))
            self.w = self.w - learning_rate * w_gradient
            self.b=self.b-learning_rate*b_gradient
            print('step:', i, 'Loss:', loss[0][0])

    def __str__(self):
        return 'weights\t:%s\n bias\t:%f\n' % (self.w, self.b)

def generate_data():
    '''generate the examples to implement linear regression by numpy'''
    x = np.linspace(1, 30, num=30).reshape(-1, 3)
    y = np.linspace(1, 10, num=10) + np.random.random(10)
    y = y.reshape(-1, 1)
    return x,y

if __name__ == "__main__":
    x, y = generate_data()
    lr = LinearRegression()
    lr.fit(x, y)
    y_hat = lr.predict(x)
    lr.plot(x, y, y_hat)

参考

Previous评价指标 Next逻辑回归

Last updated 14 days ago