需要代码的可以直接跳到最后一章，我们这里会进行一些数学推导

一、相关系数

1.1、相关系数的计算

相关系数是对变量之间关系的密切程度的度量，对两个变量之间线性相关程度的度量称为简单相关系数，若相关系数是根据全部数据计算的，称为总体相关系数记为 ρ\rhoρ，若是根据样本数据计算的则称为样本先关系数记为rrr。

样本相关系数的计算公式
r=∑(xi−xˉ)(yi−yˉ)∑(xi−xˉ)2⋅(yi−yˉ)2=lxylxx⋅lyyr = \frac{{\sum {({x_i} - \bar x)({y_i} - \bar y)} }}{{\sqrt {{{\sum {{{({x_i} - \bar x)}^2} \cdot ({y_i} - \bar y)} }^2}} }} = \frac{{{l_{xy}}}}{{\sqrt {{l_{xx}} \cdot {l_{yy}}} }} r=∑(xi​−xˉ)2⋅(yi​−yˉ​)2​∑(xi​−xˉ)(yi​−yˉ​)​=lxx​⋅lyy​​lxy​​
r的取值范围是-1到1

r的绝对值为1，完全相关

r=0，不相关（线性），但是不能说明变量之间没有线性关系

r<0，负相关。r>0，正相关

r的数值与量纲无关，r也不表示因果关系

1.2、相关系数的显著性检验

检验两个总体变量之间是否存在线性相关性，采用t检验，提出假设 H0:ρ=0;H1:ρ≠0;H_0:\rho = 0;H_1: \rho \ne 0;H0​:ρ=0;H1​:ρ​=0;计算检验的统计量：
T=rn−21−r2∼t(n−2)T=\frac {r \sqrt {n-2}}{\sqrt {1-r^2}} \sim t(n-2) T=1−r2​rn−2​​∼t(n−2)
确定显著性水平α\alphaα,并作出决策，若∣T∣>tα|T|>t_{\alpha}∣T∣>tα​则拒绝原假设。若∣T∣

二、一元线性回归

在相关分析中，变量x和变量y处于同等地位，在回归分析中，变量y称为因变量，处在被解释的地位，x称为自变量，用于预测因变量的变化。回归模型分为一元回归和多元回归，一元回归分为线性回归和非线性回归，多元回归也分为线性回归和非线性回归。一元回归只有一个自变量，多元回归有多个自变量。

对于只涉及一个自变量的简单线性回归模型可以表示为
y=β0+β1x+εy = {\beta _0} + {\beta _1}x + \varepsilon y=β0​+β1​x+ε
误差项ε\varepsilonε是随机变量，反映了除x和y之间的线性关系之外的随机因素对y的影响，是不能由x和y之间的线性关系所解释的变异性。

1.1、一元线性回归的基本假设

E(ε)=0;E(y∣x)=β0+β1x;E(\varepsilon)=0;E(y|x)=\beta _ 0 + \beta _ 1 x; E(ε)=0;E(y∣x)=β0​+β1​x;

对于所有的x值，ε\varepsilonε的方差σ2\sigma ^ 2σ2都相同，误差项ε\varepsilonε是一个服从正态分布的随机变量，且相互独立，即ε∼N(σ2)\varepsilon\sim N(\sigma^2)ε∼N(σ2)

1.2、回归方程

描述y的平均值或者期望值如何依赖于x的方程称为回归方程，回归方程的形式如下
E(y)=β0+β1xE(y)=\beta_0+\beta_1x E(y)=β0​+β1​x
总体回归参数β0,β1\beta_0 , \beta_1β0​,β1​是未知的，需要利用样本数据去估计，用样本统计量β^0{\hat\beta}_0β^​0​β^1{\hat\beta}_1β^​1​去估计β0,β1\beta_0 , \beta_1β0​,β1​就得到了估计的回归方程。
E(y)=β^0+β^1xE(y)=\hat\beta_0+\hat\beta_1x E(y)=β^​0​+β^​1​x

1.3、参数的最小二乘估计

最小二乘法是使得因变量的观察值与估计值之间的离差平方和达到最小来求参数的方法
{β1=∑i=1n(xi−xˉ)(yi−yˉ)∑i=1n(xi−xˉ)2=lxylxxβ0=yˉ−β^1xˉ\left\{ \begin{array}{l} {\beta _1} = \frac{{\sum\limits_{i = 1}^{n} {({x_i} - \bar x)({y_i} - \bar y)} }}{{\sum\limits_{i=1}^n}(x_i-\bar x) ^ 2} = \frac {l_{xy}} {l_{xx}}\\ {\beta _0} = {\bar y } - {\hat \beta} _ 1 {\bar x} \end{array} \right. ⎩⎪⎪⎨⎪⎪⎧​β1​=i=1∑n​(xi​−xˉ)2i=1∑n​(xi​−xˉ)(yi​−yˉ​)​=lxx​lxy​​β0​=yˉ​−β^​1​xˉ​

1.4、离差平方和的分解

因变量y的取值是不同的，y取值的变异性来源于两个方面，自变量x不同和其他因素的影响（测量误差等）。

通过一个具体的观测值，变异的大小可以通过该实际观测值与其均值之差y−yˉy-\bar yy−yˉ​来表示
SST=∑i=1n(yi−yˉ)2SSR=∑i=1n(y^i−yˉ)2SSE=∑i=1n(yi−y^)2SST=SSR+SSE{SST = \sum\limits_{i=1}^n (y_i-\bar y) ^ 2 }\\ {SSR = \sum\limits_{i=1}^n (\hat y_i-\bar y) ^ 2 }\\ {SSE = \sum\limits_{i=1}^n (y_i-\hat y) ^ 2 }\\ {SST = SSR + SSE} SST=i=1∑n​(yi​−yˉ​)2SSR=i=1∑n​(y^​i​−yˉ​)2SSE=i=1∑n​(yi​−y^​)2SST=SSR+SSE

SSR是总变异平方和，SSR是回归平方和，SSE是残差平方和

1.5、判定系数

判定系数是回归平方和占总离差平方和的比例
r2=SSRSST=∑i=1n(y^i−yˉ)2∑i=1n(yi−yˉ)2r^2 = \frac {SSR} {SST} = \frac {\sum\limits_{i=1}^{n} {(\hat y_i - \bar y)^2}} {\sum\limits_{i=1}^{n} {( y_i - \bar y)^2}} r2=SSTSSR​=i=1∑n​(yi​−yˉ​)2i=1∑n​(y^​i​−yˉ​)2​
反映了回归直线的拟合程度，取值范围在-1到1之间，r2→1r_2 \to 1r2​→1回归方程的拟合效果更好，判定系数等于相关系数的平方（一元线性回归）

1.6、估计标准误差

估计标准误差是实际观察值与回归估计值离差平方和的均方根，反映实际观察值在回归直线附近的分散情况，从另一个角度说明了回归直线的拟合程度，计算公式为
Sy=∑i=1n(yi−y^i)2n−2=SSEn−2S_y=\sqrt {\frac {\sum\limits_{i=1}^{n} (y_i - \hat y_i) ^ 2} {n-2}} = \sqrt {\frac {SSE} {n-2}} Sy​=n−2i=1∑n​(yi​−y^​i​)2​​=n−2SSE​​

1.7、回归方程的显著性检验

检验自变量和因变量之间的线性关系是否显著，具体方法是将SSR同SSE比较，用F检验来分析两者之间的差别是否显著，如果显著，两变量之间不存在线性关系，如果不显著，两个变量之间存在线性关系。提出假设H0:H_0:H0​:线性关系不显著；H1：H_1：H1​：线性关系显著
F=MSRMSE=SSR/1SSE/(n−2)=(∑i=1n(y^i−yˉ)2)/1(∑i=1n(yi−yˉ)2)/(n−2)∼F(1,n−2)F = \frac {MSR} {MSE} = \frac {SSR/1} {SSE/(n-2)} = \frac {(\sum\limits_{i=1}^{n}{(\hat y_i - \bar y)^2})/1 } {(\sum\limits_{i=1}^{n}{( y_i - \bar y)^2})/(n-2) } \sim F(1,n-2) F=MSEMSR​=SSE/(n−2)SSR/1​=(i=1∑n​(yi​−yˉ​)2)/(n−2)(i=1∑n​(y^​i​−yˉ​)2)/1​∼F(1,n−2)
确定显著性水平α\alphaα，并根据自由度查出临界值FαF_\alphaFα​，做出决策：若F>FαF>F_\alphaF>Fα​，拒绝H0H_0H0​;若F

1.8、回归系数的显著性检验

β^1\hat\beta_1β^​1​是根据最小二乘法求出的样本统计量，他有自己的分布。具有如下性质：

分布：正态分布

数学期望：D(β^1)=β1D(\hat\beta_1)=\beta_1D(β^​1​)=β1​

标准差：σβ^1=σ∑(si−xˉ)2\sigma_{\hat\beta_{1}} = \frac {\sigma} {\sqrt {\sum {(s_i - \bar x) ^ 2}}}σβ^​1​​=∑(si​−xˉ)2​σ​

由于σ\sigmaσ未知，需要用估计量SyS_ySy​来代替得到β^1\hat \beta _ 1β^​1​的估计的标准差σβ^1=Sy∑(si−xˉ)2\sigma_{\hat\beta_{1}} = \frac {S_y} {\sqrt {\sum {(s_i - \bar x) ^ 2}}}σβ^​1​​=∑(si​−xˉ)2​Sy​​

提出假设 H0:β1=0(没有线性关系)；H1:β1≠0(有线性关系)H_0:\beta_1 = 0(没有线性关系)；H_1:\beta_1 \ne 0 (有线性关系)H0​:β1​=0(没有线性关系)；H1​:β1​​=0(有线性关系)

计算检验的统计量 T=β^1Sβ^1∼t(n−2)T=\frac {\hat \beta _ 1} {S _ {\hat \beta _ 1}} \sim t(n-2)T=Sβ^​1​​β^​1​​∼t(n−2)

确定显著性水平α\alphaα，并根据自由度查出临界值tα/2t_{\alpha / 2}tα/2​，做出决策：若t>tα/2t>t_{\alpha / 2}t>tα/2​，拒绝H0H_0H0​;若t

1.9、点估计

平均值的点估计 y^0=β^0+β^1⋅xˉ\hat y_0 = \hat \beta _ 0 + \hat \beta_1 \cdot \bar xy^​0​=β^​0​+β^​1​⋅xˉ

个别值的点估计 y^0=β^0+β^1⋅x1\hat y_0 = \hat \beta _ 0 + \hat \beta_1 \cdot x_1y^​0​=β^​0​+β^​1​⋅x1​

1.10、区间估计

y的平均值的置信区间估计，E(y0)E(y_0)E(y0​)在1−α1-\alpha1−α置信水平下的置信区间
y^0±tα/2(n−2)Sy1n+(x0−xˉ)2∑i=1n(xi−xˉ)2\hat y _ 0 \pm t_{\alpha / 2}(n-2) S_y \sqrt{\frac{1} {n} + \frac {(x_0 - \bar x)^2} {\sum\limits_{i=1}^{n}(x_i-\bar x)^2}} y^​0​±tα/2​(n−2)Sy​n1​+i=1∑n​(xi​−xˉ)2(x0​−xˉ)2​​
y的个别值的预测区间估计，对于自变量xxx的一个给定值x0x_0x0​，求出因变量yyy的一个个别值的估计区间，这一区间成为预测区间
y^0±tα/2(n−2)Sy1+1n+(x0−xˉ)2∑i=1n(xi−xˉ)2\hat y _ 0 \pm t_{\alpha / 2}(n-2) S_y \sqrt{1 + \frac{1} {n} + \frac {(x_0 - \bar x)^2} {\sum\limits_{i=1}^{n}(x_i-\bar x)^2}} y^​0​±tα/2​(n−2)Sy​1+n1​+i=1∑n​(xi​−xˉ)2(x0​−xˉ)2​​
影响区间宽度的水平

区间水平随着置信水平的增大而增大，区间宽度随离散程度的增大而增大，区间宽度随样本容量的增大而减小，区间宽度随着x0x_0x0​与xˉ\bar xxˉ的差异程度的增大而增大。

1.11、残差分析

残差 ei=yi−y^i=yi−β^0−β^1xie_i = y_i - \hat y_i = y_i-\hat \beta_0 - \hat \beta_1 x_iei​=yi​−y^​i​=yi​−β^​0​−β^​1​xi​

误差项 εi=yi−β0−β1xi\varepsilon _ i = y_i- \beta_0 - \beta_1 x_iεi​=yi​−β0​−β1​xi​

残差是误差项的估计值

残差标准化(残差除以其标准差后得到的数值) zei=eiSe=yi−y^iMSEz_{e_i} = \frac {e_i} {S_e} = \frac {y_i - \hat y_i} {\sqrt {MSE}}zei​​=Se​ei​​=MSE​yi​−y^​i​​

如果ε\varepsilonε服从正态分布，那么标准化残差也应该服从正态分布

三、多元线性回归

多元线性回归模型可以表示为
y=β0+β1x1+β2x2+...+βpxp+εy=\beta_0 + \beta_1 x_1 + \beta_2 x_2 + ... + \beta_p x_p + \varepsilon y=β0​+β1​x1​+β2​x2​+...+βp​xp​+ε
自变量可观测，随机误差项εi\varepsilon _ iεi​的期望值为0，且方差σ2\sigma ^ 2σ2均相同，是一个服从正态分布的随机变量，且相互独立

3.1、最小二乘法

使因变量的观察值与估计值之间的离差平方和最小来求
Q(β^0,β^1,β^1,...,β^0,)=∑i=1n(yi−y^i)2=∑i=1nei2Q({\hat \beta _ 0} , {\hat \beta _ 1} ,{\hat \beta _ 1} , ... ,{\hat \beta _ 0} ,) = \sum\limits_{i=1}^{n}(y_i - \hat y_i)^2 = \sum\limits_{i=1}^{n} e_{i}^{2} Q(β^​0​,β^​1​,β^​1​,...,β^​0​,)=i=1∑n​(yi​−y^​i​)2=i=1∑n​ei2​
根据最小二乘法可以求解各回归参数
KaTeX parse error: Unknown column alignment: 1 at position 24: … \begin{array}{1̲} {{\frac {\par…

3.2、回归方程的显著性检验

多重判定系数R2R^2R2
R2=SSRSST=∑i=1n(y^i−yˉ)2∑i=1n(yi−yˉ)2R^2 = \frac {SSR} {SST} = \frac {\sum\limits_{i=1}^{n} {(\hat y_i - \bar y)^2}} {\sum\limits_{i=1}^{n} {( y_i - \bar y)^2}} R2=SSTSSR​=i=1∑n​(yi​−yˉ​)2i=1∑n​(y^​i​−yˉ​)2​
取值在0-1之间，R2→1R^2 \to 1R2→1说明回归方程拟合好

修正的多重判定系数R2R^2R2
R修2=1−(1−R2)×n−1n−p−1R_{修}^{2} = 1 - (1 - R^2) \times \frac {n-1} {n-p-1} R修2​=1−(1−R2)×n−p−1n−1​
提出假设H0:H_0:H0​:线性关系不显著；H1：H_1：H1​：线性关系显著
F=MSRMSE=SSR/pSSE/(n−p−1)=(∑i=1n(y^i−yˉ)2)/p(∑i=1n(yi−yˉ)2)/(n−p−1)∼F(p,n−p−1)F = \frac {MSR} {MSE} = \frac {SSR/p} {SSE/(n-p-1)} = \frac {(\sum\limits_{i=1}^{n}{(\hat y_i - \bar y)^2})/p } {(\sum\limits_{i=1}^{n}{( y_i - \bar y)^2})/(n-p-1) } \sim F(p,n-p-1) F=MSEMSR​=SSE/(n−p−1)SSR/p​=(i=1∑n​(yi​−yˉ​)2)/(n−p−1)(i=1∑n​(y^​i​−yˉ​)2)/p​∼F(p,n−p−1)
确定显著性水平α\alphaα，并根据自由度查出临界值FαF_\alphaFα​，做出决策：若F>FαF>F_\alphaF>Fα​，拒绝H0H_0H0​;若F

3.3、回归系数的显著性检验

计算检验的统计量 T=β^1Sβ^1∼t(n−p−1)T=\frac {\hat \beta _ 1} {S _ {\hat \beta _ 1}} \sim t(n-p-1)T=Sβ^​1​​β^​1​​∼t(n−p−1)

确定显著性水平α\alphaα，并根据自由度查出临界值tα/2t_{\alpha / 2}tα/2​，做出决策：若∣T∣>tα/2|T|>t_{\alpha / 2}∣T∣>tα/2​，拒绝H0H_0H0​;若∣T∣

四、非线性回归

几种常见的非线性模型

指数函数 y=αeβxy=\alpha e ^ {\beta x}y=αeβx 两边取对数线性化

幂函数 y=axβy=ax^{\beta}y=axβ 两边取对数线性化

双曲线函数 y=xαx+βy=\frac {x} {\alpha x + \beta}y=αx+βx​

对数函数 y=α+βlnxy=\alpha + \beta ln {x}y=α+βlnx

S型曲线y=1α+βe−xy = \frac {1} {\alpha + \beta e^{-x}}y=α+βe−x1​

代码

代码需要安装的库

linux系统下

sudo pip install numpy
sudo pip install matplotlib

windows系统下

pip install numpy
pip install matplotlib

代码如下

# 最小二乘拟合
import numpy as np
import matplotlib.pyplot as plt# 需要拟合的两个变量
x = np.array([1,2,3,4,5])
y = np.array([1,2,4,9,15])
# 需要多少阶拟合
m = 10# 生成系数矩阵A
def gen_coefficient_matrix(X):A = []# 计算每一个方程的系数for i in range(m):a = []# 计算当前方程中的每一个系数for j in range(m):a.append(sum(X ** (i+j)))A.append(a)return A# 计算方程组的右端向量b
def gen_right_vector(X, Y):b=[]for i in range(m):b.append(sum(X**i*y))return bA = gen_coefficient_matrix(x)
b = gen_right_vector(x,y)res = np.linalg.solve(A, b)
print(f"res is {res}")x_ = [x0 / 10 for x0 in range(1 ,max(x)*10+1)]
z = []for i in x_:temp = 0for j in range(len(res)):temp = temp + res[j] * (i) ** jz.append(temp) plt.figure(1)
plt.plot(x,y,'ro')
plt.plot(x_,z,'-')
plt.show()