概率统计Python计算：随机变量的线性回归

若随机向量$(X,Y)$存在相关系数${\rho }_{XY}$，由
$$\begin{gathered} E[(Y-(aX+b){)}^2]=D(Y)+a^{2}D(X)-2a\text{Cov}(Y,X)+(E(Y)-aE(X)-b{)}^2 \\ =D(Y)(1-{\rho }_{XY}^2)+D(X){\left(a-{\rho }_{XY}\sqrt{\frac{D(Y)}{D(X)}}\right)}^2+(E(Y)-aE(X)-b{)}^2 \end{gathered}$$
当
$$\{\begin{array}{l}a={\rho }_{XY}\sqrt{\frac{D(Y)}{D(X)}}\\ b=E(Y)-aE(X)\end{array}$$
时，即
$$\{\begin{array}{l}a={\rho }_{XY}\sqrt{\frac{D(Y)}{D(X)}}\\ b=E(Y)-{\rho }_{XY}\sqrt{\frac{D(Y)}{D(X)}}E(X)\end{array}$$
时$E[(Y-(aX+b){)}^2]$最小，等于$D(Y)(1-{\rho }_{XY}^2)$。此时，称$aX+b$为$Y$对于$X$的线性回归。$Y$对$X$的线性回归$aX+b$中系数$a$、$b$的取值${\rho }_{XY}\sqrt{\frac{D(Y)}{D(X)}}$和$E(Y)-{\rho }_{XY}\sqrt{\frac{D(Y)}{D(X)}}E(X)$称为$Y$对$X$的回归系数。相仿地，$X$对$Y$的线性回归为$aY+b$，其中
$$\{\begin{array}{l}a={\rho }_{XY}\sqrt{\frac{D(X)}{D(Y)}}\\ b=E(X)-{\rho }_{XY}\sqrt{\frac{D(X)}{D(Y)}}E(Y)\end{array}.$$
例1 设随机向量$(X,Y)$的联合分布律为

计算$Y$对于$X$的线性回归。
解： 先计算出$X$和$Y$的边缘分布：

和

从而算得$E(X)=\frac{8}{3}$，$E(Y)=3$，$D(X)=2.22$，$D(Y)=2$。
又$E(XY)=8.125$，于是$\text{Cov}(X,Y)=E(XY)-E(X)E(Y)=0.125$。进而${\rho }_{XY}=\frac{\text{Cov}(X,Y)}{\sqrt{D(X)}\sqrt{D(Y)}}=0.593$。最终得到
$$\{\begin{array}{l}a={\rho }_{XY}\sqrt{\frac{D(Y)}{D(X)}}=0.056\\ b=E(Y)-aE(X)=2.85\end{array}.$$
即，$Y$对于$X$的线性回归为$0.056X+2.85$。
下列代码验算本例的计算结果。

import numpy as np                                      #导入numpy
X=np.array([1, 2, 3, 4, 5])                             #X的取值
Y=np.array([1, 2, 3, 4, 5])                             #Y的取值
Pxy=np.array([[1/12, 1/24, 1/12, 1/12, 1/24],           #联合分布律的概率矩阵
             [1/24, 1/24, 1/24, 0, 1/24],
             [0, 1/24, 1/24, 1/24, 1/24],
             [1/24, 1/24, 0, 1/24, 1/24],
             [1/30, 1/30, 1/30, 1/30, 1/30]])
Ex=expect(Pxy, X)                                       #E(X)
Ex2=expect(Pxy, X, func=lambda x, y: x*x)               #E(X^2)
sigmax=np.sqrt(Ex2-Ex**2)                               #X的标准差
Ey=expect(Pxy,Yv=Y, func=lambda x,y:y)                  #E(Y)
Ey2=expect(Pxy, Yv=Y, func=lambda x, y: y*y)            #E(Y^2)
sigmay=np.sqrt(Ey2-Ey**2)                               #Y的标准差
Exy=expect(Pxy, X, Y, lambda x, y:x*y)                  #E(XY)
rho=rhoxy(Exy, Ex, Ey, sigmax, sigmay)                  #X，Y的相关系数
a=rho*sigmay/sigmax                                     #回归系数a
b=Ey-a*Ex                                               #回归系数b
print('Y=%.3f*X+%.3f'%(a, b))                           #Y对X的线性回归

借助程序中各行注释，不难理解代码意义。其中，第9、10、12、13和15调用计算离散型随机变量期望的函数expect（定义见博文《离散型自定义分布数学期望的计算》），分别计算$E(X)$，$E(X^2)$，$E(Y)$，$E(Y^2)$及$E(XY)$。第16行调用计算随机变量$(X,Y)$的相关系数${\rho }_{XY}$的函数rhoxy（定义见博文《协方差与相关系数计算》）运行程序，输出

Y=0.056*X+2.850

写博不易，敬请支持：
如果阅读本文于您有所获，敬请点赞、评论、收藏，谢谢大家的支持！
代码诚可贵，原理价更高。若为AI学，读正版书好。
返回《导引》