理论知识

理论模型 y=a+bx+ε
X是解释变量，又称为自变量，它是确定性变量，是可以控制的。是已知的。

Y是被解释变量，又称因变量，它是一个随机性变量。是已知的。

a,b是待定的参数。是未知的。

ε 是误差项。

保证最小二乘估计是最佳无偏估计。
（1）正态性假设：要求总体误差项服从正态分布
（2）零均值性假设：在自变量取一定值得条件下，其总体各误差项的条件平均值为零。
（3）零方差性假设：在自变量取一定值得条件下，其总体各误差项的条件方差为一常数。
（4）独立性假设：误差项之间相互独立，误差项和自变量之间相互独立。

三种线性回归

正相关、负相关、非线性相关（非线性相关不代表不相关，也可能有二次函数关系、指数关系等其他相关关系，只是非线性 相关）

协方差

cov(X,Y)=E[(X-E(X))(Y-E(Y))]
E 是期望值
协方差反映X和Y的相关性，及相关性程度，小于0为负相关，大于0为正相关。缺点是会受变量变化幅度的影响。

相关系数

r=cov(X,Y)/(σxσy)
X和Y的协方差除以X的标准差与Y的标准差
相关系数，是标准化后的协方差，是线性分析中比较常用的，r=1正线性相关，r=-1负线性相关，r=0完全无关。且r的绝对值越大表示相关性越强。（注：相关关系≠因果关系）

误差

image.png
image.png

似然函数： 由数据去推参数等于什么，让我们的数据跟真实数据接近的可能性越大越好

image.png
损失函数

求：为什么用最小二乘法
误差函数的由来，它是由极大似然估计得到的。我们要最小化这个误差函数，即最小化1/2倍的每个样本预测值与真值差值的平方和

最小二乘法下的参数最优解

当我们有了损失函数，希望损失函数最小的时候，有两种办法： 最小二乘法；梯度下降

最小二乘法有很多种推理方法，结果就是直接把参数求出来
梯度下降就是一点一点的用学习率迭代，找到最小值

求目标函数:最小二乘法

image.png

评估模型的好坏

image.png

梯度下降

image.png

随机梯度下降

stochastic gradient descent

小批量随机梯度下降

mini batch SGD

优缺点

优点：结果易于理解，计算上不复杂。
缺点：对非线性数据拟合不好。线性回归受噪声影响较大
适用数据类型：数值型和标称型数据。

代码实现：梯度下降

import numpy as np
import pylab

def compute_error(b,m,data):

    totalError = 0
    #Two ways to implement this
    #first way
    # for i in range(0,len(data)):
    #     x = data[i,0]
    #     y = data[i,1]
    #
    #     totalError += (y-(m*x+b))**2

    #second way
    x = data[:,0]
    y = data[:,1]
    totalError = (y-m*x-b)**2
    totalError = np.sum(totalError,axis=0)

    return totalError/float(len(data))

def optimizer(data,starting_b,starting_m,learning_rate,num_iter):
    b = starting_b
    m = starting_m

    #gradient descent
    for i in range(num_iter):
        #update b and m with the new more accurate b and m by performing
        # thie gradient step
        b,m =compute_gradient(b,m,data,learning_rate)
        if i%100==0:
            print 'iter {0}:error={1}'.format(i,compute_error(b,m,data))
    return [b,m]

def compute_gradient(b_current,m_current,data ,learning_rate):

    b_gradient = 0
    m_gradient = 0

    N = float(len(data))
    #Two ways to implement this
    #first way
    # for i in range(0,len(data)):
    #     x = data[i,0]
    #     y = data[i,1]
    #
    #     #computing partial derivations of our error function
    #     #b_gradient = -(2/N)*sum((y-(m*x+b))^2)
    #     #m_gradient = -(2/N)*sum(x*(y-(m*x+b))^2)
    #     b_gradient += -(2/N)*(y-((m_current*x)+b_current))
    #     m_gradient += -(2/N) * x * (y-((m_current*x)+b_current))

    #Vectorization implementation
    x = data[:,0]
    y = data[:,1]
    b_gradient = -(2/N)*(y-m_current*x-b_current)
    b_gradient = np.sum(b_gradient,axis=0)
    m_gradient = -(2/N)*x*(y-m_current*x-b_current)
    m_gradient = np.sum(m_gradient,axis=0)
        #update our b and m values using out partial derivations

    new_b = b_current - (learning_rate * b_gradient)
    new_m = m_current - (learning_rate * m_gradient)
    return [new_b,new_m]


def plot_data(data,b,m):

    #plottting
    x = data[:,0]
    y = data[:,1]
    y_predict = m*x+b
    pylab.plot(x,y,'o')
    pylab.plot(x,y_predict,'k-')
    pylab.show()


def Linear_regression():
    # get train data
    data =np.loadtxt('data.csv',delimiter=',')

    #define hyperparamters
    #learning_rate is used for update gradient
    #defint the number that will iteration
    # define  y =mx+b
    learning_rate = 0.001
    initial_b =0.0
    initial_m = 0.0
    num_iter = 1000

    #train model
    #print b m error
    print 'initial variables:\n initial_b = {0}\n intial_m = {1}\n error of begin = {2} \n'\
        .format(initial_b,initial_m,compute_error(initial_b,initial_m,data))

    #optimizing b and m
    [b ,m] = optimizer(data,initial_b,initial_m,learning_rate,num_iter)

    #print final b m error
    print 'final formula parmaters:\n b = {1}\n m={2}\n error of end = {3} \n'.format(num_iter,b,m,compute_error(b,m,data))

    #plot result
    plot_data(data,b,m)

if __name__ =='__main__':

    Linear_regression()

代码实现：最小二乘法，以波士顿房价为例子

import numpy as np
import copy
from sklearn.datasets import load_boston#导入博士顿房价数据集
 
 
class LinerRegression:
    M_x = []  #
    M_y = []  #
    M_theta = []  # 参数向量
    trained = False
 
    def __init__(self):
        pass
 
    def regression(self, data, target):
        self.M_x = np.mat(data)
        # 每个向量添加一个分量1，用来对应系数θ0
        fenliang = np.ones((len(data), 1))
        self.M_x = np.hstack((self.M_x, fenliang))
        self.M_y = np.mat(target)
        M_x_T = self.M_x.T  # 计算X矩阵的转置矩阵
        self.M_theta = (M_x_T * self.M_x).I * M_x_T * self.M_y.T# 通过最小二乘法计算出参数向量
        self.trained = True
 
    def predict(self, vec):
        if not self.trained:
            print("You haven't finished the regression!")
            return
        M_vec = np.mat(vec)
        fenliang = np.ones((len(vec), 1))
        M_vec = np.hstack((M_vec, fenliang))
        estimate = np.matmul(M_vec,self.M_theta)
        return estimate
 
 
if __name__ == '__main__':
    # 从sklearn的数据集中获取相关向量数据集data和房价数据集target
    data, target = load_boston(return_X_y=True)
    lr = LinerRegression()
    lr.regression(data, target)
    # 提取一批样例观察一下拟合效果
    test = data[::51]
    M_test = np.mat(test)
    real = target[::51]#实际数值real
    estimate=np.array(lr.predict(M_test))#回归预测数值estimate
    #打印结果
    for i in range(len(test)):
        print("实际值:",real[i]," 估计值:",estimate[i,0])

reference：
https://zhuanlan.zhihu.com/p/35887329
https://blog.csdn.net/sxf1061926959/article/details/66976356
https://blog.csdn.net/lc013/article/details/55002463
https://blog.csdn.net/back_to_dream/article/details/51361431
https://blog.csdn.net/qq_32864683/article/details/80368135
https://blog.csdn.net/qsczse943062710/article/details/55508314