【机器学习】基于线性回归的医疗费用预测模型

文章目录

    • 一、线性回归
      • 定义和工作原理
      • 假设表示
    • 二、导入库和数据集
      • 矩阵表示
      • 可视化
    • 三、成本函数
      • 向量的内积
    • 四、正态方程
    • 五、探索性数据分析
      • 描述性统计
      • 检查缺失值
      • 数据分布图
        • 相关性热图
        • 保险费用分布
        • 保险费用与性别和吸烟情况的关系
        • 保险费用与子女数量的关系
        • 保险费用与地区和性别的关系
        • 保险费用与年龄和BMI的关系
    • 六、数据预处理
      • 编码
      • Box-Cox 变换
    • 七、训练集和测试集划分
    • 八、模型构建
    • 九、模型评估
      • 使用正态方程进行预测
      • 使用 Scikit-Learn 进行预测
    • 十、模型验证
      • 具体验证步骤
      • 结果分析

建立一个线性回归模型来预测个人的医疗费用。该数据集包括以下特征:年龄、性别、BMI(身体质量指数)、子女数量、是否吸烟和地区。这些特征是独立变量,费用是依赖变量。模型的目标是预测健康保险收取的个人医疗费用。

一、线性回归

定义和工作原理

线性回归是一种监督学习算法,适用于目标变量(因变量)是连续实数的情况。它通过最佳拟合线来建立因变量 y y y 与一个或多个自变量 x x x 之间的关系。

线性回归基于普通最小二乘法(OLS)或均方误差(MSE)的原理。在统计学中,OLS 是一种估计线性回归函数未知参数的方法,其目标是最小化给定数据集中观察到的因变量与线性回归函数预测的因变量之间的平方差和。

假设表示

我们用 x i \mathbf{x_i} xi 表示独立变量,用 y i \mathbf{y_i} yi 表示因变量。一对训练示例表示为 ( x i , y i ) \mathbf{(x_i, y_i)} (xi,yi)。其中, i \mathbf{i} i 是训练集的索引,我们有 m \mathbf{m} m 个训练示例,即 i = 1 , 2 , 3 , . . . , m \mathbf{i = 1, 2, 3, ..., m} i=1,2,3,...,m

监督学习的目标是学习给定训练集的假设函数 h \mathbf{h} h,该函数可以根据 x \mathbf{x} x 估计 y \mathbf{y} y。假设函数表示为:
h θ ( x i ) = θ 0 + θ 1 x i \mathbf{ h_\theta(x_i) = \theta_0 + \theta_1 x_i } hθ(xi)=θ0+θ1xi
其中, θ 0 \mathbf{\theta_0} θ0 θ 1 \mathbf{\theta_1} θ1 是假设参数。这是简单/单变量线性回归的方程。

对于多元线性回归,如果存在多个独立变量,我们用 x i j \mathbf{x_{ij}} xij 表示独立变量,用 y i \mathbf{y_i} yi 表示因变量。有 n \mathbf{n} n 个独立变量,则 j = 1 , 2 , 3 , . . . , n \mathbf{j = 1, 2, 3, ..., n} j=1,2,3,...,n。假设函数表示为:
h θ ( x i ) = θ 0 + θ 1 x i 1 + θ 2 x i 2 + . . . + θ j x i j + . . . + θ n x i n \mathbf{ h_\theta(x_i) = \theta_0 + \theta_1 x_{i1} + \theta_2 x_{i2} + ... + \theta_j x_{ij} + ... + \theta_n x_{in} } hθ(xi)=θ0+θ1xi1+θ2xi2+...+θjxij+...+θnxin
其中, θ 0 , θ 1 , . . . , θ j , . . . , θ n \mathbf{\theta_0, \theta_1, ..., \theta_j, ..., \theta_n} θ0,θ1,...,θj,...,θn 为假设参数, m \mathbf{m} m 为训练样本数, n \mathbf{n} n 为独立变量数, x i j \mathbf{x_{ij}} xij 为第 j \mathbf{j} j 个特征的第 i \mathbf{i} i 个训练样本。

二、导入库和数据集

# 导入库
import pandas as pd  # 数据处理
import numpy as np  # 数据处理
import matplotlib.pyplot as plt  # 数据可视化
import seaborn as sns  # 数据可视化

# 设置图形参数
plt.rcParams['figure.figsize'] = [8, 5]
plt.rcParams['font.size'] = 14
plt.rcParams['font.weight'] = 'bold'
plt.style.use('seaborn-whitegrid')

# 导入数据集
path = '../linear-regression-tutorial/'
df = pd.read_csv(path + 'insurance.csv')

# 查看数据集的行数和列数
print('\nNumber of rows and columns in the data set: ', df.shape)
print('')

# 查看数据集的前几行
df.head()

数据集的形状为 ( 1338 , 7 ) (1338, 7) (1338,7),包含1338个样本和7个特征。具体特征如下:

在这里插入图片描述
现在我们有了导入的数据集。该数据集包含 m = 1338 \mathbf{m = 1338} m=1338 个训练示例和 n = 7 \mathbf{n = 7} n=7 个特征。目标变量是费用,其他六个变量(例如年龄、性别、BMI、子女数量、吸烟状态、地区)是独立变量。

由于有多个独立变量,我们需要拟合多元线性回归。假设函数如下:
h θ ( x i ) = θ 0 + θ 1 ⋅ a g e + θ 2 ⋅ s e x + θ 3 ⋅ b m i + θ 4 ⋅ c h i l d r e n + θ 5 ⋅ s m o k e r + θ 6 ⋅ r e g i o n \mathbf{ h_\theta(x_{i}) = \theta_0 + \theta_1 \cdot age + \theta_2 \cdot sex + \theta_3 \cdot bmi + \theta_4 \cdot children + \theta_5 \cdot smoker + \theta_6 \cdot region } hθ(xi)=θ0+θ1age+θ2sex+θ3bmi+θ4children+θ5smoker+θ6region
这是给定数据集的多元线性回归方程。例如,如果 i = 1 \mathbf{i = 1} i=1,则:
h θ ( x 1 ) = θ 0 + θ 1 ⋅ 19 + θ 2 ⋅ f e m a l e + θ 3 ⋅ 27.900 + θ 4 ⋅ 1 + θ 5 ⋅ y e s + θ 6 ⋅ s o u t h w e s t y 1 = 16884.92400 \mathbf{ h_\theta(x_{1}) = \theta_0 + \theta_1 \cdot 19 + \theta_2 \cdot female + \theta_3 \cdot 27.900 + \theta_4 \cdot 1 + \theta_5 \cdot yes + \theta_6 \cdot southwest } \\ \mathbf{ y_1 = 16884.92400 } hθ(x1)=θ0+θ119+θ2female+θ327.900+θ41+θ5yes+θ6southwesty1=16884.92400
如果 i = 3 \mathbf{i = 3} i=3,则:
h θ ( x 3 ) = θ 0 + θ 1 ⋅ 28 + θ 2 ⋅ m a l e + θ 3 ⋅ 33.000 + θ 4 ⋅ 3 + θ 5 ⋅ n o + θ 6 ⋅ n o r t h w e s t y 3 = 4449.46200 \mathbf{ h_\theta(x_{3}) = \theta_0 + \theta_1 \cdot 28 + \theta_2 \cdot male + \theta_3 \cdot 33.000 + \theta_4 \cdot 3 + \theta_5 \cdot no + \theta_6 \cdot northwest } \\ \mathbf{ y_3 = 4449.46200 } hθ(x3)=θ0+θ128+θ2male+θ333.000+θ43+θ5no+θ6northwesty3=4449.46200
注意:在 Python 中,索引从 0 开始。即 x 1 \mathbf{x_1} x1 表示为:
x 1 = ( 19 f e m a l e 27.900 1 n o n o r t h w e s t ) \mathbf{ x_1 = \left(\begin{matrix} 19 & female & 27.900 & 1 & no & northwest \end{matrix}\right) } x1=(19female27.9001nonorthwest)

矩阵表示

通常,可以将上述向量表示为:
x i j = ( x i 1 x i 2 ⋯ x i n ) \mathbf{ x_{ij} = \left( \begin{matrix} \mathbf{x_{i1}} & \mathbf{x_{i2}} & \cdots & \mathbf{x_{in}} \end{matrix} \right) } xij=(xi1xi2xin)
现在,我们将所有单个向量组合成一个 ( m , n ) (m, n) (m,n) 的输入矩阵 X \mathbf{X} X,其中包含所有训练示例:
X = ( x 11 x 12 ⋯ x 1 n x 21 x 22 ⋯ x 2 n x 31 x 32 ⋯ x 3 n ⋮ ⋮ ⋱ ⋮ x m 1 x m 2 ⋯ x m n ) ( m , n ) \mathbf{X} = \left( \begin{matrix} x_{11} & x_{12} & \cdots & x_{1n} \\ x_{21} & x_{22} & \cdots & x_{2n} \\ x_{31} & x_{32} & \cdots & x_{3n} \\ \vdots & \vdots & \ddots & \vdots \\ x_{m1} & x_{m2} & \cdots & x_{mn} \\ \end{matrix} \right)_{(m,n)} X= x11x21x31xm1x12x22x32xm2x1nx2nx3nxmn (m,n)
我们将函数参数和因变量表示为向量形式:

θ = ( θ 0 θ 1 ⋮ θ j ⋮ θ n ) ( n + 1 , 1 ) y = ( y 1 y 2 ⋮ y i ⋮ y m ) ( m , 1 ) \mathbf{\theta} = \left( \begin{matrix} \theta_0 \\ \theta_1 \\ \vdots \\ \theta_j \\ \vdots \\ \theta_n \end{matrix} \right)_{(n+1,1)} \mathbf{y} = \left( \begin{matrix} y_1 \\ y_2 \\ \vdots \\ y_i \\ \vdots \\ y_m \end{matrix} \right)_{(m,1)} θ= θ0θ1θjθn (n+1,1)y= y1y2yiym (m,1)
因此,假设函数可以表示为:
h θ ( x ) = X θ \mathbf{ h_\theta(x) = X\theta } hθ(x)=Xθ

可视化

为了可视化,我们将使用 seaborn 库,将 BMI 作为自变量,费用作为因变量进行拟合。

sns.lmplot(x='bmi', y='charges', data=df, aspect=2, height=6)
plt.xlabel('Boby Mass Index $(kg/m^2)$: as Independent variable')
plt.ylabel('Insurance Charges: as Dependent variable')
plt.title('Charge Vs BMI')
plt.show()

在这里插入图片描述

三、成本函数

成本函数用于衡量模型在估计 x x x y y y 之间关系的准确性。通过计算观察到的因变量与假设函数预测的因变量之间的平均差异,我们可以评估模型的表现。
J ( θ ) = 1 m ∑ i = 1 m ( y ^ i − y i ) 2 \mathbf{J(\theta) = \frac{1}{m} \sum_{i=1}^{m}(\hat{y}_i - y_i)^2} J(θ)=m1i=1m(y^iyi)2
为了实现线性回归,我们在训练样本中添加一个额外的列,即 x 0 x_0 x0 特征,其中 x 0 = 1 \mathbf{x_0=1} x0=1。于是输入矩阵 X \mathbf{X} X 变为:
X = ( x 10 x 11 x 12 ⋯ x 1 n x 20 x 21 x 22 ⋯ x 2 n x 30 x 31 x 32 ⋯ x 3 n ⋮ ⋮ ⋮ ⋱ ⋮ x m 0 x m 1 x m 2 ⋯ x m n ) ( m , n + 1 ) \mathbf{X} = \begin{pmatrix} x_{10} & x_{11} & x_{12} & \cdots & x_{1n} \\ x_{20} & x_{21} & x_{22} & \cdots & x_{2n} \\ x_{30} & x_{31} & x_{32} & \cdots & x_{3n} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ x_{m0} & x_{m1} & x_{m2} & \cdots & x_{mn} \\ \end{pmatrix}_{(m,n+1)} X= x10x20x30xm0x11x21x31xm1x12x22x32xm2x1nx2nx3nxmn (m,n+1)
其中 x i 0 = 1 \mathbf{x_{i0} = 1} xi0=1。现在我们可以将普通最小二乘成本函数以矩阵形式重写为:
J ( θ ) = 1 m ( X θ − y ) T ( X θ − y ) \mathbf{J(\theta) = \frac{1}{m} (X\theta - y)^T (X\theta - y)} J(θ)=m1(Xθy)T(Xθy)
输入矩阵 X \mathbf{X} X 的大小为 ( m , n + 1 ) \mathbf{(m,n+1)} (m,n+1),函数参数 θ \mathbf{\theta} θ 的大小为 ( n + 1 , 1 ) \mathbf{(n+1,1)} (n+1,1),因变量向量 y \mathbf{y} y 的大小为 ( m , 1 ) \mathbf{(m,1)} (m,1)。矩阵 X ( m , n + 1 ) θ ( n + 1 , 1 ) \mathbf{X_{(m,n+1)} \theta_{(n+1,1)}} X(m,n+1)θ(n+1,1) 的乘积将返回一个大小为 ( m , 1 ) \mathbf{(m,1)} (m,1) 的向量,然后 ( X θ − y ) ( 1 , m ) T ( X θ − y ) ( m , 1 ) \mathbf{(X\theta - y)^T_{(1,m)} (X\theta - y)_{(m,1)}} (Xθy)(1,m)T(Xθy)(m,1) 的乘积将返回一个标量。

向量的内积

为了说明为什么 ∑ i = 1 m ( y ^ i − y i ) 2 \sum_{i=1}^{m} (\hat{y}_i - y_i)^2 i=1m(y^iyi)2 等于 ( X θ − y ) T ( X θ − y ) (X\theta - y)^T (X\theta - y) (y)T(y),我们需要详细解释一下符号和操作。

首先,定义一些符号:

  • y \mathbf{y} y 是实际值的向量,维度为 m × 1 m \times 1 m×1
  • y ^ \mathbf{\hat{y}} y^ 是预测值的向量,维度为 m × 1 m \times 1 m×1,可以表示为 y ^ = X θ \mathbf{\hat{y} = X\theta} y^=Xθ,其中 X \mathbf{X} X 是特征矩阵,维度为 m × n m \times n m×n
  • e \mathbf{e} e 是误差向量,表示为 e = y ^ − y \mathbf{e = \hat{y} - y} e=y^y,维度为 m × 1 m \times 1 m×1
  • $\mathbf{e}^T \mathbf{e} = (X\theta - y)^T (X\theta - y) $

那么,误差向量的每个元素就是对应的预测值和实际值之间的差异:
e i = y ^ i − y i \mathbf{e_i} = \hat{y}_i - y_i ei=y^iyi
我们感兴趣的是误差的平方和:
∑ i = 1 m ( y ^ i − y i ) 2 \sum_{i=1}^{m} (\mathbf{\hat{y}_i} - \mathbf{y_i})^2 i=1m(y^iyi)2
现在来看矩阵运算。误差向量 e \mathbf{e} e 的转置是一个 1 × m 1 \times m 1×m 的行向量:
e T = [ e 1 , e 2 , … , e m ] \mathbf{e}^T = [\mathbf{e_1}, \mathbf{e_2}, \ldots, \mathbf{e_m}] eT=[e1,e2,,em]
e \mathbf{e} e 转置与 e \mathbf{e} e 相乘,得到一个标量:
e T e = [ e 1 e 2 … e m ] [ e 1 e 2 ⋮ e m ] = e 1 2 + e 2 2 + ⋯ + e m 2 \mathbf{e}^T \mathbf{e} = \begin{bmatrix} \mathbf{e_1} & \mathbf{e_2} & \ldots & \mathbf{e_m} \end{bmatrix} \begin{bmatrix} \mathbf{e_1} \\ \mathbf{e_2} \\ \vdots \\ \mathbf{e_m} \end{bmatrix} = \mathbf{e_1}^2 + \mathbf{e_2}^2 + \cdots + \mathbf{e_m}^2 eTe=[e1e2em] e1e2em =e12+e22++em2
这是因为矩阵乘法中,行向量和列向量的内积就是对应元素的乘积和。展开来看:
e T e = ∑ i = 1 m e i 2 = ∑ i = 1 m ( y ^ i − y i ) 2 \mathbf{e}^T \mathbf{e} = \sum_{i=1}^{m} \mathbf{e_i}^2 = \sum_{i=1}^{m} (\mathbf{\hat{y}_i} - \mathbf{y_i})^2 eTe=i=1mei2=i=1m(y^iyi)2
因此,可以看到 e T e \mathbf{e}^T \mathbf{e} eTe ∑ i = 1 m ( y ^ i − y i ) 2 \sum_{i=1}^{m} (\mathbf{\hat{y}_i} - \mathbf{y_i})^2 i=1m(y^iyi)2 是相等的。这表明通过矩阵运算得到的误差平方和与逐元素计算的误差平方和是一致的。

四、正态方程

正态方程是具有普通最小二乘成本函数的线性回归问题的解析解。为了最小化成本函数,对 J ( θ ) \mathbf{J(\theta)} J(θ) θ \mathbf{\theta} θ 的偏导数,并使其等于零。函数的导数表示当输入发生微小变化时,函数输出的变化情况。
m i n θ 0 , θ 1 , … , θ n J ( θ 0 , θ 1 , … , θ n ) \mathbf{min_{\theta_0, \theta_1, \ldots, \theta_n} J(\theta_0, \theta_1, \ldots, \theta_n)} minθ0,θ1,,θnJ(θ0,θ1,,θn)
对每个 θ j \mathbf{\theta_j} θj 求导数,使其等于零:
∂ J ( θ j ) ∂ θ j = 0 \mathbf{\frac{\partial J(\theta_j)}{\partial \theta_j} = 0} θj∂J(θj)=0
其中 j = 0 , 1 , 2 , … , n \mathbf{j = 0,1,2, \ldots, n} j=0,1,2,,n

现在我们将应用成本函数的偏导数公式:
∂ J ( θ j ) ∂ θ j = ∂ ∂ θ 1 m ( X θ − y ) T ( X θ − y ) \mathbf{\frac{\partial J(\theta_j)}{\partial \theta_j} = \frac{\partial}{\partial \theta} \frac{1}{m} (X\theta - y)^T (X\theta - y)} θj∂J(θj)=θm1(Xθy)T(Xθy)
我们丢弃 1 m \mathbf{\frac {1}{m}} m1 部分,因为我们将导数与零进行比较。接着我们求解 J ( θ ) \mathbf{J(\theta)} J(θ)
J ( θ ) = ( X θ − y ) T ( X θ − y ) = ( ( X θ ) T − y T ) ( X θ − y ) = ( θ T X T − y T ) ( X θ − y ) = θ T X T X θ − y T X θ − θ T X T y + y T y = θ T X T X θ − 2 θ T X T y + y T y \mathbf{J(\theta) = (X\theta - y)^T (X\theta - y)}\\ \mathbf{= ((X\theta)^T - y^T) (X\theta - y)}\\ \mathbf{= (\theta^T X^T - y^T) (X\theta - y)}\\ \mathbf{= \theta^T X^T X \theta - y^T X \theta - \theta^T X^T y + y^T y}\\ \mathbf{= \theta^T X^T X \theta - 2\theta^T X^T y + y^T y} J(θ)=(Xθy)T(Xθy)=((Xθ)TyT)(Xθy)=(θTXTyT)(Xθy)=θTXTXθyTXθθTXTy+yTy=θTXTXθ2θTXTy+yTy
补充:矩阵运算的基本性质

  1. 向量的转置运算满足对称性,即:

( a T b ) = ( b T a ) \begin{equation} (\mathbf{a}^T \mathbf{b}) = (\mathbf{b}^T \mathbf{a}) \end{equation} (aTb)=(bTa)

  1. 矩阵乘法的转置满足:

( A B ) T = B T A T \begin{equation} (\mathbf{AB})^T = \mathbf{B}^T \mathbf{A}^T \end{equation} (AB)T=BTAT

所以
y T ( X θ ) = ( X θ ) T y = ( θ T X T ) y = θ T X T y \begin{align*} \mathbf{y}^T (\mathbf{X\theta}) & = (\mathbf{X\theta})^T \mathbf{y} \\ & = (\mathbf{\theta}^T \mathbf{X}^T) \mathbf{y} \\ & = \mathbf{\theta}^T \mathbf{X}^T \mathbf{y} \end{align*} yT(Xθ)=(Xθ)Ty=(θTXT)y=θTXTy
现在我们对 θ \mathbf{\theta} θ 求导数:
∂ J ( θ ) ∂ θ = ∂ ∂ θ ( θ T X T X θ − 2 θ T X T y + y T y ) = X T X ∂ θ T θ ∂ θ − 2 X T y ∂ θ T ∂ θ + ∂ y T y ∂ θ \mathbf{\frac{\partial J(\theta)}{\partial \theta} = \frac{\partial}{\partial \theta} (\theta^T X^T X \theta - 2\theta^T X^T y + y^T y)} \mathbf{= X^T X \frac {\partial \theta^T \theta}{\partial \theta} - 2 X^T y \frac{\partial \theta^T}{\partial \theta} + \frac{\partial y^T y}{\partial \theta}} θ∂J(θ)=θ(θTXTXθ2θTXTy+yTy)=XTXθθTθ2XTyθθT+θyTy
其中, ∂ x 2 ∂ x = 2 x \mathbf{\frac {\partial x^2}{\partial x} = 2x} ∂xx2=2x ∂ k x 2 ∂ x = k x \mathbf{\frac {\partial kx^2}{\partial x} = kx} ∂x∂kx2=kx ∂ C o n s t a c t ∂ x = 0 \mathbf{\frac {\partial Constact}{\partial x} = 0} ∂x∂Constact=0

  • x 2 x^2 x2 求导数得到 2 x 2x 2x
  • k x 2 kx^2 kx2 求导数得到 2 k x 2kx 2kx。(上面是把2k整合到k里面了)
  • 对常数函数求导数得到 0。

∂ J ( θ ) ∂ θ = X T X 2 θ − 2 X T y + 0 0 = 2 X T X θ − 2 X T y X T X θ = X T y θ = ( X T X ) − 1 X T y \mathbf{\frac{\partial J(\theta)}{\partial \theta} = X^T X 2\theta - 2X^T y + 0} \\ \mathbf{0 = 2X^T X \theta - 2X^T y}\\ \mathbf{X^T X \theta = X^T y}\\ \mathbf{\theta = (X^T X)^{-1} X^T y} θ∂J(θ)=XTX2θ2XTy+00=2XTXθ2XTyXTXθ=XTyθ=(XTX)1XTy

这就是线性回归的正态方程。

五、探索性数据分析

描述性统计

以下表格提供了数据集的描述性统计信息,包括每列的计数、均值、标准差、最小值和各个百分位数:

在这里插入图片描述

检查缺失值

使用热图可视化数据集中的缺失值:

plt.figure(figsize=(12,4))
sns.heatmap(df.isnull(), cbar=False, cmap='viridis', yticklabels=False)
plt.title('数据集中缺失值');

在这里插入图片描述

通过热图可以看出数据集中没有缺失值。

数据分布图

相关性热图

生成变量之间的相关性热图:

corr = df.corr()
sns.heatmap(corr, cmap='Wistia', annot=True)

在这里插入图片描述

热图显示变量之间没有显著相关性。

保险费用分布

生成保险费用的分布图和对数分布图:

f= plt.figure(figsize=(12,4))

ax=f.add_subplot(121)
sns.distplot(df['charges'],bins=50,color='r',ax=ax)
ax.set_title('Distribution of insurance charges')

ax=f.add_subplot(122)
sns.distplot(np.log10(df['charges']),bins=40,color='b',ax=ax)
ax.set_title('Distribution of insurance charges in $log$ sacle')
ax.set_xscale('log');

在这里插入图片描述

从左图可以看出,保险费用从 1120 到 63500 不等,图是右偏的。在右边的图中,我们将应用对数变换,数据分布大致趋于正常。为了进一步分析,我们将对目标变量保险费用应用对数变换。

保险费用与性别和吸烟情况的关系

生成性别与保险费用的提琴图以及吸烟情况与保险费用的提琴图:

f = plt.figure(figsize=(14,6))
ax = f.add_subplot(121)
sns.violinplot(x='sex', y='charges',data=df,palette='Wistia',ax=ax)
ax.set_title('Violin plot of Charges vs sex')

ax = f.add_subplot(122)
sns.violinplot(x='smoker', y='charges',data=df,palette='magma',ax=ax)
ax.set_title('Violin plot of Charges vs smoker');

在这里插入图片描述

从左图可以看出,男性和女性的保险费用大致相同,平均约为 5000 美元。从右图可以看出,吸烟者的保险费用明显高于非吸烟者,非吸烟者的平均费用约为 5000 美元,而吸烟者的最低保险费用就已经是 5000 美元。

保险费用与子女数量的关系

生成子女数量与保险费用的箱线图:

plt.figure(figsize=(14,6))
sns.boxplot(x='children', y='charges',hue='sex',data=df,palette='rainbow')
plt.title('Box plot of charges vs children');

在这里插入图片描述

生成子女数量与保险费用的聚合表:

df.groupby('children').agg(['mean','min','max'])['charges']

在这里插入图片描述

保险费用与地区和性别的关系

生成地区与保险费用的提琴图:

plt.figure(figsize=(14,6))
sns.violinplot(x='region', y='charges',hue='sex',data=df,palette='rainbow',split=True)
plt.title('Violin plot of charges vs children');

在这里插入图片描述

保险费用与年龄和BMI的关系

生成年龄与保险费用的散点图以及BMI与保险费用的散点图:

f = plt.figure(figsize=(14,6))
ax = f.add_subplot(121)
sns.scatterplot(x='age',y='charges',data=df,palette='magma',hue='smoker',ax=ax)
ax.set_title('Scatter plot of Charges vs age')

ax = f.add_subplot(122)
sns.scatterplot(x='bmi',y='charges',data=df,palette='viridis',hue='smoker')
ax.set_title('Scatter plot of Charges vs bmi')
plt.savefig('sc.png');

在这里插入图片描述

从左图可以看出,投保人的最低年龄为 18 岁。保单中有多个等级,大多数非吸烟者选择 1 st 1^{\text{st}} 1st 2 nd 2^{\text{nd}} 2nd 等级,吸烟者保单从 2 nd 2^{\text{nd}} 2nd 3 rd 3^{\text{rd}} 3rd 等级开始。

身体质量指数 (BMI) 是基于身高和体重的身体脂肪测量指标,适用于成年男性和女性。最低 BMI 为 16 kg/m 2 16 \text{kg/m}^2 16kg/m2,最高可达 54 kg/m 2 54 \text{kg/m}^2 54kg/m2

六、数据预处理

编码

机器学习算法不能直接处理分类数据,必须将分类数据转换为数字。编码的主要方法包括标签编码、独热编码和虚拟变量陷阱。

标签编码是将分类标签转换为数字形式,使算法能够理解和处理这些标签。

独热编码是将分类变量表示为二进制向量。这种方法首先将分类值映射到整数值(标签编码),然后每个整数值都表示为一个二进制向量,除该整数的索引位置为1外,其余位置均为0。

虚拟变量陷阱是指当一个变量可以由其他变量预测出来时的多重共线性问题。为了避免虚拟变量陷阱,可以通过删除一个变量和原始变量来处理。

使用 pandas 库的 get_dummies 函数,可以在一行代码中完成上述所有编码步骤。下面的代码演示了如何获取性别、子女、吸烟者和地区特征的虚拟变量,并通过设置 drop_first=True 来避免虚拟变量陷阱:

# 虚拟变量编码
categorical_columns = ['sex', 'children', 'smoker', 'region']
df_encode = pd.get_dummies(data=df, prefix='OHE', prefix_sep='_',
                           columns=categorical_columns, drop_first=True, dtype='int8')

# 验证虚拟变量过程
print('原始数据框的列:\n', df.columns.values)
print('\n数据集的行和列数:', df.shape)
print('\n编码后数据框的列:\n', df_encode.columns.values)
print('\n数据集的行和列数:', df_encode.shape)

Box-Cox 变换

Box-Cox 变换是一种将非正态因变量转换为正态形状的方法。正态性是许多统计技术的重要假设;如果数据不正态,可以应用 Box-Cox 变换以满足正态性假设。Box-Cox 变换的公式如下:
{ y λ − 1 λ , y i ¬ = 0 l o g ( y i ) λ = 0 \mathbf{ \begin {cases}\frac {y^\lambda - 1}{\lambda},& y_i\neg=0 \\ log(y_i) & \lambda = 0 \end{cases}} {λyλ1,log(yi)yi¬=0λ=0
Box-Cox 变换的关键在于找到合适的 λ \mathbf{\lambda} λ 值。以下代码展示了如何进行 Box-Cox 变换,并返回转换后的变量、 λ \mathbf{\lambda} λ 值和置信区间:

from scipy.stats import boxcox
y_bc, lam, ci = boxcox(df_encode['charges'], alpha=0.05)

# 输出置信区间和lambda值
ci, lam

Box-Cox 变换后,置信区间为 ( ( − 0.01140290617294196 , 0.0988096859767545 ) , λ = 0.043649053770664956 ) ((-0.01140290617294196, 0.0988096859767545), \mathbf{\lambda} = 0.043649053770664956) ((0.01140290617294196,0.0988096859767545),λ=0.043649053770664956)

由于 Box-Cox 变换并未在本模型中表现更好,因此采用对数变换:

# 对数变换
df_encode['charges'] = np.log(df_encode['charges'])

原始分类变量被删除,并且特定分类变量的独热编码变量列之一也被删除。因此,通过使用 get_dummies 函数,我们完成了所有三个编码步骤。

七、训练集和测试集划分

在数据预处理完成后,需要将数据划分为训练集和测试集。下面的代码展示了如何使用 scikit-learn 库中的 train_test_split 函数进行数据划分:

from sklearn.model_selection import train_test_split

X = df_encode.drop('charges', axis=1)  # 独立变量
y = df_encode['charges']  # 因变量

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=23)

八、模型构建

在此步骤中,使用我们的线性回归方程 θ = ( X T X ) − 1 X T y \mathbf{\theta = (X^T X)^{-1} X^Ty} θ=(XTX)1XTy 构建模型。第一步,我们需要向原始数据集添加一个特征 x 0 = 1 \mathbf{x_0 =1} x0=1

添加特征 x 0 = 1 x_0=1 x0=1

我们首先为训练集和测试集添加一个常数特征 x 0 = 1 x_0=1 x0=1

# Step 1: add x0 =1 to dataset
X_train_0 = np.c_[np.ones((X_train.shape[0],1)),X_train]
X_test_0 = np.c_[np.ones((X_test.shape[0],1)),X_test]

构建模型

接下来,我们使用正态方程来计算模型参数:

# Step 2: build model
theta = np.matmul(np.linalg.inv(np.matmul(X_train_0.T, X_train_0)), np.matmul(X_train_0.T, y_train))

将参数和特征列整理成一个数据框,以便查看:

# The parameters for linear regression model
parameter = ['theta_'+str(i) for i in range(X_train_0.shape[1])]
columns = ['intersect:x_0=1'] + list(X.columns.values)
parameter_df = pd.DataFrame({'Parameter': parameter, 'Columns': columns, 'theta': theta})

使用 Scikit-Learn 构建模型

为了验证正态方程求解的参数,我们还使用 Scikit-Learn 的线性回归模块来构建模型:

# Scikit Learn module
from sklearn.linear_model import LinearRegression
lin_reg = LinearRegression()
lin_reg.fit(X_train, y_train)  # Note: x_0 =1 is no need to add, sklearn will take care of it.

获取 Scikit-Learn 模型的参数并与正态方程的参数进行比较:

# Parameter
sk_theta = [lin_reg.intercept_] + list(lin_reg.coef_)
parameter_df = parameter_df.join(pd.Series(sk_theta, name='Sklearn_theta'))
parameter_df

结果如下:

外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传

从上述结果可以看出,两个模型获得的参数相同。因此,我们成功地使用正态方程建立了模型,并使用 Scikit-Learn 线性回归模块进行了验证。下一步是预测和模型评估。

九、模型评估

我们将使用模型参数对测试数据集进行预测,然后将预测值与测试集中的实际值进行比较。我们使用公式计算均方误差
J ( θ ) = 1 m ∑ i = 1 m ( y ^ i − y i ) 2 \mathbf{J(\theta) = \frac{1}{m} \sum_{i=1}^{m}(\hat{y}_i - y_i)^2} J(θ)=m1i=1m(y^iyi)2
R 2 \mathbf{R^2} R2 是数据与拟合回归线接近程度的统计度量。 R 2 \mathbf{R^2} R2 始终介于 0 到 100% 之间。0% 表示模型无法解释响应数据围绕其平均值的任何变化。100% 表示模型可以解释响应数据围绕平均值的所有变化。
R 2 = 1 − S S E S S T \mathbf{R^2 = 1 - \frac{SSE}{SST}} R2=1SSTSSE
SSE 是误差平方和:
S S E = ∑ i = 1 m ( y ^ i − y i ) 2 \mathbf{SSE = \sum_{i=1}^{m}(\hat{y}_i - y_i)^2} SSE=i=1m(y^iyi)2
SST 是总平方和:
S S T = ∑ i = 1 m ( y i − y ˉ i ) 2 \mathbf{SST = \sum_{i=1}^{m}(y_i - \bar{y}_i)^2} SST=i=1m(yiyˉi)2
其中, y ^ \mathbf{\hat{y}} y^ 是预测值, y ˉ \mathbf{\bar{y}} yˉ y \mathbf{y} y 的平均值。

使用正态方程进行预测

# Normal equation
y_pred_norm =  np.matmul(X_test_0, theta)

# Evaluation: MSE
J_mse = np.sum((y_pred_norm - y_test)**2) / X_test_0.shape[0]

# R_square 
sse = np.sum((y_pred_norm - y_test)**2)
sst = np.sum((y_test - y_test.mean())**2)
R_square = 1 - (sse / sst)

print('The Mean Square Error(MSE) or J(theta) is: ', J_mse)
print('R square obtained for normal equation method is :', R_square)

结果:

The Mean Square Error(MSE) or J(theta) is:  0.1872962232298182
R square obtained for normal equation method is : 0.7795687545055328

使用 Scikit-Learn 进行预测

# sklearn regression module
y_pred_sk = lin_reg.predict(X_test)

# Evaluation: MSE
from sklearn.metrics import mean_squared_error
J_mse_sk = mean_squared_error(y_pred_sk, y_test)

# R_square
R_square_sk = lin_reg.score(X_test, y_test)

print('The Mean Square Error(MSE) or J(theta) is: ', J_mse_sk)
print('R square obtained for scikit learn library is :', R_square_sk)

结果:

The Mean Square Error(MSE) or J(theta) is:  0.18729622322981887
R square obtained for scikit learn library is : 0.7795687545055319

该模型返回的 R 2 R^2 R2 值为 77.95%,因此它非常适合我们的数据测试,但我们仍然可以通过不同的技术提高其性能。请注意,我们通过应用自然对数变换来处理输出变量。当我们将模型投入生产时,对方程应用反对数变换。

十、模型验证

为了验证模型,我们需要检查线性回归模型的一些假设。线性回归模型的常见假设如下:

  1. 线性关系

在线性回归中,因变量和自变量之间的关系是线性的。这可以通过实际值与预测值的散点图来检查。

  1. 残差正态分布

残差误差图应呈正态分布。

  1. 残差均值

残差误差的平均值应尽可能为 0 或接近 0。

  1. 多元正态性

线性回归要求所有变量都是多元正态的。这个假设可以用 Q-Q 图来检查。

  1. 多重共线性

线性回归假设数据中几乎没有或没有多重共线性。当独立变量彼此高度相关时,就会出现多重共线性。方差膨胀因子 (VIF) 确定独立变量之间的相关性以及该相关性的强度:
V I F = 1 1 − R 2 \mathbf{VIF = \frac {1}{1-R^2}} VIF=1R21
如果 VIF > 1 且 VIF < 5 为中等相关性,VIF > 5 为多重共线性的临界水平。

  1. 同方差性

数据是同方差的,这意味着残差在回归线上相等。我们可以通过残差与拟合值的散点图来检查。如果存在异方差,图将呈现漏斗形状。

具体验证步骤

检查线性关系

# Check for Linearity
f = plt.figure(figsize=(14,5))
ax = f.add_subplot(121)
sns.scatterplot(y_test,y_pred_sk,ax=ax,color='r')
ax.set_title('Check for Linearity:\n Actual Vs Predicted value')

检查残差正态性和均值

# Check for Residual normality & mean
ax = f.add_subplot(122)
sns.distplot((y_test - y_pred_sk),ax=ax,color='b')
ax.axvline((y_test - y_pred_sk).mean(),color='k',linestyle='--')
ax.set_title('Check for Residual normality & mean: \n Residual eror');

在这里插入图片描述

检查多元正态性

# Check for Multivariate Normality
# Quantile-Quantile plot 
f,ax = plt.subplots(1,2,figsize=(14,6))
import scipy as sp
_,(_,_,r)= sp.stats.probplot((y_test - y_pred_sk),fit=True,plot=ax[0])
ax[0].set_title('Check for Multivariate Normality: \nQ-Q Plot')

检查同方差性

#Check for Homoscedasticity
sns.scatterplot(y = (y_test - y_pred_sk), x= y_pred_sk, ax = ax[1],color='r') 
ax[1].set_title('Check for Homoscedasticity: \nResidual Vs Predicted');

在这里插入图片描述

检查多重共线性

# 计算方差膨胀因子
VIF = 1 / (1 - R_square_sk)
VIF

结果为 4.536561945911138,表明不存在多重共线性。

结果分析

在我们的模型中,实际值与预测值的图是曲线,因此线性假设失败。残差均值为零,但残差误差图右偏。Q-Q 图显示对数值大于 1.5 的趋势增加。该图表现出异方差,误差在某些点之后加剧。方差膨胀因子值小于 5,因此不存在多重共线性。


参考: Linear Regression Tutorial
推荐:

  • python 错误记录
  • python 笔记
  • 数据结构

在这里插入图片描述

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:/a/781944.html

如若内容造成侵权/违法违规/事实不符,请联系我们进行投诉反馈qq邮箱809451989@qq.com,一经查实,立即删除!

相关文章

Halcon 铣刀刀口破损缺陷检测

一 OTSU OTSU&#xff0c;是一种自适应阈值确定的方法,又叫大津法&#xff0c;简称OTSU&#xff0c;是一种基于全局的二值化算法,它是根据图像的灰度特性,将图像分为前景和背景两个部分。当取最佳阈值时&#xff0c;两部分之间的差别应该是最大的&#xff0c;在OTSU算法中所采…

张量分解(2)——张量运算(内积、外积、直积、范数)

&#x1f345; 写在前面 &#x1f468;‍&#x1f393; 博主介绍&#xff1a;大家好&#xff0c;这里是hyk写算法了吗&#xff0c;一枚致力于学习算法和人工智能领域的小菜鸟。 &#x1f50e;个人主页&#xff1a;主页链接&#xff08;欢迎各位大佬光临指导&#xff09; ⭐️近…

Stream流真的很好,但答应我别用toMap()

你可能会想&#xff0c;toList 和 toSet 都这么便捷顺手了&#xff0c;当又怎么能少得了 toMap() 呢。 答应我&#xff0c;一定打消你的这个想法&#xff0c;否则这将成为你噩梦的开端。 让我们先准备一个用户实体类。 Data AllArgsConstructor public class User { priv…

【C#】函数方法、属性分文件编写

1.思想 分文件编写是面向对象编程的重要思想&#xff0c;没有实际项目作为支撑很难理解该思想的精髓&#xff0c;换言之&#xff0c;一两个函数代码量因为太少无法体现分文件编写减少大量重复代码的优势。 2.项目结构介绍 整项目的名称叫AutoMetadata&#xff0c;是一个基于W…

【第三版 系统集成项目管理工程师】第4章 信息系统架构

持续更新。。。。。。。。。。。。。。。 【第三版】系统集成项目管理工程师 考情分析4.1架构基础4.1.1指导思想&#xff08;非重点&#xff09; P1364.1.2设计原则&#xff08;非重点&#xff09; P1364.1.3建设目标&#xff08;非重点&#xff09; P1374.1.4总体框架 P138练习…

【web前端HTML+CSS+JS】--- CSS学习笔记02

一、CSS&#xff08;层叠样式表&#xff09;介绍 1.优势 2.定义解释 如果有多个选择器共同作用的话&#xff0c;只有优先级最高那层样式决定最终的效果 二、无语义化标签 div和span&#xff1a;只起到描述的作用&#xff0c;不带任何样式 三、标签选择器 1.标签/元素选择器…

什么牌子的头戴式蓝牙耳机好性价比高?

说起性价比高的头戴式蓝牙耳机,就不得不提倍思H1s,作为倍思最新推出的新款,在各项功能上都实现了不错的升级,二字开头的价格,配置却毫不含糊, 倍思H1s的音质表现堪称一流。它采用了40mm天然生物纤维振膜,这种振膜柔韧而有弹性,能够显著提升低音的量感。无论是深沉的低音还是清…

Android 10年,35岁,该往哪个方向发力

网上看到个网友发的帖子&#xff0c;觉的这个可能是很多开发人员都会面临和需要思考的问题。 不管怎样&#xff0c; 要对生活保持乐观&#xff0c;生活还是有很多的选择和出路的。 &#xff08;内容来自网络&#xff0c;不代表个人观点&#xff09; 《Android Camera开发入门》…

机器人动力学模型及其线性化阻抗控制模型

机器人动力学模型 机器人动力学模型描述了机器人的运动与所受力和力矩之间的关系。这个模型考虑了机器人的质量、惯性、关节摩擦、重力等多种因素&#xff0c;用于预测和解释机器人在给定输入下的动态行为。动力学模型是设计机器人控制器的基础&#xff0c;它可以帮助我们理解…

element-plus的文件上传组件el-upload

el-upload组件 支持多种风格&#xff0c;如文件列表&#xff0c;图片&#xff0c;图片卡片&#xff0c;支持多种事件&#xff0c;预览&#xff0c;删除&#xff0c;上传成功&#xff0c;上传中等钩子。 file-list&#xff1a;上传的文件集合&#xff0c;一定要用v-model:file-…

基于B/S模式和Java技术的生鲜交易系统

你好呀&#xff0c;我是计算机学姐码农小野&#xff01;如果有相关需求&#xff0c;可以私信联系我。 开发语言&#xff1a;Java 数据库&#xff1a;MySQL 技术&#xff1a;B/S模式、Java技术 工具&#xff1a;Visual Studio、MySQL数据库开发工具 系统展示 首页 用户注册…

RAG综述汇总

第一篇&#xff1a;Retrieval-Augmented Generation for Large Language Models: A Survey(同济/复旦) 论文链接 1.简介 这篇全面的综述论文详细研究了 RAG 范式的发展&#xff0c;包括 Naive RAG、Advanced RAG 和 Modular RAG。介绍了 RAG 框架的三个基础技术&#xff0c;…

Python28-7.4 独立成分分析ICA分离混合音频

独立成分分析&#xff08;Independent Component Analysis&#xff0c;ICA&#xff09;是一种统计与计算技术&#xff0c;主要用于信号分离&#xff0c;即从多种混合信号中提取出独立的信号源。ICA在处理盲源分离&#xff08;Blind Source Separation&#xff0c;BSS&#xff0…

CANopen协议开发梳理总结笔记教程

0、提醒 CANOpen使用时&#xff0c;需要清楚什么是大端和小端&#xff0c;这对于CANOpen数据发送及解析时&#xff0c;有很大的帮助。且学习开发CANOpen时&#xff0c;需要具备一定的CAN基础。 1、CANOpen协议介绍 ①、什么是CANOpen协议 CANOpen协议是一种架构在控制局域网络…

yaml格式转换成json格式

yaml格式转换成json格式 ①postman生成的结果是yaml格式 ps&#xff1a;postman输出的格式是没有自动换行的&#xff0c;需要将内容换行 ②复制到Python的脚本跑一趟&#xff1a;自动换行并去掉/n&#xff1b; str " "//(postman输出的内容&#xff09; print(st…

【python技巧】parser传入参数

参考网址: https://lightning.ai/docs/pytorch/LTS/api/pytorch_lightning.utilities.argparse.html#pytorch_lightning.utilities.argparse.add_argparse_args 1. 简单传入参数. parse_known_args()方法的作用就是把不在预设属性里的参数也返回,比如下面这个例子, 执行pytho…

2024年信息系统项目管理师1批次上午客观题参考答案及解析(1)

1、新型基础设施建设是以新发展理念为引领&#xff0c;以()为驱动&#xff0c;以信息网络为基础&#xff0c;面向高质量发展需要&#xff0c;提供数字转型、智能升级、融合创新等服务的基础设施体系。 A&#xff0e;技术创新 B&#xff0e;人工智能 C&#xff0e;区块链 D&…

代码随想录算法训练营第二十七天|452. 用最少数量的箭引爆气球、435. 无重叠区间、763.划分字母区间

452. 用最少数量的箭引爆气球 如何使用最少的弓箭呢&#xff1f; 直觉上来看&#xff0c;貌似只射重叠最多的气球&#xff0c;用的弓箭一定最少&#xff0c;那么有没有当前重叠了三个气球&#xff0c;我射两个&#xff0c;留下一个和后面的一起射这样弓箭用的更少的情况呢&am…

STM32-输入捕获IC和编码器接口

本内容基于江协科技STM32视频学习之后整理而得。 文章目录 1. 输入捕获IC1.1 输入捕获IC简介1.2 频率测量1.3 输入捕获通道1.4 主从触发模式1.5 输入捕获基本结构1.6 PWMI基本结构 2. 输入捕获库函数及代码2.1 输入捕获库函数2.2 6-6 输入捕获模式测频率2.2.1 硬件连接2.2.2 硬…

曹操的五色棋布阵 - 工厂方法模式

定场诗 “兵无常势&#xff0c;水无常形&#xff0c;能因敌变化而取胜者&#xff0c;谓之神。” 在三国的战场上&#xff0c;兵法如棋&#xff0c;布阵如画。曹操的五色棋布阵&#xff0c;不正是今日软件设计中工厂方法模式的绝妙写照吗&#xff1f;让我们从这个神奇的布阵之…