精华内容
下载资源
问答
  • import math x_train_set = x[: math.floor(len(x) * 0.8), :] y_train_set = y[: math.floor(len(y) * 0.8), :] x_validation = x[math.floor(len(x) * 0.8): , :] y_validation = y[math.floor(len(y) * 0.8): , :...
  • 多元线性回归求解过程 解析解求解

    万次阅读 多人点赞 2018-10-30 17:21:00
    常用算法一 多元线性回归详解2(解析解求解多元线性回归)  上一篇讲到什么是多元线性回归以及多元线性回归的推导过程详解,本章我们一起来看如何求得最优解,就是我们得到了多元线性回归到损失函数就是最小二乘公式...

    常用算法一 多元线性回归详解2(解析解求解多元线性回归)

            上一篇讲到什么是多元线性回归以及多元线性回归的推导过程详解,本章我们一起来看如何求得最优解,就是我们得到了多元线性回归到损失函数就是最小二乘公式,那么如何利用最小二乘公式求得最优解。对多元线性回归到推导过程没有概念的同学欢迎查看上一篇文章::常用算法一 多元线性回归详解1(推导过程)

     

    求解多元线性回归        

            多元线性回归常用的求解方法有两种:

                      1-解析解求解法

                      2-梯度下降法求解

             本章我们来看多元线性回归的解析解求解法。

    解析解求解法

             说到解析解求解,很多同学都已经忘记了什么事解析解。解析解就是指通过公式就可以求得到方程的解。我们只需要方程的参数带入到公式中,计算公式结果就可以得到方程的解,而不用一步一步化简求解。比如我们初中学的一元二次方程的解细节是。是不是豁然开朗,原来就是你小子。

             想要用解析解来求解最小二乘函数,那我们首先得知道他的解析解是啥。

    a.求得最小二乘公式的解析解。

             这里要用到上一章讲到的知识点,求一个函数在某一点上的导数,就是求在这个函数的图像上,过这一点所做切线的斜率。这一点的导函数就是切线的函数。一个二次函数的图像是一个抛物线,那想想一下,通过图像的顶点所做的切线是一条怎样的直线。应该是一条与x轴平行的直线,此时这条直线的斜率为0.函数图像的顶点就是函数的解,也就是说,我们通过函数的解这一点来做切线,切线的斜率就是0.

             那我们反过来利用一下刚刚总结出的结论。如果我找到了函数图像上切线为0的点,是不是找到了函数的解?切线是什么?对函数上某一点求导就等于通过这一点在函数图像上做切线,作出的切线就是求导得到的导函数的图像,切线的函数就是对函数求导所得到的导函数,那我们只要找到导函数为0对点,是不是就得到了图像的解?(这一段一定要理解。多读几遍)

             所以,我们可以通过对最小二乘函数求导,让导函数为0时的结果,就是最小二乘的解。求导过程如下:

            首先对最小二乘进行变形,变为矩阵表达形式:

            

           展开矩阵函数:

                  

        展开之后我们对J(θ)求导并令导数等于0:

                           

       最终求的解析解为:θ=(X^{T}X)^{-1}X^{T}Y

    2-解析解的代码实现

      手动实现: 

    
    import numpy as np
    import matplotlib.pyplot as plt
    from bz2 import __author__
    
    #设置随机种子
    seed = np.random.seed(100)
    
    #构造一个100行1列到矩阵。矩阵数值生成用rand,得到到数字是0-1到均匀分布到小数。
    X = 2 * np.random.rand(100,1)   #最终得到到是0-2均匀分布到小数组成到100行1列到矩阵。这一步构建列X1(训练集数据)
    #构建y和x的关系。 np.random.randn(100,1)是构建的符合高斯分布(正态分布)的100行一列的随机数。相当于给每个y增加列一个波动值。
    y= 4 + 3 * X + np.random.randn(100,1)
    
    #将两个矩阵组合成一个矩阵。得到的X_b是100行2列的矩阵。其中第一列全都是1.
    X_b = np.c_[np.ones((100,1)),X]
    
    #解析解求theta到最优解
    theta_best = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)
    # print(theta_best)
    # 生成两个新的数据点,得到的是两个x1的值
    X_new = np.array([[0],[2]])
    
    # 填充x0的值,两个1
    X_new_b = np.c_[(np.ones((2,1))),X_new]
    
    print (X_new_b)
    
    # 用求得的theata和构建的预测点X_new_b相乘,得到yhat
    y_predice = X_new_b.dot(theta_best)
    print(y_predice)
    # 画出预测函数的图像,r-表示为用红色的线
    plt.plot(X_new,y_predice,'r-')
    # 画出已知数据X和掺杂了误差的y,用蓝色的点表示
    plt.plot(X,y,'b.')
    # 建立坐标轴
    plt.axis([0,2,0,15,])
    
    plt.show()
    

      利用sklearn包实现

    from sklearn.linear_model import  LinearRegression
    import numpy as np
    import matplotlib.pyplot as plt
    # 解析解求线性回归
    
    
    # 手动构建数据集和y与x的对应关系
    x = 2 * np.random.rand(100,1)
    y= 4 + 3*x + np.random.randn(100,1)
    
    line_reg = LinearRegression()
    # 训练数据集,训练完成后,参数会保存在对象line_reg中。
    line_reg.fit(x,y)
    
    # line_reg.intercept为截距,就是w0,line_reg.coef_为其他参数,coef的全拼为coefficient
    print(line_reg.intercept_,line_reg.coef_)
    
    x_new = np.array([[0],[2]])
    # line_reg.predict(x_new) 为预测结果
    print(line_reg.predict(x_new))
    
    plt.plot(x_new,line_reg.predict(x_new),'r-')
    # 画出已知数据X和掺杂了误差的y,用蓝色的点表示
    plt.plot(x,y,'b.')
    # 建立坐标轴
    plt.axis([0,2,0,15,])
    
    plt.show()

      运行结果如下:

      

      如果有的同学想要运行代码,需要安装pycharm和anaconda,将python的interrupt设置为anaconda的bin目录下的python就可以了。网上有很多教程,请原谅这里不再赘述了。

     此处需要说明,因为在使用解析解求解最小二乘的过程中,出现了矩阵求逆的步骤。因为有些矩阵没有逆矩阵,只能使用近似矩阵来代替,所以结果的精度会降低。二则矩阵求逆随着维度的增加,计算量也大大增加,求解速度变慢。所以一般情况下我们都会使用第二种求解办法:梯度下降。

    文章链接:梯度下降和随机梯度下降

          

    展开全文
  • 多元线性回归原理代码实现 原理 多元线性回归是一元线性回归的升级版吧,都是线性模型。 线性回归就是试图学到一个线性模型,尽可能的准确的预测出真实值。就是给机器数据集,其中包括x特征值对应的y值,通过训练...

    多元线性回归

    原理

    多元线性回归是一元线性回归的升级版吧,都是线性模型
    线性回归就是试图学到一个线性模型,尽可能的准确的预测出真实值。就是给机器数据集,其中包括x特征值和对应的y值,通过训练得出一个模型,再只拿一些x特征值给它,这个模型给你预测出较为精准的y值。
    线性回归试图学到的模型是: f ( x ) = ω x i + b f(x)=\omega x_{i}+b f(x)=ωxi+b,使得预测值f(x)跟真实值y结果相似。看着眼熟不?其实本质就有点像我们的 f ( x ) = k x + b f(x)=kx+b f(x)=kx+b这条直线。上面的 ω \omega ω是权重,b是偏置。
    特征值单一的时候用一元线性回归,如果是多个特征决定一个y值,有多个权重,此时用多元线性回归,多元线性回归试图学得模型为: f ( x ) = ω T x i + b f(x)=\omega^{T} x_{i}+b f(x)=ωTxi+b

    在多元线性回归中,使用最小二乘法对参数 ω \omega ω和b进行参数估计,求解 ω \omega ω和b。在此为方便计算,把这两个参数吸收为向量形式进行处理,用 ω ^ \hat{\omega} ω^来表示: ω ^ = ( ω ; b ) \hat{\omega}=(\omega;b) ω^=ω;b,相应的把数据集表示为一个 m ∗ ( d + 1 ) m*(d+1) m(d+1)大小矩阵的X,其中每一行都应一个示例,该行前d个元素对应于示例的d个属性值,把最后一个元素恒置为1便于计算(m个数据,d+1个特征)
    再把y也写为向量形式:y=(y1;y2;y3;…ym)

    由上通过推导公式可得 ω ^ \hat{\omega} ω^的解
    最终学的多元线性回归模型是 f ( x ^ i ) = x ^ i ( X T X ) − 1 X T y f(\hat{x}_{i})=\hat{x}_{i}(X^{T}X)^{-1}X^{T}y f(x^i)=x^i(XTX)1XTy

    具体公式推导有些繁琐,想看看的可私我发手写推导。

    代码实现

    有两种方法进行实现:
    第一种,解析解形式,通过公式代码进行模型训练。
    第二种用sklearn库,简单粗暴的解决,它提供模型的训练函数,不用敲上面推导出来的公式代码了。

    这里使用混凝土抗压强度数据集进行测试,预测值y为混凝土抗压强度定量MPa

    第一种解析解形式:

    import numpy as np
    from sklearn.model_selection import train_test_split
    #混凝土抗压强度数据测试,预测值为混凝土抗压强度定量MPa
    import matplotlib.pyplot as plt
    
    data=np.loadtxt("F:/comdata/Concrete_Data.csv",delimiter=",",skiprows=1,dtype=np.float32);
    #print(data)
    #print(data.shape[0]) 1030行数据
    data1=np.ones((data.shape[0],1)) #全1矩阵,m=1030行n=1列
    #print(data1)
    data=np.hstack((data,data1))
    #print(data)
    y=data[:,0]
    X=data[:,1:]
    #print(y)
    X_train,X_test,y_train,y_test=train_test_split(X,y,test_size=0.8,random_state=42)
    w_hat=np.dot(np.dot(np.linalg.inv(np.dot(X_train.T,X_train)),X_train.T),y_train)
    y_predict=np.dot(X_test,w_hat)
    print(y_predict) #预测的mpa抗压强度值
    
    plt.rcParams['font.sans-serif']=['SimHei']
    plt.title("解析解求解多元线性回归模型得到预测值与真实值对比图")
    plt.xlabel("样本ID")
    plt.ylabel("混凝土抗压强度MPa")
    plt.plot(range(len(X_test)),y_test,alpha=0.4,c="red") #s=2,c="red",
    plt.plot(range(len(X_test)),y_predict,ls="-",alpha=0.5,c="green")
    plt.show()
    
    

    结果:
    在这里插入图片描述

    第二种sklearn实现

    import numpy as np
    from sklearn import linear_model
    from sklearn.model_selection import train_test_split
    import matplotlib.pyplot as plt
    #混凝土抗压强度数据测试,预测值为混凝土抗压强度定量MPa
    #读入数据
    data=np.loadtxt("F:/comdata/Concrete_Data.csv",delimiter=",",skiprows=1,dtype=np.float32)
    #切片,分离数据,X,Y
    X=data[:,1:]
    y=data[:,0]
    #分割数据,分为训练集和测试集
    X_train,X_test,y_train,y_test=train_test_split(X,y,test_size=0.8)
    
    model=linear_model.LinearRegression();#获得线性模型
    #训练模型,用sklearn库里面的
    model.fit(X_train,y_train)
    #查看训练好的参数,就是w
    #print(model.coef_)
    
    #预测模型
    y_predict=model.predict(X_test)
    print(y_predict)
    
    #可视化展示
    plt.rcParams['font.sans-serif']=['SimHei']
    plt.title("sklearn线性回归模型预测值与真实值对比图")
    plt.xlabel("样本ID")
    plt.ylabel("混凝土抗压强度MPa")
    plt.scatter(np.arange(300),y_test[:300],c="red")
    plt.scatter(np.arange(300),y_predict[:300],c="black")
    plt.show()
    

    结果:
    在这里插入图片描述

    数据集若需要可私我取。

    展开全文
  • 线性回归及python代码实现

    千次阅读 2020-02-26 19:52:51
    什么是线性回归 ...回归分析中,只包括一个自变量一个因变量,且二者的关系可用一条直线近似表示,这种回归分析称为一元线性回归分析。如果回归分析中包括两个或两个以上的自变量,且因变量自变量之间是...

    回归分析

    监督学习中,如果预测的变量是离散的,我们称其为分类(如决策树,支持向量机等),如果预测的变量是连续的,我们称其为回归。

    • 离散:连续的对应(就是反义词)就是离散 。离散就是不连续。例如像整数1,2,3,4,5,...这种就是离散的

    在统计学中,回归分析regression analysis)指的是确定两种或两种以上变量间相互依赖的定量关系的一种统计分析方法。

    在大数据分析中,回归分析是一种预测性的建模技术,它研究的是因变量(目标)和自变量(预测器)之间的关系。这种技术通常用于预测分析,时间序列模型以及发现变量之间的因果关系。例如,司机的鲁莽驾驶与道路交通事故数量之间的关系,最好的研究方法就是回归。

    1. 回归分析按照涉及的变量的多少,分为一元回归和多元回归分析;
    2. 按照因变量的多少,可分为简单回归分析和多重回归分析;
    3. 按照自变量和因变量之间的关系类型,可分为线性回归分析和非线性回归分析

    线性回归

    线性回归是利用数理统计中回归分析,来确定两种或两种以上变量间相互依赖的定量关系的一种统计分析方法,运用十分广泛。其表达形式为 y = w ′ x + e y = w'x+e y=wx+e,其中 e e e 为误差服从均值为 0 的正态分布。

    • 回归分析中,只包括一个自变量和一个因变量,且二者的关系可用一条直线近似表示,这种回归分析称为一元线性回归分析。
    • 如果回归分析中包括两个或两个以上的自变量,且因变量和自变量之间是线性关系,则称为多元线性回归分析。

    例子

    假如现在我想去找银行贷款,那我能不能在贷款之前先预测一下我能贷到多少钱呢?下面有这样一张表:

    工资年龄贷款额度
    40002520000
    80003070000
    50002835000
    75003350000
    120004085000

    通过这张表我们可以比较明显的看到工资、年龄与贷款额度之间是有一个关系的,大概是工资越高、年龄越小贷款的额度就相应的越高。那么我能不能根据自己的工资和年龄去预测一个具体的贷款值?

    这里我们就可以引入回归分析了,我们可以把工资年龄当作两个自变量贷款额度当作因变量,通过求解自变量和因变量的线性关系,从而得到一个线性回归方程。利用得到的线性回归方程我们就可以进行预测了,那这种统计分析中利用线性回归方程对自变量和应变量建模的回归分析就叫做线性回归。

    函数拟合

    继续上面的问题,现在我们把工资和年龄设为 x 1 , x 2 x_1,x_2 x1,x2 (特征值),贷款额度设为 h h h (目标值),去求一个线性函数,那么我们可以假设方程为: h = α 0 + α 1 x 1 + α 2 x 2 h=\alpha_0+\alpha_1x_1+\alpha_2x_2 h=α0+α1x1+α2x2,接着假设定义 x 0 x_0 x0的值为 1,上面的式子可以表示成: h = ∑ i = 0 n α i x i = α T X h = \displaystyle\sum_{i=0}^n \alpha_i x_i=\alpha^ T X h=i=0nαixi=αTX

    实际情况中,变量之间不一定有线性关系,这时我们就要通过拟合的方式来解决。也就是找到一组最佳的参数 α 0 , α 1 , α 2 \alpha_0,\alpha_1,\alpha_2 α0,α1,α2 使得预测值最接近真实值。最后我们可能得到这样的结果:
    dk-1
    线性回归模型经常用最小二乘逼近来拟合,当然也可能用别的方法来拟合。

    高斯分布

    通过上面的图,可以看出真实值和预测值之间一般是存在误差(用 ϵ \epsilon ϵ 来表示误差),于是(真实值) y = ∑ i = 0 n α i x i = α T X + ϵ y =\displaystyle\sum_{i=0}^n \alpha_i x_i=\alpha^ T X+\epsilon y=i=0nαixi=αTX+ϵ,对于每个样本有: y ( i ) = α T x ( i ) + ϵ ( i ) y^{(i)} = \alpha^Tx^{(i)}+\epsilon^{(i)} y(i)=αTx(i)+ϵ(i),误差 ϵ ( i ) \epsilon^{(i)} ϵ(i)是独立并且具有相同的分布,并且服从均值为0,方差为 σ 2 \sigma^2 σ2的高斯分布。

    通常我们认为误差往往是很微小的,根据“中心极限定理” :在自然界与生产中,一些现象受到许多相互独立的随机因素的影响,如果每个因素所产生的影响都很微小时,总的影响可以看作是服从正态分布的。中心极限定理就是从数学上证明了这一现象 。所以高斯分布对误差假设来说是一种很好的模型。。

    正态分布(Normal distribution),也称“常态分布”,又名高斯分布(Gaussian distribution);正态曲线呈钟型,两头低,中间高,左右对称因其曲线呈钟形,因此人们又经常称之为钟形曲线。

    若随机变量 X X X服从一个数学期望为 μ \mu μ、方差为 σ 2 \sigma^2 σ2的正态分布,记为 N ( μ , σ 2 ) N(\mu,\sigma^2) N(μ,σ2)。其概率密度函数为正态分布的期望值 μ \mu μ决定了其位置,其标准差 σ \sigma σ决定了分布的幅度。当 μ = 0 , σ = 1 \mu=0,\sigma=1 μ=0,σ=1时的正态分布是标准正态分布。
    通过概率密度函数来定义高斯分布:
    在这里插入图片描述
    高斯分布的概率密度函数: p ( y ) = 1 σ 2 π e − ( y − μ ) 2 2 σ 2 p(y)=\frac{1}{\sigma\sqrt{2\pi}}e^{-\frac{(y-\mu)^2}{2\sigma^2}} p(y)=σ2π 1e2σ2(yμ)2

    由于误差服从高斯分布,将其带入正太分布函数中( μ = 0 \mu=0 μ=0),可以得到: p ( y ( i ) ∣ x ( i ) ; α ) = 1 σ 2 π e − ( y ( i ) − α T x ( i ) ) 2 2 σ 2 p(y^{(i)}|x^{(i)};\alpha)=\frac{1}{\sigma\sqrt{2\pi}}e^{-\frac{(y^{(i)} - \alpha ^Tx^{(i)})^2}{2\sigma^2}} p(y(i)x(i);α)=σ2π 1e2σ2(y(i)αTx(i))2

    似然函数

    在数理统计学中,似然函数是一种关于统计模型中的参数的函数,表示模型参数中的似然性。似然性是用于在已知某些观测所得到的结果时,对有关事物的性质的参数进行估计。给定输出x时,关于参数 θ \theta θ 的似然函数 L ( θ ∣ x ) L(\theta|x) L(θx)(在数值上)等于给定参数 θ \theta θ 后变量X的概率:
    L ( θ ∣ x ) = p ( X = x ∣ θ ) L(\theta|x)=p(X=x|\theta) L(θx)=p(X=xθ)
    似然函数的主要用法在于比较它相对取值,虽然这个数值本身不具备任何含义。例如,考虑一组样本,当其输出固定时,这组样本的某个未知参数往往会倾向于等于某个特定值,而不是随便的其他数,此时,似然函数是最大化的。

    最大似然估计:对同一个似然函数,如果存在一个参数值,使得它的函数值达到最大的话,那么这个值就是最为“合理”的参数值。是似然函数最初也是最自然的应用。在 θ \theta θ 的所有取值上,使这个函数最大化。这个使可能性最大的值即被称为 θ \theta θ 的最大似然估计。

    现在我们要求解的就是最佳参数,也就是满足最大似然估计的 α \alpha α ,所以将 p ( y ( i ) ∣ x ( i ) ; α ) = 1 σ 2 π e − ( y ( i ) − α T x ( i ) ) 2 2 σ 2 p(y^{(i)}|x^{(i)};\alpha)=\frac{1}{\sigma\sqrt{2\pi}}e^{-\frac{(y^{(i)} - \alpha ^Tx^{(i)})^2}{2\sigma^2}} p(y(i)x(i);α)=σ2π 1e2σ2(y(i)αTx(i))2带入:
    L ( α ) = ∏ i = 1 m p ( y ( i ) ∣ x ( i ) ; α ) = ∏ i = 1 m 1 σ 2 π e − ( y ( i ) − α T x ( i ) ) 2 2 σ 2 L(\alpha)=\displaystyle\prod_{i=1}^mp(y^{(i)}|x^{(i)};\alpha)=\displaystyle\prod_{i=1}^m\frac{1}{\sigma\sqrt{2\pi}}e^{-\frac{(y^{(i)} - \alpha ^Tx^{(i)})^2}{2\sigma^2}} L(α)=i=1mp(y(i)x(i);α)=i=1mσ2π 1e2σ2(y(i)αTx(i))2

    对数似然:涉及到似然函数的许多应用中,更方便的是使用似然函数的自然对数形式,即“对数似然函数”。求解一个函数的极大化往往需要求解该函数的关于未知参数的偏导数。由于对数函数是单调递增的,而且对数似然函数在极大化求解时较为方便,所以对数似然函数常用在最大似然估计及相关领域中。
    log ⁡ L ( α ) = log ⁡ ∏ i = 1 m p ( y ( i ) ∣ x ( i ) ; α ) = log ⁡ ∏ i = 1 m 1 σ 2 π e − ( y ( i ) − α T x ( i ) ) 2 2 σ 2 \log L(\alpha)=\log \displaystyle\prod_{i=1}^mp(y^{(i)}|x^{(i)};\alpha)=\log \displaystyle\prod_{i=1}^m\frac{1}{\sigma\sqrt{2\pi}}e^{-\frac{(y^{(i)} - \alpha^T x^{(i)})^2}{2\sigma^2}} logL(α)=logi=1mp(y(i)x(i);α)=logi=1mσ2π 1e2σ2(y(i)αTx(i))2
    = ∑ i = 1 m log ⁡ 1 σ 2 π e − ( y i − α T x i ) 2 2 σ 2 = \displaystyle\sum_{i=1}^m\log\frac{1}{\sigma\sqrt{2\pi}}e^{-\frac{(y^i - \alpha^T x^i)^2}{2\sigma^2}} =i=1mlogσ2π 1e2σ2(yiαTxi)2
    = m log ⁡ 1 σ 2 π − 1 σ 2 ⋅ 1 2 ∑ i = 1 m ( y ( i ) − α T x ( i ) ) 2 \quad \quad =m\log \frac{1}{\sigma\sqrt{2\pi}}-\frac{1}{\sigma^2}\cdot\frac{1}{2}\displaystyle\sum_{i=1}^m(y^{(i)} - \alpha^T x^{(i)})^2 =mlogσ2π 1σ2121i=1m(y(i)αTx(i))2
    目标是要使的对数似然函数最大, m log ⁡ 1 σ 2 π m\log \frac{1}{\sigma\sqrt{2\pi}} mlogσ2π 1是一个常数,所以就得让后面减去的最小,得到目标函数:
    J ( α ) = 1 2 ∑ i = 1 m ( y ( i ) − α T x ( i ) ) 2 = 1 2 ∑ i = 1 m ( Y ( x ( i ) ) − y ( i ) ) 2 = 1 2 ( X α − y ) T ( X α − y ) J(\alpha)=\frac{1}{2}\displaystyle\sum_{i=1}^m(y^{(i)} - \alpha^T x^{(i)})^2=\frac{1}{2}\displaystyle\sum_{i=1}^m(Y(x^{(i)})-y^{(i)})^2=\frac{1}{2}(X\alpha-y)^T(X\alpha-y) J(α)=21i=1m(y(i)αTx(i))2=21i=1m(Y(x(i))y(i))2=21(Xαy)T(Xαy)
    α \alpha α 求偏导:
    ∇ α J ( α ) = ∇ α ( 1 2 ( X α − y ) T ( X α − y ) ) = ∇ α ( 1 2 ( α T X T − y T ) ( X α − y ) ) \nabla_\alpha J(\alpha)=\nabla_\alpha(\frac{1}{2}(X\alpha-y)^T(X\alpha-y))=\nabla_\alpha(\frac{1}{2}(\alpha^TX^T-y^T)(X\alpha-y)) αJ(α)=α(21(Xαy)T(Xαy))=α(21(αTXTyT)(Xαy))
    = ∇ α 1 2 ( α T X T X α − α T X T y − y T X α − y + y T y ) ) =\nabla_\alpha\frac{1}{2}(\alpha^TX^TX\alpha-\alpha^TX^Ty-y^TX\alpha-y+y^Ty)) =α21(αTXTXααTXTyyTXαy+yTy))
    = 1 2 ( 2 X T X α − X T y − ( y T X ) T ) = X T X α − X T y =\frac{1}{2}(2X^TX\alpha-X^Ty-(y^TX)^T)=X^TX\alpha-X^Ty =21(2XTXαXTy(yTX)T)=XTXαXTy
    令偏导等于 0:
    α = ( X T X ) − 1 X T y ( 最 小 二 乘 ) \alpha=(X^TX)^{-1}X^Ty\quad\quad(最小二乘) α=(XTX)1XTy

    最小二乘法

    假如我们在研究两个变量 ( x , y ) (x,y) x,y之间的相互关系时,通常可以得到一系列成对的数据 ( x 1 , y 1. x 2 , y 2... x m , y m ) (x1,y1.x2,y2... xm,ym) x1,y1.x2,y2...xm,ym;将这些数据描绘在 x-y 直角坐标系中,若发现这些点在一条直线附近,可以设这条直线方程为: y i = a 0 + a 1 x ( 其 中 a 0 、 a 1 任 意 实 数 ) y_i=a_0+a_1x(其中 a_0、a_1任意实数) yi=a0+a1x(a0a1)

    为了建立这条直线,那么我们就要确定 a 0 a_0 a0 a 1 a_1 a1的值,我们将实际值 y i y_i yi与计算值 y j ( y j = a 0 + a 1 x ) [ 1 ] y_j(y_j=a_0+a_1x )\quad[1] yj(yj=a0+a1x)[1] 的差 ( y i − y j ) (y_i-y_j) (yiyj)的平方和 φ = ∑ i = 1 n ( y i − y j ) 2 [ 2 ] \varphi=\displaystyle\sum_{i=1}^{n} (y_i-y_j)^2\quad[2] φ=i=1n(yiyj)2[2] 最小作为优化判据
    我们把[1]式带入[2]式中,得到 φ = ∑ i = 1 n ( y i − a 0 − a 1 x i ) 2 [ 3 ] \varphi=\displaystyle\sum_{i=1}^{n} (y_i-a_0-a_1x_i )^2 \quad[3] φ=i=1n(yia0a1xi)2[3]
    为了使 φ \varphi φ 的值最小,我们可以通过求解最小值的方法,让函数 φ \varphi φ a 0 、 a 1 a_0、a_1 a0a1 求偏导,令两个偏导数等于零。这样我们就的到了式[4]和式[5]。
    { ∑ i = 1 n 2 ( a 0 + a 1 x i − y i ) = 0 [ 4 ] ∑ i = 1 n 2 x i ( a 0 + a 1 x i − y i ) = 0 [ 5 ] \begin{cases}\displaystyle\sum_{i=1}^{n} 2(a_0+a_1x_i-y_i )=0 \quad[4]\\\displaystyle\sum_{i=1}^{n} 2x_i(a_0+a_1x_i-y_i )=0 \quad[5]\end{cases} i=1n2(a0+a1xiyi)=0[4]i=1n2xi(a0+a1xiyi)=0[5]
    求解方程组,可以得到 a 0 = 1 n ∑ i = 1 n ( y i − a 1 x i ) , a 1 = n ∑ i = 1 n x i y i − ∑ i = 1 n x i ∑ i = 1 n y i n ∑ i = 1 n x i 2 − ∑ i = 1 n x i ∑ i = 1 n x i a_0=\frac{1}{n}\displaystyle\sum_{i=1}^n(y_i-a_1x_i) , a_1=\frac{n\displaystyle\sum_{i=1}^nx_iy_i-\displaystyle\sum_{i=1}^nx_i\displaystyle\sum_{i=1}^ny_i}{n\displaystyle\sum_{i=1}^nx_i^2-\displaystyle\sum_{i=1}^nx_i\displaystyle\sum_{i=1}^nx_i} a0=n1i=1n(yia1xi),a1=ni=1nxi2i=1nxii=1nxini=1nxiyii=1nxii=1nyi ,将 a 0 a_0 a0 a 1 a_1 a1 带入[1]式中就的到了我们的一元线性方程,即:数学模型。

    这里以最简单的一元线性模型来解释最小二乘法。回归分析中,如果只包括一个自变量和一个因变量,且二者的关系可用一条直线近似表示,这种回归分析称为一元线性回归分析。如果回归分析中包括两个或两个以上的自变量,且因变量和自变量之间是线性关系,则称为多元线性回归分析。对于二维空间线性是一条直线;对于三维空间线性是一个平面,对于多维空间线性是一个超平面。

    梯度下降

    梯度下降是迭代法的一种,沿着这个函数下降的方向找,逼近最小值,可以用于求解最小二乘问题。在求解机器学习算法的模型参数,即无约束优化问题时,梯度下降(Gradient Descent)是最常采用的方法之一,另一种常用的方法是最小二乘法。

    在求解损失函数的最小值时,可以通过梯度下降法来一步步的迭代求解,得到最小化的损失函数和模型参数值。

    其迭代公式为 a k + 1 = a k + ρ k s ˉ ( k ) a_{k+1}=a_k+\rho k^{\bar{s}^{(k)}} ak+1=ak+ρksˉ(k) ,其中 s ˉ ( k ) \bar{s}^{(k)} sˉ(k) 代表梯度负方向, ρ k \rho k ρk 表示梯度方向上的搜索步长。梯度方向我们可以通过对函数求导得到,步长的确定比较麻烦,太大了的话可能会发散,太小收敛速度又太慢。

    评估

    通常使用RSE(残差标准差)和 R 来评估模型

    R 2 = 1 − R S S T S S = 1 − ∑ i = 1 m ( y i ^ − y i ) 2 ∑ i = 1 m ( y i − y ˉ ) 2 残 差 平 方 和 类 似 方 差 项 R^2=1-\frac{RSS}{TSS}=1-\frac{\displaystyle\sum_{i=1}^m(\hat{y_i}-y_i)^2}{\displaystyle\sum_{i=1}^m(y_i-\bar{y})^2}\quad\quad\quad\quad \frac{残差平方和}{类似方差项} R2=1TSSRSS=1i=1m(yiyˉ)2i=1m(yi^yi)2
    R 2 R^2 R2 可以衡量因变量的变化比例,并用自变量 x x x 表述。因此,假设在一个线性方程中,自变量 x x x 可以解释因变量,那么变化比例就高, R 2 R^2 R2 将接近 1。反之,则接近0。所以 R 2 R^2 R2的取值越接近 1 ,我们就认为模型拟合的与好

    简单线性回归python案例

    导入用到的包

    import numpy as np
    import matplotlib.pyplot as plt
    

    生成数据

    # 设置随机生成算法的初始值
    np.random.seed(300)
    data_size = 100
    #  生出size个符合均分布的浮点数,取值范围为[low, high),默认取值范围为[0, 1.0)
    x = np.random.uniform(low=1.0, high=10.0, size=data_size)
    y = x * 20 + 10 + np.random.normal(loc=0.0, scale=10.0, size=data_size)
    

    查看数据

    print(x)
    print(y.shape)
    plt.scatter(x, y)
    plt.show()
    

    [5.06010312 2.98920108 4.32173127 3.61698133 3.16970135 8.0384179 8.29984929 4.60711487 1.23399282 6.48435743 6.66179289 6.86728789 8.72731994 6.28122472 7.33751513 5.37311508 7.15974709 2.01131411 2.64152172 7.27343316 4.31517483 4.84567523 4.59836046 6.83840906 4.7201361 2.47241009 5.57532693 5.65184176 4.13384671 3.72027221 3.04603181 1.49130937 7.68238581 3.78143155 9.79455216 4.27478265 5.42937585 7.76009388 5.53962905 2.12052548 8.35322843 4.32938243 3.18827513 8.28250053 3.09192118 3.2946656 4.92652317 8.48630166 6.78218775 7.45073433 7.36576515 7.96927182 4.13018028 7.03576017 6.43589797 7.3106855 9.21539573 9.0112781 9.5537087 7.9801062 5.44528565 6.62427759 4.46709413 5.58802991 8.67279532 3.2084603 9.61530432 8.82526244 2.18656259 5.71768665 7.5517145 3.89302569 8.74359262 2.08295081 5.91832925 5.46471639 9.43168746 1.58599839 7.63747752 4.87551191 8.56518562 7.88432931 3.48748529 1.72525194 6.83057399 2.76268027 6.48435105 8.16047282 2.30107294 2.63178986 8.74522455 5.82819158 6.58947199 5.36626558 4.85639013 6.27787997 3.56740106 7.88396527 1.66851535 4.97195607]
    (100,)
    在这里插入图片描述

    区分训练集和测试集

    # 返回一个随机排列
    shuffled_index = np.random.random(data_size)
    x = x[shuffled_index]
    y = y[shuffled_index]
    split_index = int(data_size * 0.75)
    x_train = x[:split_index]
    y_train = y[:split_index]
    x_test = x[split_index:]
    y_test = y[split_index:]
    

    线性回归类

    class LinerRegression():
    
        def __init__(self, learning_rate=0.01, max_iter=100, seed=None):
            np.random.seed(seed)
            self.lr = learning_rate
            self.max_iter = max_iter
            self.a = np.random.normal(1, 0.1)
            self.b = np.random.normal(1, 0.1)
            self.loss_arr = []
    
        def fit(self, x, y):
            self.x = x
            self.y = y
            for i in range(self.max_iter):
                self._train_step()
                self.loss_arr.append(self.loss())
    
        def _f(self, x, a, b):
            """一元线性函数"""
            return x * a + b
    
        def predict(self, x=None):
            """预测"""
            if x is None:
                x = self.x
            y_pred = self._f(x, self.a, self.b)
            return y_pred
    
        def loss(self, y_true=None, y_pred=None):
            """损失"""
            if y_true is None or y_pred is None:
                y_true = self.y
                y_pred = self.predict(self.x)
            return np.mean((y_true - y_pred)**2)
    
        def _calc_gradient(self):
            """梯度"""
            d_a = np.mean((self.x * self.a + self.b - self.y) * self.x)
            d_b = np.mean(self.x * self.a + self.b - self.y)
            print(d_a, d_b)
            return d_a, d_b
    
        def _train_step(self):
            """训练频度"""
            d_a, d_b = self._calc_gradient()
            self.a = self.a - self.lr * d_a
            self.b = self.b - self.lr * d_b
            return self.a, self.b
    

    训练数据

    regr = LinerRegression(learning_rate=0.01, max_iter=10, seed=314)
    regr.fit(x_train, y_train)
    

    展示

    print(f'cost: \t{regr.loss():.3}')
    print(f"f={regr.a:.2} x + {regr.b:.2}")
    plt.scatter(np.arange(len(regr.loss_arr)), regr.loss_arr, marker='o', c='green')
    plt.show()
    

    cost: 1.11e+02
    f=2.1e+01 x+4.2
    在这里插入图片描述
    在这里插入图片描述

    评估

    y_pred = regr.predict(x_test)
    def sqrt_R(y_pred, y_test):
        y_c = y_pred - y_test
        rss = np.mean(y_c ** 2)
        y_t = y_test - np.mean(y_test)
        tss = np.mean(y_t ** 2)
        return (1 - rss / tss)
    r = sqrt_R(y_pred, y_test)
    print(f"R: {r}")
    

    R: 0.9327511073300032


    sklearn案例

    导包

    from sklearn.datasets import load_boston
    from sklearn.externals import joblib
    from sklearn.linear_model import LinearRegression, SGDRegressor, Ridge
    from sklearn.metrics import mean_squared_error
    from sklearn.model_selection import train_test_split
    from sklearn.preprocessing import StandardScaler
    

    载入数据

    """线性回归预测房子价格"""
    # 获取数据
    lb = load_boston()
    print(lb)
    

    处理数据

    # 分割数据集 (训练集、测试集)
    x_train, x_text, y_train, y_test = train_test_split(lb.data, lb.target, test_size=0.25)
    # 进行标准化
    # 特征值 目标值 都标准化
    std_x = StandardScaler()
    x_train = std_x.fit_transform(x_train)
    x_text = std_x.fit_transform(x_text)
    
    std_y = StandardScaler()
    y_train = std_y.fit_transform(y_train.reshape(-1, 1))  # 传入数据必须二维
    y_test = std_y.fit_transform(y_test.reshape(-1, 1))
    

    方程求解

    # estimator预测
    # 正规方程求解方式
    lr = LinearRegression()
    lr.fit(x_train, y_train)
    
    # print(lr.coef_)
    
    y_lr_predict = std_y.inverse_transform(lr.predict(x_text))
    # print("预测:\n", y_lr_predict)
    print("均方误差:", mean_squared_error(std_y.inverse_transform(y_test), y_lr_predict))
    print(f"R : {lr.score(x_test, y_test)}")
    

    [[-0.10633041 0.10911523 0.00409775 0.06205962 -0.23187842 0.29700813 -0.00083156 -0.32863065 0.2570794 -0.18807137 -0.22742163 0.10321854 -0.39041034]
    均方误差: 21.92038221080247
    R : 0.6956028823910707

    梯度下降

    # 梯度下降方式
    sgd = SGDRegressor()
    sgd.fit(x_train, y_train)
    
    y_sgd_predict = std_y.inverse_transform(sgd.predict(x_text))
    print(y_sgd_predict)
    print("均方误差:", mean_squared_error(std_y.inverse_transform(y_test), y_sgd_predict))
    print(f"R : {sgd.score(x_test, y_test)}")
    

    均方误差: 17285.95206027733
    R : 0.6907551016044395

    展开全文
  • IDL 多元线性回归求解

    2010-04-01 08:46:36
    求解多元线性回归方程解的IDL代码,其中包括矩阵的高斯消元的步骤,txt文件保存的是导入的矩阵数据,矩阵的行列数根据数据确定
  • ref. 《机器学习》周志华 P53一、线性模型线性模型(linear model)试图学得一个通过...二、(多元)线性回归1. 问题的提出给定数据集D如下:则(多元)线性回归试图学得,使得,其中={w1;w2;……;wd;b}。为...

    ref. 《机器学习》周志华 P53

    一、线性模型

    线性模型(linear model)试图学得一个通过属性的线性组合来进行预测的函数,即


    线性模型形式简单,易于建模。w直观表达了各个属性在预测任务中的重要性,因此线性模型有很好的可解释性(comprehensibility)。

    二、(多元)线性回归

    1. 问题的提出

    给定数据集D如下:

    则(多元)线性回归试图学得,使得,其中={w1;w2;……;wd;b}。

    为进行向量化计算,将数据表示为X,每行对应一个样本,如下所示。在每个样本最后增加元素“1”,是为了用于拟合f()里的偏置b。


    则现在要拟合的函数为f(X)=

    2. 问题求解

    ① 目标

    确定的关键在于如何衡量f(x)与y之间的差别,通常使用均方误差作为性能度量,使其最小化作为目标,即,


    这种基于均方误差最小化来进行模型求解的方法称为“最小二乘法”。

    ② XTX是满秩矩阵的情况:

     同样为了向量化计算,需要把y也表示为向量的形式,。则问题求解的目标如(3.9)所示。


    当XTX为满秩矩阵或正定矩阵时,式(3.10)推导出式(3.11)略。

    式(3.11)是我们求解参数w,b的关键公式。

    ③ XTX不是满秩矩阵的情况:

    在现实任务中的往往遇到X的属性个数大于样本个数的情况,即XTX不是满秩矩阵,此时用上述方法会得到多个,且它们都能使均方误差最小化,选择哪一个解作为输出,将由学习算法的归纳偏好决定。常见的做法是引入正则化项,可参考岭回归、LASSO、Elastic Net等(我上篇博文有简单介绍)。

    三、我的笔记

    我最近的工作中,领导给出的任务是拟合一个包含7个样本,每个样本7个属性的数据集,其实就是一个7个变量的7个方程组,要求方程组没有b。接到这个任务,立刻就想到了是用线性回归问题。在机器学习算法中,对待拟合的函数普遍都加了偏置b,在我这个要求不设置b的任务中,将“二、(多元)线性回归”的求解问题稍改变一下即可,如下:

    X无需在最后一列增加元素“1”(因为无需拟合b),仍然使用式(3.11)进行的求解,这个中不再含有b。

    四、python代码

    1. “二、(多元)线性回归”部分涉及的代码

    # 线性回归 ref.《机器学习》.周志华.P55
    def fitLM(X,Y):
        samplesize = np.size(X, 0)
        # 1. 样本最后一个元素置1
        X_b = np.ones([1, samplesize]).tolist()
        X = np.transpose(X).tolist()
        X.extend(X_b)
        X = np.transpose(X)
    
        # 2. 根据公式(3.11)求w*
        XT=np.transpose(X)
        equ1=np.linalg.inv(np.dot(XT,X))
        equ2=np.dot(XT,Y)
        w=np.dot(equ1,equ2)
        return w
    
    

    2.  “三、我的笔记”部分涉及的代码

    def fitLM(X,Y):
        samplesize = np.size(X, 0)
        
        # 2. 根据公式(3.11)求w*
        XT=np.transpose(X)
        equ1=np.linalg.inv(np.dot(XT,X))
        equ2=np.dot(XT,Y)
        w=np.dot(equ1,equ2)
        return w

    展开全文
  • 简单线性回归(最小二乘法) # 0. 引入依赖 import numpy as np import matplotlib.pyplot as plt # 1. 导入数据 points = np.genfromtxt('D:\学习资料\推荐系统\代码\练习代码\自己手打\data.csv', delimiter = ',...
  • 线性回归模型代码实现

    千次阅读 2020-02-13 19:30:05
    线性回归 线性回归基本要素 1.模型 根据房屋面积(平方米)以及房龄(年)预测房屋价格: 2.数据集 我们通常收集一系列的真实数据,例如多栋房屋的真实售出价格它们对应的面积房龄。对模型训练以及测试,一栋...
  • 梯度下降法求解一元线性回归 依然是这个房价预测的任务,这是一个一元线性回归问题,这次我们采用梯度下降法来求解它可以分为5步 第1步加载样本数据x,y 第2步设置超参数,在这个例子中,超参数包括学习率迭代...
  • 目录1、线性回归的原理基础定义公式推导简单理解2、最小二乘法PYTHON实现0. 导入相关库1. 导入数据2. 定义损失函数3. 定义算法拟合函数4. 测试定义的函数5. 画出拟合曲线3、最小二乘简单例子 1、线性回归的原理 基础...
  • 本篇内容本来想在写在上篇博客中的,那样篇幅过长,就单独提出来了。 本篇文章采用两种方式实现线性回归,一种是使用scikit...而通过上篇博客,我们已经知道了最小二乘法求解线性回归参数,所以完全可以自己手动实现
  • 文章目录原理以及公式【1】一元线性回归问题【2】多元线性回归问题【3】学习率【4】流程分析(一元线性回归)【5】流程分析(多元线性回归)归一化原理以及每种归一化适用的场合一元线性回归代码以及可视化结果多元...
  • matlab 多元线性回归模型求解 实现对多影响因素的分析
  • 梯度下降法求解多元线性回归 我们来用编程实现一个多元线性回归问题。 这是样本数据 其中的属性包括房屋面积房间数。可以看到面积房间数的取值范围相差很大,如果直接使用这样的数据来训练模型,这个面积的贡献...
  • 线性模型基本形式:f(x)=w1*x1+w2*x2+w3*x3+...+wd*xd+b 向量形式:f(x)=w'x+b(w'指w转置w'=(w1,w2,w3,...,wd)) 回归任务最常用均方误差作为性能度量,见下图 ...如对数线性回归:lny=w'x+b,让e^(w'x+b)逼近y ...
  • 一元线性回归的Python代码测试数据,主要是包含一个属性值label的测试数据
  • 深入理解线性回归的数学推导过程 能够使用原生代码实现线性回归模型(包括Ridge,LASSO,EasticNet) 能够使用skLearn实现线性回归模型 能够使用线性回归分析实际数据
  • 最小二乘法线性回归作为最基础的线性回归,在统计机器学习中都有重要的地位。在机器学习中,线性回归用来从数据中获得启示来帮助预测,因此如何得到最拟合数据的函数防止过拟合是研究重点。 假设我们的拟合函数...
  • matlab实现一元线性回归和多元线性回归

    万次阅读 多人点赞 2018-01-30 10:58:46
    回归分析中,如果有两个或两个以上的自变量,就称为多元回归。事实上,一种现象常常是与多个因素相联系的,由多个自变量的最优组合共同来预测或估计因变量,比只用一个自变量进行预测或估计更有效,更符合实际。 ...
  • 前面利用梯度下降算法求解线性回归参数,我们设置的最大迭代次数为100,此时计算出的theta参数值为:matrix([[0.18629876, 3.48668984]]) 利用正规方程 # 定义正规方程函数: def normalEqn(X, y): theta ...
  • 代码实现简单线性回归(单变量)

    万次阅读 多人点赞 2019-03-02 18:23:12
    1. 使用自己生成数据(代码实现) import numpy as np m = 1000000 test_x = np.random.random(size=m) test_y = test_x * 2 + ...2. 算术法求解代码实现) class LinearRegressionArithmetical: #算术法实现 d...
  • 多元线性回归 及其Python实现

    万次阅读 多人点赞 2019-04-05 21:14:11
    多元线性回归求解过程 多元线性回归的形式: 目标函数: 将一个样本的向量化: 将所有样本的向量化: 向量化后的目标函数及求解结果: ps.上述多元线性回归的正规方程解问题是:事件复杂度高;优点是:不需要对...
  • 线性回归其实就是寻找一条直线拟合数据点,使得损失函数最小。直线的表达式为: yi=ω1xi,1+ω2xi,2+ωjxi,j+...+by_i = \omega_1x_{i,1}+\omega_2x_{i,2}+\omega_jx_{i,j}+...+b 损失函数的表达式为: J=12∑i=0...
  • 线性回归原理和实现基本认识

    万次阅读 多人点赞 2017-04-28 17:04:57
     线性回归在假设特证满足线性关系,根据给定的训练数据训练一个模型,并用此模型进行预测。先举个简单的例子;我们假设一个线性方程 Y=2x+1, x变量为商品的大小,y代表为销售量;当月份x =5时,我们就能根据线性...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 16,033
精华内容 6,413
关键字:

线性回归的求解和代码实现