精华内容
下载资源
问答
  • 2019-11-26 13:01:26

    【更新日志】
    3/3/2020 对部分公式中出现的错误进行了修正
    4/5/2020 修改了文章标题


    1. 前文回顾

    在上一篇文章中,我们建立了多元线性回归模型,以及模型的相关假设,并给出了对应的样本模型(详情请见:【统计学习系列】多元线性回归模型(一)——模型的建立与基本假设)(别问我为什么点开404,因为这一章的内容太多我还在补充中,暂未开放):
    y i = β 0 + ∑ j = 1 p x i j β j + ϵ i ,   i = 1 , . . . , N y_i= \beta_0 + \sum_{j=1}^{p} x_{ij} \beta_j + \epsilon_i , \ i=1,...,N yi=β0+j=1pxijβj+ϵi, i=1,...,N其中:
    ϵ = ( ϵ i ) N × 1 ∼ N ( 0 , I n σ 2 ) \bm{\epsilon} = (\epsilon_i)_{N \times 1} \thicksim N(0,\bm{I_n}\sigma^2) ϵ=(ϵi)N×1N(0,Inσ2)
    因此,需要估计的参数有 p + 2 个:β0, β1, …, βp, σ。我们如何利用样本来估计模型参数呢?我们的一般思路是:1)制定评价标准,2)在给定标准下寻找最优参数。这篇文章主要介绍最小二乘估计法极大似然估计法这两种方法。话不多说,让我们开始探寻参数估计的秘密吧!

    :本篇涉及大量数学定理推导。尽管作者力求内容通俗易懂,但同时也希望保证证明的严谨性。因此在参数估计的推导中,我会把证明思路与过程尽可能清晰、完整地展示出来,这可能需要读者具有一定程度的凸优化问题求解和线性代数的基础。实操应用类读者可直接跳过2、3两章,而直接阅读第4章结论。


    2. 最小二乘法估计(Ordinary Least Squared Estimate, OLS)

    如果将 N 组样本对 (yi, xi) 看做是一个 p+1 维实空间中的N个点,那么我们现在要做的就是在空间中找到一个 p 维超平面,来尽可能“好”的拟合空间中的这N个样本的点。什么样的指标可以衡量这种拟合的好坏呢?样本点到拟合平面的距离则是一个衡量拟合好坏的测量工具:当点到平面距离和越小,说明估计量与真实值之间的“距离”越小(离得越近),模型对样本数据的拟合情况越好;而点到平面距离和越大,说明估计量与真实值之间的“距离”越大(离得越远),模型对样本数据的拟合情况越差。而对于距离的度量,我们可以使用两点差值的平方这一指标:
    d i s t a n c e ( y i , y ^ i ) = ( y i − y ^ i ) 2 distance(y_i, \hat{y}_i) = (y_i - \hat{y}_i)^2 distance(yi,y^i)=(yiy^i)2
    注1:使用平方而非绝对值是为了后边方便求导。
    注2:这里的表述并不严谨。作为距离的测度,我们应该使用范数(例如2-范数)。其本质是因为定义范数作为距离的欧式空间是一个赋范线性空间。

    最小二乘估计法的核心思想是:找到一组参数 β ,使得“样本点到平面的距离和最小”,或者说最小化残差平方和。用数学语言可以表示为:

    min ⁡ β 0 , β 1 , . . . , β p R S S ( β 0 , β 1 , . . . , β p ) = ∑ i = 1 N ( y i − f ( x i ) ) 2 = ∑ i = 1 N ( y i − β 0 − ∑ j = 1 p x i j β j ) 2 \min_{\beta_0, \beta_1,...,\beta_p} RSS(\beta_0, \beta_1,...,\beta_p) \\ \hskip{1.5em} = \sum_{i=1}^{N}(y_i - f(\bm{x_i} ))^2 \\ \hskip{5em} = \sum_{i=1}^{N}(y_i - \beta_0 - \sum_{j=1}^{p} x_{ij} \beta_j)^2 β0,β1,...,βpminRSS(β0,β1,...,βp)=i=1N(yif(xi))2=i=1N(yiβ0j=1pxijβj)2

    注1:利用这一准则成立的合理前提是:(xi, yi)是从总体一种独立抽取的随机样本。若样本不满足随机性,yi与xi也应该条件独立。

    注2:在利用这一准则进行参数估计时,其过程本身并不蕴含任何假设条件(即该过程并不蕴含模型的有效性假设)

    为方便计算,我们首先将RSS的表达式改写为矩阵运算形式:
    R S S ( β ) = ( y − X β ) T ( y − X β ) RSS( \bm{\beta} ) = (\bm{y} - \bm{X} \bm{\beta} )^T (\bm{y} - \bm{X} \bm{\beta} ) RSS(β)=(yXβ)T(yXβ)

    其中:
    β = [ β 0 β 1 ⋮ β p ] ( p + 1 ) × 1 ,       y = [ y 0 y 1 ⋮ y N ] N × 1 ,       X = [   1    x 11    …    x 1 p   1    x 21    …    x 2 p ⋮       ⋮       ⋱      ⋮   1    x N 1    …    x N p ] N × ( p + 1 ) \bm{\beta} = \begin{bmatrix} \beta_0 \\ \beta_1 \\ \vdots\\ \beta_p \end{bmatrix}_{(p+1) \times1}, \space \space \space \space \space \bm{y} = \begin{bmatrix} y_0 \\ y_1 \\ \vdots\\ y_N \end{bmatrix}_{N \times1}, \space \space \space \space \space \bm{X} = \begin{bmatrix} \space 1 \space\space x_{11} \space\space \dots \space\space x_{1p} \\ \space 1 \space\space x_{21} \space\space \dots \space\space x_{2p} \\ \vdots \space \space\space \space \space \vdots \space \space \space \space \space ⋱ \space \space \space \space \vdots\\ \space 1 \space\space x_{N1} \space\space \dots \space\space x_{Np} \end{bmatrix}_{N \times (p+1)} β=β0β1βp(p+1)×1,     y=y0y1yNN×1,     X= 1  x11    x1p 1  x21    x2p               1  xN1    xNpN×(p+1)
    由于问题(I)可以视为无约束优化问题,分别对目标函数求其一阶和二阶导数,则有:
    ∂ R S S ( β ) ∂ β = − 2 X T ( y − X β ) \frac{ \partial RSS( \bm{\beta} ) } { \partial \bm{\beta} } = -2 \bm{X}^T ( \bm{y} - \bm{X} \bm{\beta}) βRSS(β)=2XT(yXβ) ∂ 2 R S S ( β ) ∂ β 2 = 2 X T X ≽ 0 \frac{ \partial^2 RSS( \bm{\beta} ) } { \partial \bm{\beta}^2 } = 2 \bm{X}^T \bm{X} \succcurlyeq 0 β22RSS(β)=2XTX0

    由线性代数和凸优化问题(或多元函数极值问题)基础知识可知:RSS(β)为向量 β半正定二次型,同时RSS(β)是定义在 Rp+1 上的凸函数

    XTX 满秩(亦称非奇异)时,RSS(β)为向量 β正定二次型,因此RSS(β)是定义在 Rp+1 上的严格凸函数。由严格凸函数性质可知,RSS(β)最小值点存在且唯一,且其取得最小值满足的充分必要条件是:RSS(β)关于向量 β一阶导函数为0

    根据上述结论,该最优化问题的解满足: X T ( y − X β ^ ) = 0 \bm{X}^T ( \bm{y} - \bm{X} \bm{ \hat\beta}) = 0 XT(yXβ^)=0

    移项得: X T X β ^ = X T y \bm{X}^T \bm{X} \bm{\hat\beta} = \bm{X}^T \bm{y} XTXβ^=XTy

    当自变量不存在多重共线性时, XTX 矩阵非奇异,因此其逆矩阵存在,因此解得优化问题的解:
    β ^ = ( X T X ) − 1 X T y \bm{\hat\beta} = ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y} β^=(XTX)1XTy
    注1矩阵 XTX 非奇异, 等价于 样本资料阵 X 列满秩,等价于模型自变量不存在多重共线性。这部分证明请参考文献[1]《高等代数(第二版)上册 》
    注2:关于矩阵函数的代数运算,例如乘法运算、逆运算、求导运算等更多细节,请参考文献[1]《高等代数(第二版)上册 》
    注3:关于凸规划问题的相关基础知识,请参考文献[2]《数值最优化》

    将估计参数 β 带入模型中,可以得到学习样本的拟合值:
    y ^ = X β ^ = X ( X T X ) − 1 X T y = d e f H y \bm{\hat y} = \bm{X} \bm{\hat\beta} = \bm{X} ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y} \xlongequal{def} \bm{H} \bm{y} y^=Xβ^=X(XTX)1XTydef Hy其中,矩阵 H H = X ( X T X ) − 1 X T \bm{H} = \bm{X} ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T H=X(XTX)1XT被称为“帽子矩阵”(仿佛给 y 带了一个帽子)。

    参数 σ 表示误差项的波动率。我们可以把残差序列作为误差项的估计量,因此可以使用残差序列的样本标准差对 σ 进行估计:
    σ ^ = ( y − y ^ ) T ( y − y ^ ) N − p − 1 = ∑ i = 1 N ( y i − y ^ i ) 2 N − p − 1 \bm{\hat \sigma} = \sqrt{ \frac {( \bm{y} - \bm{\hat{y}} )^T ( \bm{y} - \bm{\hat{y}} ) } {N-p-1} } = \sqrt{ \frac { \sum _{i=1}^{N} (y_i - \hat{y}_i)^2} {N-p-1} } σ^=Np1(yy^)T(yy^) =Np1i=1N(yiy^i)2

    :这里分母使用 N-p-1 是为了确保该方差估计量的无偏性,其具体论证我将在下一篇文章中进行说明。


    3. 极大似然估计(Maximum Likelihood Estimate, ML)

    在上一章使用最小二乘法对模型参数进行估计时,我们采用了用“残差平方和”作为评判模型拟合的好坏,其实质是希望找到一个拟合超平面(即为一个线性模型),使得所有样本到超平面的距离和最短(也就是最接近)。因此,这种评估方法是站在几何学的角度上进行的。除此之外,我们还有没有其他衡量拟合好坏的角度呢?

    首先我们知道,样本具有“二元性” :当在抽样之前,样本可以视为一个 随机变量(或随机向量);而抽样之后,样本的取值被固定,因此又可以视为一个常数(或常向量)。因此,在抽样之前,N 个来自于总体的独立随机样本应该有联合分布,这是随机变量的特征决定的。若多元线性回归模型的正态性假设成立,那么N 个样本的因变量 y 的条件分布应满足:
    y ∼ N ( X β , I n σ 2 ) \bm{y} \thicksim N(\bm{X} \bm{\beta},\bm{I_n}\sigma^2) yN(Xβ,Inσ2)

    联合概率密度函数存在,有:

    f ( y 1 , . . . , y N ; x 1 , . . . , x N , β , σ ) = ∏ i = 1 N 1 2 π σ exp ⁡ { − 1 2 σ 2 ( y i − β 0 − β 1 x i 1 − ⋯ − β p x i p ) 2 } = 1 ( 2 π σ ) N × exp ⁡ { − ∑ i = 1 N 1 2 σ 2 ( y i − β 0 − β 1 x i 1 − ⋯ − β p x i p ) 2 } = ( 2 π σ ) − N × exp ⁡ { − 1 2 σ 2 ( y − X β ) T ( y − X β ) } f ( y_1, ..., y_N; \bm{x_1} , ..., \bm{x_N}, \bm{\beta }, \sigma ) \\ \hspace{4.5em} = \prod_{i=1}^{N} \frac {1} {\sqrt{2\pi}\sigma} \exp\{ - \frac{1}{2\sigma^2} (y_i - \beta_0 - \beta_1 x_{i1}- \dots - \beta_p x_{ip})^2 \} \\ \hspace{4.5em} = \frac {1} { ( \sqrt{2\pi}\sigma ) ^ N} \times \exp\{ - \sum_{i=1}^N \frac{1}{2\sigma^2}(y_i - \beta_0 - \beta_1 x_{i1}- \dots - \beta_p x_{ip} ) ^2 \} \\ \hspace{4.5em} =( \sqrt{2\pi}\sigma ) ^ {-N} \times \exp \{ - \frac{1}{2\sigma^2} (\bm{y} - \bm{X} \bm{\beta} )^ T(\bm{y} - \bm{X} \bm{\beta} ) \} f(y1,...,yN;x1,...,xN,β,σ)=i=1N2π σ1exp{2σ21(yiβ0β1xi1βpxip)2}=(2π σ)N1×exp{i=1N2σ21(yiβ0β1xi1βpxip)2}=(2π σ)N×exp{2σ21(yXβ)T(yXβ)}其中: φ(x) 为一元标准正态分布的分布函数。

    :关于随机变量与分布更多的资料请参考文献[3]《测度论与概率论基础》

    在抽样之后,样本被固定,联合概率密度函数变成了关于总体参数 β 的函数。我们重新定义这个函数,称这一函数为“似然函数”(Likelihood Function),并记为 L (β),其表达式为:
    L ( β , σ ;   ( y 1 , x 1 ) , . . . , ( y N , x N ) ) = ( 2 π σ ) − N × exp ⁡ { − 1 2 σ 2 ( y − X β ) T ( y − X β ) } L( \bm{\beta } , \sigma; \space (y_1, \bm{x_1}), ..., (y_N, \bm{x_N}) ) =( \sqrt{2\pi}\sigma ) ^ {-N} \times \exp \{ - \frac{1}{2\sigma^2} (\bm{y} - \bm{X} \bm{\beta} )^ T(\bm{y} - \bm{X} \bm{\beta} ) \} L(β,σ; (y1,x1),...,(yN,xN))=(2π σ)N×exp{2σ21(yXβ)T(yXβ)}

    :似然函数与联合概率密度函数在表达形式上一致,但是似然函数是关于总体参数 β 的函数,而概率密度函数是关于随机变量序列(即样本序列)( yi, xi ), i = 1, …, N 的函数。

    为什么称他为似然函数呢?我们知道,概率密度函数展现了的随机变量发生的可能性的大小。因此似然函数的意义是:若我们已经抽出来了一组样本,那么这组样本来自参数为 β 的总体的可能性是多少

    :这个说法不够严谨,但这是我能想到的最容易理解的表达了。。。

    那么我们如何利用似然函数来估计参数呢?我们期望:“在一次试验中,若某一事件发生了,那么这一事件最有可能来自于发生可能性最大的哪一种情形”,或者说“发生概率最高的事件在一次实验中最有可能发生”。

    这听起来完全就是废话啊!我们不妨举一个简单的例子:

    例:一个箱子里有十个球,这十个球有可能是下面三种情况:
    情况 1:十个球中,一个红球,九个白球;
    情况 2:十个球中,五个红球,五个白球;
    情况 3:十个球中,九个红球,一个白球;
    已知有放回地随机抽,每次抽一个,结果三次都抽中了红球。请问,箱子里最有可能是这三种情况中的哪一种?

    分别计算这三种情况所发生“有放回三次都抽中红球”这一事件的概率,可以知道情况3发生的可能性最大。因此我们可以认为箱子里装有九个红球和一个白球。

    :对一思想进行严格说明需要涉及贝叶斯决策,有兴趣的读者可自行查阅更多资料

    通过这个例子的说明可知:我们需要找到一组参数,使得在这组参数下,样本的联合概率密度达到最大。因此,这个问题就变成了:在所有有可能的参数取值中(专业说法:在参数空间中),所抽中的样本来自哪一组参数的可能性最高

    至此,我们将问题抽象成了一个无约束优化问题,这与最小二乘估计的手段十分相似。但是由于似然函数是连乘运算,这在求导的过程中不够方便,因此我们对似然函数其进行对数化处理,得到“对数似然函数”(Logarithm Likelihood Function ),并记为 l (β)。因此,这一优化问题转化为了:

    min ⁡ β l ( β , σ ;   ( y 1 , x 1 ) , ( y 2 , x 2 ) , . . . , ( y N , x N ) ) = − N ln ⁡ 2 π σ − 1 2 σ 2 ( y − X β ) T ( y − X β ) \min_ \bm{\beta} l( \bm{\beta }, \sigma;\space (y_1, \bm{x_1}), ( y_2, \bm{x_2}), ... , (y_N, \bm{x_N}) ) = -N \ln { \sqrt{2\pi}\sigma } - \frac{1}{2\sigma^2} ( \bm{y} - \bm{X} \bm{\beta} )^ T(\bm{y} - \bm{X} \bm{\beta} ) βminl(β,σ; (y1,x1),(y2,x2),...,(yN,xN))=Nln2π σ2σ21(yXβ)T(yXβ)

    :不难证明:对数似然函数与似然函数取最小值时,β 相等。

    与第2部分过程相仿,其一阶导函数有:
    ∂ l ( β , σ ) ∂ β = − 2 X T ( y − X β ) = 0 \frac{ \partial l( \bm{\beta},\sigma ) } { \partial \bm{\beta} } = -2 \bm{X}^T ( \bm{y} - \bm{X} \bm{\beta}) = 0 βl(β,σ)=2XT(yXβ)=0 ∂ l ( β , σ ) ∂ σ = − N σ + 1 σ 3 ( y − X β ) T ( y − X β ) = 0 \frac{ \partial l( \bm{\beta},\sigma) } { \partial \sigma } = -\frac{N}{\sigma} + \frac{1}{\sigma^3} ( \bm{y} - \bm{X} \bm{\beta} )^ T(\bm{y} - \bm{X} \bm{\beta} ) = 0 σl(β,σ)=σN+σ31(yXβ)T(yXβ)=0

    解得:
    β ^ = ( X T X ) − 1 X T y \bm{\hat\beta} = ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y} β^=(XTX)1XTy σ ^ = 1 N ( y − X β ) T ( y − X β ) \hat{\sigma} = \sqrt { \frac{1}{N} ( \bm{y} - \bm{X} \bm{\beta} )^ T(\bm{y} - \bm{X} \bm{\beta} ) } σ^=N1(yXβ)T(yXβ)


    4. 结论

    基于第2部分和第3部分的讨论,我们知道:基于最小二乘法极大似然估计法两种方法得到的参数估计结果是一致的有没有一种殊途同归的艺术感!),估计值为:
    β ^ = ( X T X ) − 1 X T y \bm{\hat\beta} = ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y} β^=(XTX)1XTy

    其中:
    β = [ β 0 β 1 ⋮ β p ] ( p + 1 ) × 1 ,       y = [ y 0 y 1 ⋮ y N ] N × 1 ,       X = [   1    x 11    …    x 1 p   1    x 21    …    x 2 p ⋮       ⋮       ⋱      ⋮   1    x N 1    …    x N p ] N × ( p + 1 ) \bm{\beta} = \begin{bmatrix} \beta_0 \\ \beta_1 \\ \vdots\\ \beta_p \end{bmatrix}_{(p+1) \times1}, \space \space \space \space \space \bm{y} = \begin{bmatrix} y_0 \\ y_1 \\ \vdots\\ y_N \end{bmatrix}_{N \times1}, \space \space \space \space \space \bm{X} = \begin{bmatrix} \space 1 \space\space x_{11} \space\space \dots \space\space x_{1p} \\ \space 1 \space\space x_{21} \space\space \dots \space\space x_{2p} \\ \vdots \space \space\space \space \space \vdots \space \space \space \space \space ⋱ \space \space \space \space \vdots\\ \space 1 \space\space x_{N1} \space\space \dots \space\space x_{Np} \end{bmatrix}_{N \times (p+1)} β=β0β1βp(p+1)×1,     y=y0y1yNN×1,     X= 1  x11    x1p 1  x21    x2p               1  xN1    xNpN×(p+1)
    不过这两种方法对误差项方差的估计不同。最小二乘法对误差项标准差的估计量为:
    σ ^ O L S = 1 N − p − 1 ( y − y ^ ) T ( y − y ^ ) \bm{\hat \sigma} _{OLS} = \sqrt{ \frac {1 } {N-p-1} ( \bm{y} - \bm{\hat{y}} )^T ( \bm{y} - \bm{\hat{y}} ) } σ^OLS=Np11(yy^)T(yy^)

    极大似然估计对误差项标准差的估计量为: σ ^ M L = 1 N ( y − y ^ ) T ( y − y ^ ) \bm{\hat \sigma} _{ML} = \sqrt{ \frac { 1 } {N} ( \bm{y} - \bm{\hat{y}} )^T ( \bm{y} - \bm{\hat{y}} ) } σ^ML=N1(yy^)T(yy^) 其中: y ^ = X β ^ = X ( X T X ) − 1 X T y \bm{\hat y} = \bm{X} \bm{\hat\beta} = \bm{X} ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y} y^=Xβ^=X(XTX)1XTy
    :在机器学习领域中,我们首先设置目标函数(或称代价方程),通过寻找最优的参数使得目标函数达到最小值,从而实现参数估计的目的。因此从这个意义上讲,多元线性回归模型可以算机器学习模型中的一种。


    5. 新问题的提出

    在估计出了参数之后,小伙伴们 (无中生有) 一定会顺理成章地提出如下疑问:

    (1)这些估计出来的参数好不好?准不准确呢?准确率有多高呢?
    (2)估计出来的模型好不好?有没有效呢?
    (3)估计出来的模型是否满足模型的基本假设条件呢?不满足可咋办呢?
    (4)如何利用估计出来的模型进行预测呢?预测的准确性有多高呢?
    (5)你还得拖多长时间才往后更新啊?

    面对小伙伴们的灵魂拷问,让我们下回分解。


    参考文献

    [1]《高等代数(第二版)上册 》丘维声 著
    [2]《数值最优化》李董辉等 著
    [3]《测度论与概率论基础》程士宏 著


    写在最后

    欢迎感兴趣的小伙伴来跟作者一起挑刺儿~ 包括但不限于语言上的、排版上的和内容上的不足和疏漏~ 一起进步呀!
    有任何问题,欢迎在本文下方留言,或者将问题发送至勘误邮箱: mikeysun_bugfix@163.com
    谢谢大家!XD

    更多相关内容
  • 多元线性回归参数估计方法,吴仕勋,赵东方,本文依据高斯—马尔可夫定理,通过对最小二乘估计方法得出的参数估计值的分析,从另外两个角度出发得出了参数估计的值与最小二乘
  • 基于杂交粒子群算法的多元线性回归参数估计及预测区间研究.pdf
  • 5.4多元线性回归中的参数估计.pdf5.4多元线性回归中的参数估计.pdf5.4多元线性回归中的参数估计.pdf5.4多元线性回归中的参数估计.pdf5.4多元线性回归中的参数估计.pdf5.4多元线性回归中的参数估计.pdf5.4多元线性...
  • 5.4多元线性回归中的参数估计.docx5.4多元线性回归中的参数估计.docx5.4多元线性回归中的参数估计.docx5.4多元线性回归中的参数估计.docx5.4多元线性回归中的参数估计.docx5.4多元线性回归中的参数估计.docx5.4多元...
  • 其9个相关影响变量数据(包括风力,机动车保有量,火电厂、炼钢厂、炼焦厂平均各排口每小时各主要污染物的排放量),在MATALB中采用多元线性回归方法建立了模型、参数估计和模型检验,并在已得模型的基础上剔除不...
  • 一元线性回归:梯度下降法 一元线性回归线性回归的最简单的一种,即只有一个特征变量。首先是梯度下降法,这是比较经典的求法。一元线性回归通俗易懂地说,就是一元一次方程。只不过这里的斜率和截距要通过最小...

    一元线性回归:梯度下降法

    一元线性回归是线性回归的最简单的一种,即只有一个特征变量。首先是梯度下降法,这是比较经典的求法。一元线性回归通俗易懂地说,就是一元一次方程。只不过这里的斜率和截距要通过最小二乘和梯度下降法不断迭代找到最优解。我们先来看看运用到的代价函数:最小二乘法。
    在这里插入图片描述
    这其实和高中学的最小二乘法一模一样,不过值得注意的是,这里的2其实是可以消去的,这对结果的影响不大,之所以保留是因为,方便与之后求导所得的2消掉。
    梯度下降法:
    开始就说过,斜率和截距是要通过迭代计算求得的。因为初始化的斜率和截距所求得的代价函数(最小二乘法)的值很大,所以要通过梯度下降法来不断改变斜率和截距的值,从而达到全局最优解局部最优解。全局最优解和局部最优解,这两者是有分别的,如下图:
    在这里插入图片描述
    θ0和θ1分别对应斜率和截距,J(θ0,θ1)就是代价函数的值。假设红色的点代价函数的值较大,蓝色的点较小,θ0和θ1的初始值的不同会导致所位于的红色峰值不同。如果J(θ0,θ1)的值在左边的红色山峰上(黑色十角星),那么经过迭代得到的就是左边的蓝色峰底,也就是全局最优解(因为J(θ0,θ1)最小)。但是,如果J(θ0,θ1)的值在右边的红色山峰上,那么那么经过迭代得到的就是右边的蓝色峰底,即局部最优解(J(θ0,θ1)不是最小)。
    重点来了,怎么进行θ0和θ1的更新呢?用到下面的公式:
    在这里插入图片描述
    新的值=新的值—学习率*代价函数对旧的值的求导。这就是梯度下降法的更新公式, 学习率是步长的意思,这里不过多解释。这里要注意学习率的取值:
    在这里插入图片描述
    梯度下降法适用于多种场景,对于线性回归的解法是这样的——
    首先是求导的结果:
    在这里插入图片描述
    在这里插入图片描述
    然后带入公式得到:
    在这里插入图片描述

    线性回归的原理大概就是这样了,接下来用python进行实战训练——
    这里用到的数据是这样的,很多行的由逗号分隔的两列数据
    在这里插入图片描述
    接下来是代码部分,首先是载入数据,
    在这里插入图片描述
    让我们看看数据的分布图像是怎么样的?
    在这里插入图片描述
    定义参数
    在这里插入图片描述
    程序主体部分:
    在这里插入图片描述在这里插入图片描述
    开始测试
    在这里插入图片描述
    测试结果
    在这里插入图片描述
    在这里插入图片描述
    从测试结果我们可以看到代价函数的值从最开始的5565经过迭代后变为112,而且最后显示的图像也比较符合,这就是梯度下降法在一元线性回归的应用。

    一元线性回归:sklearn

    其实在Python中,很多算法别人已经帮写好了,只需要直接调用即可,方便快捷,结果也是一样的。

    开始都是一样的,注意改变数据的维度
    在这里插入图片描述
    两句话就可以完成模型的创建和拟合
    在这里插入图片描述
    最后输出一样的结果
    在这里插入图片描述

    多元线性回归:梯度下降法

    现实生活中,很多事情不仅仅是一个特征,所以一元线性回归不怎么常用,多元线性回归比较常见。多元线性回归和一元线性回归的本质都是一样的,不同的是,多元线性回归的特征较多,但是核心算法不变,多元线性回归显示的图像是3d立体图,如下面的例子:
    数据是两个特征(左边两个),然后得到一个结果(最右边)
    在这里插入图片描述
    python实现:
    这里导入了画3d图的包
    在这里插入图片描述
    因为有两个特征,所以这里的最小二乘法要有两个参数
    在这里插入图片描述
    梯度更新的代码也是多加了一个参数而已:
    在这里插入图片描述
    运行程序:
    在这里插入图片描述
    在这里插入图片描述
    结果显示:
    在这里插入图片描述
    在这里插入图片描述
    从结果我们可以看到,代价函数的值从开始的47点多降到了1左右,而且图像也比较吻合。这就是梯度下降法在多元线性回归中的应用。

    多元线性回归:sklearn

    和一元线性回归一样,多元线性回归也可以使用sklearn进行计算。
    python的实现,开始都是一样的
    在这里插入图片描述
    开始创建并拟合模型,这里是有两个系数
    在这里插入图片描述
    运行:
    在这里插入图片描述
    这里的结果和梯度下降法不太一样,因为sklearn的算法是标准方程法
    在这里插入图片描述

    以上就是一元线性回归和多元线性回归在机器学习里的一般解决方法。

    展开全文
  • 多元线性回归参数估计,介绍多元线性回归参数估计
  • 衡量参数估计的指标2.1 无偏性2.2 一致性2.3 有效性3. 一些引理3.1 期望运算的线性性3.2 期望运算的线性性4. *β*^~OLS~ 的性质4.1 *β*^~OLS~ 服从的分布4.2 *β*^~OLS~ 与误差项之间的关系4.3 *β*^~OLS~ 的无偏...

    【更新日志】
    4/5/2020 对文章中公式与排版的部分错误进行修正


    1. 前文回顾

    在前面的文章中,我们介绍了多元线性回归模型的两种参数估计。对于模型:
    y i = β 0 + ∑ j = 1 p x i j β j + ϵ i ,   i = 1 , . . . , N y_i= \beta_0 + \sum_{j=1}^{p} x_{ij} \beta_j + \epsilon_i , \ i=1,...,N yi=β0+j=1pxijβj+ϵi, i=1,...,N ϵ = ( ϵ i ) N × 1 ∼ N ( 0 , I n σ 2 ) \bm{\epsilon} = (\epsilon_i)_{N \times 1} \thicksim N(0,\bm{I_n}\sigma^2) ϵ=(ϵi)N×1N(0,Inσ2)
    利用最小二乘估计法(OLS)得到的参数估计量为:
    β ^ O L S = ( X T X ) − 1 X T y \bm{\hat\beta}_{OLS} = ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y} β^OLS=(XTX)1XTy σ ^ O L S = 1 N − p − 1 ( y − X β ^ O L S ) T ( y − X β ^ O L S ) \bm{\hat \sigma} _{OLS} = \sqrt{ \frac {1 } {N-p-1} ( \bm{y} - \bm{X} \bm{\hat{\beta}}_{OLS} )^T ( \bm{y} -\bm{X} \bm{\hat{\beta}}_{OLS} ) } σ^OLS=Np11(yXβ^OLS)T(yXβ^OLS)
    而利用极大似然估计法(ML)得到的参数估计量为:
    β ^ M L = ( X T X ) − 1 X T y \bm{\hat\beta}_{ML} = ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y} β^ML=(XTX)1XTy σ ^ M L = 1 N ( y − X β ^ M L ) T ( y − X β ^ M L ) \bm{\hat \sigma} _{ML} = \sqrt{ \frac { 1 } {N} ( \bm{y} - \bm{X} \bm{\hat{\beta}}_{ML} )^T ( \bm{y} - \bm{X} \bm{\hat{\beta}}_{ML} ) } σ^ML=N1(yXβ^ML)T(yXβ^ML) 其中:
    β = [ β 0 β 1 ⋮ β p ] ( p + 1 ) × 1 ,       y = [ y 0 y 1 ⋮ y N ] N × 1 ,       X = [   1    x 11    …    x 1 p   1    x 21    …    x 2 p ⋮       ⋮       ⋱      ⋮   1    x N 1    …    x N p ] N × ( p + 1 ) \bm{\beta} = \begin{bmatrix} \beta_0 \\ \beta_1 \\ \vdots\\ \beta_p \end{bmatrix}_{(p+1) \times1}, \space \space \space \space \space \bm{y} = \begin{bmatrix} y_0 \\ y_1 \\ \vdots\\ y_N \end{bmatrix}_{N \times1}, \space \space \space \space \space \bm{X} = \begin{bmatrix} \space 1 \space\space x_{11} \space\space \dots \space\space x_{1p} \\ \space 1 \space\space x_{21} \space\space \dots \space\space x_{2p} \\ \vdots \space \space\space \space \space \vdots \space \space \space \space \space ⋱ \space \space \space \space \vdots\\ \space 1 \space\space x_{N1} \space\space \dots \space\space x_{Np} \end{bmatrix}_{N \times (p+1)} β=β0β1βp(p+1)×1,     y=y0y1yNN×1,     X= 1  x11    x1p 1  x21    x2p               1  xN1    xNpN×(p+1)
    在这两种方法中,最小二乘估计法是我们最为常用的参数估计方法,原因有二:其一,以残差平方和作为评判超平面拟合好坏的指标比较直观,也易于理解;其二,最小二乘估计参数具有良好的估计量性质。此外,以残差平方和构建代价函数(Cost Function)的方法在机器学习领域十分常见,我们会在后面的Lasso回归和Ridge回归之中继续体会其中的奥妙。

    在本章中,我们来继续探索利用最小二乘估计量在性质上究竟具有什么样的优势吧~


    2. 衡量参数估计量好坏的指标

    由于参数估计量是总体参数的估计值,因此估计量必然与总体参数之间存在一定的“误差”。如果,我们如何衡量一个参数的估计量是好是坏呢?我们可以从以下这几个性质入手:

    这里想编一个射击小游戏来说明,但是有点费脑筋,等编好了再补充上来吧。

    2.1 无偏性

    参数估计量的 无偏性(unbiasedness) 是指,在多次试验中,用总体的某参数估计值的平均值与该总体参数的真实值“没有偏差”。用数学语言来描述则可以表达为:若一个总体参数 β 的估计量 β^ 是无偏估计量,则该估计量应满足:
    E [ β ^ ] = β E[\bm{\hat\beta}] = \bm\beta E[β^]=β

    2.2 一致性

    参数估计量的 一致性(Consistency) 是指,当样本数量足够大的时候,总体参数的估计值以某种意义收敛到该参数真值,即:
    β ^ → n → ∞ 在 某 种 意 义 下 β \bm{\hat\beta} \xrightarrow[n\rarr\infin]{在某种意义下} \bm\beta β^ nβ
    (1)若 β^ 依概率收敛β ,则称 β^β弱一致估计
    (2)若 β^ 以概率1收敛β ,则称 β^β强一致估计

    注1无偏性是一种小样本性质,而一致性则是一种大样本性质
    注2一致估计又称为相合估计相容估计

    2.3 有效性

    参数估计量的 有效性(Validness) 是指,在参数的任意一无偏估计量中,该无偏估计量的方差最小,即:


    E [ β ^ ] = β E[\bm{\hat\beta}]= \bm\beta E[β^]=β ∀ β ~ ∈ { β ~ : E [ β ~ ] = β } var [ β ^ ] ≤ var [ β ~ ] \forall \bm{\tilde\beta} \in \{ \bm{\tilde\beta}: E[\bm{\tilde\beta}]= \bm\beta \} \\ \text{var}[\bm{\hat\beta}] \le \text{var}[\bm{\tilde\beta}] β~{β~:E[β~]=β}var[β^]var[β~]

    则称 β^β有效估计量(Valid Estimator)


    3. 一些引理(可略)

    为了保证后续证明的严谨性,本文列出一些比较重要的引理。在实际应用过程中可以忽略这些引理的证明过程而直接使用其结论。

    3.1 期望运算的线性性

    【引理1 期望线性性】 对于任一 n 阶随机向量 ym × n 线性变换矩阵 Am 阶随机向量 Ay 有:
    E [ A y ] = A E [ y ] E[\bm{A}\bm{y}]=\bm{A}E[\bm{y}] E[Ay]=AE[y]
    Proof:
    由于随机变量(向量)的期望由Riemann-Stieltjes积分定义:
    E [ y ] = ∫ y   d F ( y ) E[\bm{y}]=\int\bm{y}\ d\bm{F}( \bm{y} ) E[y]=y dF(y)

    其中,F(·)为随机变量(向量)的分布集函数。

    由Riemann-Stieltjes积分的线性性,可以证明:
    E [ A y ] = ∫ A y   d F ( y ) = A ∫ y   d F ( y ) = A E [ y ] E[\bm{Ay}]=\int\bm{Ay}\ d\bm{F}( \bm{y} )= \bm{A}\int\bm{y}\ d\bm{F}( \bm{y} ) = \bm{A}E[\bm{y}] E[Ay]=Ay dF(y)=Ay dF(y)=AE[y]
    Q.E.D.

    注:有关Riemann-Stieltjes积分的定义与运算性质可以参考相关文献

    3.2 协方差运算的半线性性

    【引理2 协方差半线性性】 对于任一 n 阶随机向量 xy,与两 m × n 线性变换矩阵 AB ,有:
    cov ( A x , B y ) = A cov ( x , y ) B T \text{cov} (\bm{Ax},\bm{By})= \bm{A}\text{cov} (\bm{x},\bm{y})\bm{B}^T cov(Ax,By)=Acov(x,y)BT
    Proof:
    cov ( A x , B y ) = E [ ( A x − E [ A x ] ) ( B y − E [ B y ] ) T ] = E [ ( A x − A E [ x ] ) ( B y − B E [ y ] ) T ] = A E [ ( x − E [ x ] ) ( y − E [ y ] ) T ] B T = A cov ( x , y ) B T \text{cov} (\bm{Ax},\bm{By}) \\ =E[(\bm{Ax}-E[\bm{Ax}])(\bm{By}-E[\bm{By}])^T] \\ =E[(\bm{Ax}-\bm{A}E[\bm{x}])(\bm{By}-\bm{B}E[\bm{y}])^T] \\ =\bm{A}E[(\bm{x}-E[\bm{x}])(\bm{y}-E[\bm{y}])^T]\bm{B}^T \\ =\bm{A}\text{cov} (\bm{x},\bm{y})\bm{B}^T cov(Ax,By)=E[(AxE[Ax])(ByE[By])T]=E[(AxAE[x])(ByBE[y])T]=AE[(xE[x])(yE[y])T]BT=Acov(x,y)BT
    Q.E.D.

    3.3 矩阵迹运算的性质

    【引理3 矩阵迹运算的性质】 对于任意 m × nn × m 阶实矩阵 AB ,其迹运算(trace)满足:
    t r ( A B ) = t r ( B A ) tr(AB)=tr(BA) tr(AB)=tr(BA)
    Proof:
    根据矩阵乘法与迹运算的定义:
    t r ( A B ) = ∑ i = 1 m ∑ j = 1 n a i , j b j , i tr(AB) = \sum_{i=1}^{m} \sum_{j=1}^{n} a_{i,j} b_{j,i} tr(AB)=i=1mj=1nai,jbj,i t r ( B A ) = ∑ j = 1 n ∑ i = 1 m b j , i a i , j tr(BA) = \sum_{j=1}^{n} \sum_{i=1}^{m} b_{j,i} a_{i,j} tr(BA)=j=1ni=1mbj,iai,j
    由乘法交换律和加法分配律易得:
    t r ( A B ) = ∑ i = 1 m ∑ j = 1 n a i , j b j , i = ∑ j = 1 n ∑ i = 1 m b j , i a i , j = t r ( B A ) tr(AB) = \sum_{i=1}^{m} \sum_{j=1}^{n} a_{i,j} b_{j,i} = \sum_{j=1}^{n} \sum_{i=1}^{m} b_{j,i} a_{i,j} = tr(BA) tr(AB)=i=1mj=1nai,jbj,i=j=1ni=1mbj,iai,j=tr(BA)
    Q.E.D.


    4. β^OLS 的性质

    4.1 β^OLS 服从的分布

    若模型的正态性假设成立,即:

    ϵ ∼ N ( 0 , I n σ 2 ) \bm{\epsilon} \thicksim N(0,\bm{I_n}\sigma^2) ϵN(0,Inσ2)

    则有:
    y = X β + ϵ ∼ N ( X β , I n σ 2 ) \bm{y} = \bm{X} \bm\beta + \bm\epsilon \thicksim N( \bm{X} \bm\beta, \bm{I_n}\sigma^2) y=Xβ+ϵN(Xβ,Inσ2)
    因此,在给定自变量的条件下,y 服从于均值为 ,协方差矩阵为 Inσ2 的条件正态分布。

    又因为 β^OLS 满足:
    β ^ O L S = ( X T X ) − 1 X T y \bm{\hat\beta}_{OLS} = ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y} β^OLS=(XTX)1XTy
    β^OLS 关于因变量 y线性变换(Linear Transformation)。由正态分布的性质可知, β^OLS 服从正态分布,且其均值有:
    E [ β ^ O L S ] = E [ ( X T X ) − 1 X T y ] = ( X T X ) − 1 X T E [ y ] = ( X T X ) − 1 X T X β = β E[\bm{\hat\beta}_{OLS} ] = E[ ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y} ] \\ \hspace{4em} = ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T E[\bm{y} ] \\ \hspace{3.75em} = ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{X} \bm\beta \\ \hspace{-3.5em} = \bm\beta E[β^OLS]=E[(XTX)1XTy]=(XTX)1XTE[y]=(XTX)1XTXβ=β
    其方差有:
    var [ β ^ O L S ] = var [ ( X T X ) − 1 X T y ] = cov [ ( X T X ) − 1 X T y , ( X T X ) − 1 X T y ] = ( X T X ) − 1 X T cov [ y , y ] X ( X T X ) − 1 = σ 2 ( X T X ) − 1 X T I n X ( X T X ) − 1 = σ 2 ( X T X ) − 1 \text{var} [\bm{\hat\beta}_{OLS}] \\ \hspace{3.25em} = \text{var} [( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y}] \\ \hspace{7em} = \text{cov} [( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y}, ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y}] \\ \hspace{7.25em} = ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \text{cov} [\bm{y}, \bm{y}] \bm{X} ( \bm{X}^T \bm{X} )^{-1} \\ \hspace{7.25em} = \sigma^2 ( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{I}_n \bm{X} ( \bm{X}^T \bm{X} )^{-1} \\ \hspace{-0.25em} = \sigma^2 ( \bm{X}^T \bm{X} )^{-1} var[β^OLS]=var[(XTX)1XTy]=cov[(XTX)1XTy,(XTX)1XTy]=(XTX)1XTcov[y,y]X(XTX)1=σ2(XTX)1XTInX(XTX)1=σ2(XTX)1
    至此,可以说明:在方差正态性满足的前提下,β^OLS 满足:
    β ^ O L S ∼ N ( β , σ 2 ( X T X ) − 1 ) \bm{\hat\beta}_{OLS} \thicksim N(\bm\beta, \sigma^2 ( \bm{X}^T \bm{X} )^{-1} ) β^OLSN(β,σ2(XTX)1)

    4.2 β^OLS 与误差项之间的关系

    由4.1中的论述,我们知道 β^OLS 是因变量 y 的线性表示,而 y 又是误差项 ϵ 的线性表示。实际上:
    cov ( y , ϵ ) = cov ( X β + ϵ , ϵ ) = I n σ 2 \text{cov} (\bm{y}, \bm\epsilon) = \text{cov} (\bm{X}\bm\beta + \bm\epsilon, \bm\epsilon) =\bm{I}_n\sigma^2 cov(y,ϵ)=cov(Xβ+ϵ,ϵ)=Inσ2
    因此,β^OLS 与误差项 ϵ 存在相关关系,其协防矩阵有:
    cov ( β ^ , ϵ ) = cov [ ( X T X ) − 1 X T y , ϵ ] = ( X T X ) − 1 X T cov ( y , ϵ ) = σ 2 ( X T X ) − 1 X T \text{cov} (\bm{\hat\beta},\bm\epsilon) = \text{cov} [( \bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y}, \epsilon] \\ = (\bm{X}^T \bm{X} )^{-1} \bm{X}^T \text{cov} (\bm{y}, \bm\epsilon) \\ = \sigma ^ 2( \bm{X}^T \bm{X} )^{-1} \bm{X}^T cov(β^,ϵ)=cov[(XTX)1XTy,ϵ]=(XTX)1XTcov(y,ϵ)=σ2(XTX)1XT

    4.3 β^OLS 的无偏性

    首先,我们来证明,β^OLSβ的无偏估计量。
    实际上,由4.1中的证明,我们已经得到:
    E [ β ^ O L S ] = β E[\bm{\hat\beta}_{OLS} ] = \bm\beta E[β^OLS]=β
    因此,最小二乘估计量 β^OLSβ无偏估计量

    :由于 β 的最小二乘法估计量与极大似然估计量在表达形式上相等,因此极大似然估计量 β^ML 也是 β 的无偏估计量

    4.4 β^OLS 的一致性

    下面,让我们继续证明, β^OLSβ弱一致估计量,在某些特定条件下,β^OLSβ强一致估计量。由于估计量的一致性证明过程较为艰深冗长,详细的证明过程与结论可以参考文献[1]:线性回归估计相合性问题的新进展.

    4.5 β^OLS 的有效性

    βOLS 的有效性可以由 高斯-马尔科夫定理(Gauss-Markov Theorm) 证得。

    【Gauss-Markov定理】β 的所有线性无偏估计量中, β^OLS方差最小的线性无偏估计量

    Proof:

    假设 cTyβ 的一个线性无偏估计量,其中, cT = (XTX)-1XT + DD 为一常矩阵。则有:
    E [ c T y ] = [ ( X T X ) − 1 X T + D ] X β = ( I + D X ) β = β E[\bm{c}^T\bm{y}] = [(\bm{X^TX})^{-1} \bm{X}^T + \bm{D}] \bm{X \beta} \\ =(\bm{I+DX} )\bm\beta =\bm\beta E[cTy]=[(XTX)1XT+D]Xβ=(I+DX)β=β

    因此可知:
    D X = 0 \bm{DX} =\bm{0} DX=0

    而:
    var ( c T y ) = c T var ( y ) c = σ 2 c T c = σ 2 [ ( X T X ) − 1 X T + D ] [ ( X T X ) − 1 X T + D ] T = σ 2 [ ( X T X ) − 1 + D D T ] ≥ σ 2 ( X T X ) − 1 = var ( β ^ ) \text{var}(\bm{c}^T\bm{y}) =\bm{c}^T \text{var}(\bm{y}) \bm{c} =\sigma^2 \bm{c}^T \bm{c} \\ \\ \hspace{5em} = \sigma^2 [(\bm{X^TX})^{-1} \bm{X}^T + \bm{D}] [(\bm{X^TX})^{-1} \bm{X}^T + \bm{D}]^T \\ =\sigma^2 [(\bm{X^TX})^{-1} +\bm{D}\bm{D}^T] \\ \ge \sigma^2 (\bm{X^TX})^{-1} = \text{var}(\hat{\bm{\beta}}) var(cTy)=cTvar(y)c=σ2cTc=σ2[(XTX)1XT+D][(XTX)1XT+D]T=σ2[(XTX)1+DDT]σ2(XTX)1=var(β^)
    Q.E.D.


    5. σ^OLS 的性质

    5.1 σ^2OLS 的无偏性

    我们先将σ^2OLS的表达式进行变型:
    ( N − p − 1 ) σ ^ O L S 2 = ( y − X β ^ O L S ) T ( y − X β ^ O L S ) = [ y − X ( X T X ) − 1 X T y ] T [ y − X ( X T X ) − 1 X T y ] = y T [ I N − X ( X T X ) − 1 X T ] T [ I N − X ( X T X ) − 1 X T ] y = y T [ I N − X ( X T X ) − 1 X T ] y = ( X β + ϵ ) T [ I N − X ( X T X ) − 1 X T ] ( X β + ϵ ) = β T X T [ I N − X ( X T X ) − 1 X T ] X β ( 1 ) + ϵ T [ I N − X ( X T X ) − 1 X T ] X β ( 2 ) + β T X T [ I N − X ( X T X ) − 1 X T ] ϵ ( 3 ) + ϵ T [ I N − X ( X T X ) − 1 X T ] ϵ ( 4 ) (N-p-1)\hat \sigma _{OLS}^2 \\ = ( \bm{y} - \bm{X} \bm{\hat{\beta}}_{OLS} )^T ( \bm{y} -\bm{X} \bm{\hat{\beta}}_{OLS} ) \\ = [\bm{y} - \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y}]^T [\bm{y} - \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T \bm{y}] \\ = \bm{y}^T [\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ]^T [\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ] \bm{y} \\ =\bm{y}^T [\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ] \bm{y} \\ =(\bm{\bm{X}\bm{\beta} + \bm{\epsilon}} )^T [\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ] (\bm{\bm{X}\bm{\beta} + \bm{\epsilon}} ) \\ = \bm{\beta}^T\bm{X}^T[\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ] \bm{X}\bm{\beta} \hspace{3em}(1) \\ + \bm{\epsilon}^T[\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ] \bm{X}\bm{\beta} \hspace{3em}(2) \\ + \bm{\beta}^T\bm{X}^T[\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ] \bm{\epsilon}\hspace{3em}(3) \\ +\bm{\epsilon}^T[\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ] \bm{\epsilon}\hspace{3em}(4) (Np1)σ^OLS2=(yXβ^OLS)T(yXβ^OLS)=[yX(XTX)1XTy]T[yX(XTX)1XTy]=yT[INX(XTX)1XT]T[INX(XTX)1XT]y=yT[INX(XTX)1XT]y=(Xβ+ϵ)T[INX(XTX)1XT](Xβ+ϵ)=βTXT[INX(XTX)1XT]Xβ(1)+ϵT[INX(XTX)1XT]Xβ(2)+βTXT[INX(XTX)1XT]ϵ(3)+ϵT[INX(XTX)1XT]ϵ(4)
    至此,我们将的表达式分成了四个部分。容易计算:第(1)项、第(2)项和第(3)项均 恒为0 。因此:

    σ ^ O L S 2 = 1 ( N − p − 1 ) ϵ T ( I N − X ( X T X ) − 1 X T ) ϵ \hat \sigma _{OLS}^2 = \frac{1}{(N-p-1)} \bm{\epsilon}^T (\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ) \bm{\epsilon} σ^OLS2=(Np1)1ϵT(INX(XTX)1XT)ϵ

    根据引理3关于迹运算的性质,以及迹运算与期望运算的相容性:

    E [ ϵ T ( I N − X ( X T X ) − 1 X T ) ϵ ] = E [ t r { ( I N − X ( X T X ) − 1 X T ) ϵ ϵ T } ] = t r { ( I N − X ( X T X ) − 1 X T ) E [ ϵ ϵ T ] } = t r { ( I N − X ( X T X ) − 1 X T ) I N σ 2 } = σ 2 [ t r { I N } − t r { X ( X T X ) − 1 X T } ] = σ 2 [ t r { I N } − t r { X T X ( X T X ) − 1 } ] = σ 2 [ t r { I N } − t r { I ( p + 1 ) } ] = ( N − p − 1 ) σ 2 E[\bm{\epsilon}^T (\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ) \bm{\epsilon}] \\ = E[tr\{(\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ) \bm{\epsilon} \bm{\epsilon}^T \}] \\ = tr\{(\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ) E[ \bm{\epsilon} \bm{\epsilon}^T ] \} \\ = tr\{(\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ) \bm{I}_N \sigma^2 \} \\ =\sigma^2 [ tr\{\bm{I}_N \} - tr\{\bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T \} ] \\ =\sigma^2 [ tr\{\bm{I}_N \} - tr\{\bm{X}^T\bm{X} (\bm{X}^T \bm{X} )^{-1} \} ] \\ =\sigma^2 [ tr\{\bm{I}_N \} - tr\{\bm{I}_{(p+1)} \} ] \\ = (N-p-1)\sigma^2 E[ϵT(INX(XTX)1XT)ϵ]=E[tr{(INX(XTX)1XT)ϵϵT}]=tr{(INX(XTX)1XT)E[ϵϵT]}=tr{(INX(XTX)1XT)INσ2}=σ2[tr{IN}tr{X(XTX)1XT}]=σ2[tr{IN}tr{XTX(XTX)1}]=σ2[tr{IN}tr{I(p+1)}]=(Np1)σ2
    因此:
    E [ σ ^ O L S 2 ] = σ 2 E[\hat \sigma _{OLS}^2] =\sigma^2 E[σ^OLS2]=σ2

    至此,我们证得:σ^2OLSσ2无偏估计量

    注1:显然,σ^2ML 不是 σ2 的无偏估计 ;
    注2:需要特别注意的是,σ^OLS 不是 σ 的无偏估计量(这里的平方号不能少)。


    5.2 σ^2OLS 所服从的分布

    σ^2OLS 的表达式:
    σ ^ O L S 2 = 1 ( N − p − 1 ) ϵ T ( I N − X ( X T X ) − 1 X T ) ϵ \hat \sigma _{OLS}^2 = \frac{1}{(N-p-1)} \bm{\epsilon}^T (\bm{I}_N- \bm{X} (\bm{X}^T \bm{X} )^{-1} \bm{X}^T ) \bm{\epsilon} σ^OLS2=(Np1)1ϵT(INX(XTX)1XT)ϵ
    可以看出,σ^2OLS 是关于误差项 ϵ 的二次型。因此,容易证明: σ^2OLS/ σ2 服从自由度为 N - p - 1的卡方分布,即:
    ( N − p − 1 ) σ ^ O L S 2 σ 2 ∼ χ N − p − 1 2 \frac {(N-p-1)\hat \sigma _{OLS}^2} {\sigma^2} \thicksim \chi^2_{N-p-1} σ2(Np1)σ^OLS2χNp12


    6. 结论

    至此,通过一系列冗长但富有启发性的证明,我们从 无偏性一致性有效性 这三个角度,对多元线性回归模型的最小二乘估计量 β^OLSσ^2OLS 的质量进行了评判,得到的结论如下:

    (1)β^OLSβ最优线性无偏估计量 (Best Linear Unbiased Estimator,BLUE),即其具有 无偏性、一致性与有效性

    (2)σ^2OLSσ2无偏估计量(Unbiased Estimator)。

    同时,在误差项满足正态性假设的条件下,我们推导出了 β^OLSσ^2OLS 所服从的分布:
    β ^ O L S ∼ N ( β , σ 2 ( X T X ) − 1 ) \bm{\hat\beta}_{OLS} \thicksim N(\bm\beta, \sigma^2 ( \bm{X}^T \bm{X} )^{-1} ) β^OLSN(β,σ2(XTX)1) ( N − p − 1 ) σ ^ O L S 2 σ 2 ∼ χ 2 ( N − p − 1 ) \frac {(N-p-1)\hat \sigma _{OLS}^2} {\sigma^2} \thicksim \chi^2(N-p-1) σ2(Np1)σ^OLS2χ2(Np1)
    但是,得到这两个分布有啥用啊?那就且听我们下回分解。


    7. 新问题的提出

    至此,我们解决了模型的参数估计,以及所估计参数的相关性质,并给出了参数估计量的分布。但是,我们还有一堆问题还没有解决:

    (1)估计出来的模型好不好?有没有效呢?
    (2)估计出来的模型是否满足模型的基本假设条件呢?不满足可咋办呢?
    (3)如何利用估计出来的模型进行预测呢?预测的准确性有多高呢?
    (4)下次更新是不是又得一年后啊?

    在下一篇文章中,就让我们利用参数的区间估计与假设检验,来看看我们拟合出来的模型到底有没有用。


    参考文献

    [1] 线性回归估计相合性问题的新进展 .


    写在最后

    欢迎感兴趣的小伙伴来跟作者一起挑刺儿~ 包括但不限于语言上的、排版上的和内容上的不足和疏漏~ 一起进步呀!
    有任何问题,欢迎在本文下方留言,或者将问题发送至勘误邮箱: mikeysun_bugfix@163.com
    谢谢大家!

    展开全文
  • 针对一般带约束的最小二乘估计(OCLSE)在参数估计中处理复共线性的不足,引入随机线性约束,提出了约束-d估计方法。在均方误差(MSE)下,讨论了它的性质,得到了3个主要结果,与带约束的最小二乘估计OCLSE、约束岭估计...
  • 基于改进的粒子群算法的多元线性回归模型参数估计.pdf
  • 多元线性回归模型的参数估计PPT课件.pptx
  • 讨论多元线性回归模型由极小化问题的解定义的M-估计β_n的强相合性,其中X_i为m×p矩阵.证明了:无论为随机向量((VecX) ̄′,Y′)的独立同分布观察向量还是X_i为已知的m×p设计阵,在适当的条件下β_n都是参数真值β_0的...
  • 多元线性回归模型及参数估计精PPT学习教案.pptx
  • 多元线性回归参数求解是一个矩阵求导的过程,所以需要知道一些矩阵运算、求导运算的公式: 然后对多元线性回归的损失函数进行求导,公式如下:(其中w、y、X都是矩阵) 令其为0: 其中倒数最后一步...

    前言:本文主要内容:
    1. 最小二乘法损失函数求解推导;
    2. sklearn中linear_model.LinearRegression参数介绍+案例

    最小二乘法损失函数求解推导

    最小二乘法的思路:对损失函数求导,令其为0,求得损失函数最小值时的参数,但前提条件:导数为凸函数。


    多元线性回归参数求解是一个矩阵求导的过程,所以需要知道一些矩阵运算、求导运算的公式:

    (A-B)^T=A^T-B^T;(AB)^T=B^T*A^T

    \frac{\partial{a}}{\partial{A}}=0;\frac{\partial{A^TB^TC}}{\partial{A}}=B^TC;\frac{\partial{C^TBA}}{\partial{A}}=B^TC;\frac{\partial{A^TBA}}{\partial{A}}=(B+B^T)A


     然后对多元线性回归的损失函数进行求导,公式如下:(其中w、y、X都是矩阵) 

    \begin{aligned} \frac{\partial{RSS}}{\partial{w}}&=\frac{\partial{||y-Xw||_2^2}}{\partial{w}}\\ &=\frac{\partial{(y-Xw)^T(y-Xw)}}{\partial{w}}\\ &=\frac{\partial{(y^T-w^TX^T)(y-Xw)}}{\partial{w}}\\ &=\frac{\partial{(y^Ty-w^TX^Ty-y^TXw+w^TX^TXw)}}{\partial{w}}\\ &=0-X^Ty-X^Ty+2X^TXw\\ &=X^TXw-X^Ty \end{aligned}

    令其为0: 

    \bg_white \begin{aligned} X^TXw-X^Ty&=0\\ X^TXw&=X^Ty\\ w&=(X^TX)^{-1}X^Ty \end{aligned}

    其中倒数最后一步需要左乘(X^TX)^-1。


     补充点:

    1. 逆矩阵存在的充分必要条件:特征矩阵不能存在多重共线性。

    2. 因为X^TX随着特征、数据量的增加,计算量很大,不能直接使用numpy计算,需要使用奇异值分解来计算。

    奇异值分解在PCA中也涉及使用。

    在PCA中,为了找到新的特征空间,需要使用奇异值分解来计算。 在sklearn.decomposition.PCA中,使用参数svd_solver来控制奇异值计算方法,可以使用full(精准)/arpack(截断)/randomized(随机)/auto。

    3. 最小二乘法求解线性回归属于“无偏估计”,即要求标签必须是正态分布,所以需要对y进行正态处理,可以使用quantileTransformer或者PowerTransformer。不过因为机器学习为“后验方式”,如果不对y进行处理的效果好,那么也不必正态处理。


    sklearn中linear_model.LinearRegression参数介绍+案例:加利福尼亚房价数据集

    • fit_intercept:是否有截距项,默认True;
    • normalize:对X进行标准化,默认True;
    • copy_X:是否复制一份,以保证不在原始数据上操作,默认True;
    • n_jobs:设定运算性能,填整数,默认None;
    • model.coef_ :查看系数,即w
    • model.intercept_ :查看截距项
    from sklearn.datasets import fetch_california_housing as fch
    
    housevalue=fch()
    
    X=pd.DataFrame(housevalue.data,columns=housevalue.feature_names)
    X.head()
    y=housevalue.target
    
    from sklearn.model_selection import train_test_split
    X_train,X_test,y_train,y_test=train_test_split(X,y,test_size=0.3,random_state=9)
    
    from sklearn.linear_model import LinearRegression
    model=LinearRegression().fit(X_train,y_train)
    y_hat=model.predict(X_test)
    y_hat
    
    model.coef_ # 返回模型的参数
    [*zip(X_train.columns,model.coef_)]
    
    model.intercept_ # 返回模型的截距项

    展开全文
  • 在上一篇文章中,我们分别研究了最小二乘估计量 β^OLS 和 σ^OLS 的相关性质,证明了 β^OLS 是 β 的一个最优线性无偏估计量(BLUE), σ^2OLS 是 σ2 的一个无偏估计量,并得到了其在正态性误差假设下所对应的...
  • 如何求解多元线性回归残差平方和RSS最小化的参数向量?这种通过最小化真实值和预测值之间的RSS来求解 参数方法叫做最小二乘法。 求解极值的第一步往往是求解一阶导数并让一阶导数等于0,最小二乘法也不能免俗。因此...
  • 最近忙于毕设的事情,有很长一段时间没有写笔记了,近段时间学习上需要用到一些回归模型的知识,此条笔记用来记录学习笔记,声明:参考视频来源于李进华博士,大家可以去搜他的视频,讲解深入浅出,非常到位。...
  • 多元线性回归算法

    千次阅读 2021-10-25 18:39:30
    文章目录一、概念二、EXCEL的多元线性回归三、代码实现多元线性回归1.sklearn包实现2.线性回归模型的统计学库实现四、总结参考链接 一、概念 在回归分析中,如果有两个或两个以上的自变量,就称为多元回归。事实上...
  • matlab多元线性回归

    2014-09-22 13:41:39
    matlab多元线性回归
  • 多元线性回归

    千次阅读 2019-12-29 21:45:18
    一、多元线性回归 所谓的多元线性回归就是指在回归分析中,如果有两个或两个以上的自变量,就称为多元回归。 二、多元线性回归模型 1.建立模型 以二元线性回归模型为例 ,二元线性回归模型如下: 类似的使用最小...
  • 多元线性回归分析

    2018-06-27 20:55:33
    2、掌握一元线性回归的最小二乘法参数估计的计算公式、性质和应用; 3、理解拟合优度指标:决定系数R2的含义和作用; 4、掌握解释变量 和被解释变量 之间线性关系检验,回归参数 和 的显著性检验 5、了解利用回归...
  • 设随机变量yyy与一般变量x1,x2,...,xpx_{1},x_{2},...,x_{p}x1​,x2​,...,xp​的线性回归模型为: y= y= y= 样本(x,y)(x,y)(x,y)可由y=β0+β1x+εy=\beta _{0}+\beta _{1}x+\varepsilony=β0​+β1​x+ε 表示,...
  • 多元线性回归模型检验方法

    万次阅读 2019-08-10 22:07:21
    终于找到一篇全面而又简洁的讲多元线性回归模型检验方法的文章 PDF下载地址 链接:https://pan.baidu.com/s/1UbyZcMC1VRTmlCEaX4Vybg 提取码:g481 具体内容 一、经济意义检验 经济意义检验主要检验模型参数估计量在...
  • 每个人的生活中都有一个点,那就是该人希望购买或出售房屋。 首先考虑一个人需要买房的情况。... 用户将能够输入他们想要购买的房屋类型,并在机器学习的帮助下,房屋价格预测器将显示所需房屋的估计价格。
  • Bootstrapreg 是一个 MATLAB 程序,用于引导多元线性回归模型 y = beta0 + beta1*x1 + ... + betan*xn + error。 重采样残差方法用于生成随机数以反复更新预测。 这种通过重新采样学生化残差的方法比引导对象方法更...
  • 方法采用传统的信号—传播模型,利用多元线性回归法求出该模型中的相应参数,从而利用该模型求出初始各移动终端的位置;再采用快速迭代法和各移动终端相互校正从而提高定位的精度。实验结果表明,该算法时间复杂度...
  • 分别采用最小二乘法和梯度下降算法估计线性回归模型参数

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 10,919
精华内容 4,367
关键字:

多元线性回归参数估计方法