逻辑回归模型的建立与求解——以泰坦尼克号上的幸存者为例

逻辑回归模型在数学建模中的应用

模型引出

• 逻辑回归（Logistic Regression）是一种用于分类问题的统计模型，特别适用于二分类问题。在数学建模中，逻辑回归被广泛应用于预测和分类问题，如预测疾病的发生、客户的购买行为、信用违约风险等。下面将详细介绍逻辑回归模型的建立与求解过程。

• 简而言之，逻辑回归与线性回归的主要区别在于：线性回归处理的是连续数值型数据，而逻辑回归适用于处理包含分类数据（如性别、等级等）的数据集，输出的是概率而非连续值。

模型建立

假设我们有一个响应变量$Y$和一组自变量$X_1,X_2,...,X_n$。逻辑回归模型的基本形式为：
$$P(Y=1|X)=\frac{1}{1+e^{-({\beta }_0+{\beta }_1{X}_1+{\beta }_2{X}_2+\cdots +{\beta }_n{X}_n)}}$$
其中：

• $P(Y=1|X)$是影响变量$Y$等于1得到概率，给定自变量$X$。
• ${\beta }_0$是截距项，${\beta }_1,{\beta }_2,...,{\beta }_n$是回归系数，
• $e$是自然对数的底数

为了简化计算，我们使用对数几率（Logit）函数，模型可以转化为：
$$\mathrm{Logit}(P)=\log \left(\frac{P}{1-P}\right)={\beta }_0+{\beta }_1{X}_1+{\beta }_2{X}_2+\cdots +{\beta }_n{X}_n$$
其中：

• $\mathrm{Logit}(P)$是$P$的对数几率。
• log是自然对数函数。

模型求解

逻辑回归模型的回归系数$\beta$是通过最大似然估计（MLE）来求解的。目标是找到使观察到的数据在模型下的似然函数最大的回归系数。

给定样本$(x_i,y_i)(i=1,2,...,m)$，似然函数定义为：
$$L(\beta )=\prod _{i=1}^{m}P(Y_i=1|X_i{)}^{y_i}\cdot {\left[1-P(Y_i=1|X_i)\right]}^{1-y_i}$$
其中$P(Y_i=1|X_i)$是第$i$个样本属于类别1的概率。

为了简化计算，我们使用对数似然函数（Log-Likelihood Function）：
$$\ell (\beta )=\sum _{i=1}^m\left[y_i\log (P(Y_i=1|X_i))+(1-y_i)\log \left[1-P(Y_i=1|X_i)\right]\right]$$
接下来，通过最大化对数似然函数 $\ell (\beta )$ 来估计回归系数 $\beta$。由于对数似然函数是关于回归系数 $\beta$ 的凸函数，可以使用梯度上升法、牛顿法等数值优化算法进行求解。最大化 $\ell (\beta )$ 的过程实际上是在寻找最优的回归系数，使得模型预测的概率与实际观察结果之间的差异最小。

结果解释

在模型求解之后，得到的回归系数 $\beta$ 可以解释为自变量对响应变量的影响大小。具体来说，每个回归系数 ${\beta }_j$ 表示的是在其他变量保持不变的情况下，自变量 $X_j$ 增加一个单位时，响应变量 $Y$ 为 1 的对数几率变化值。如果 ${\beta }_j$ 为正数，说明 $X_j$ 对 $Y$ 为 1 的概率有正向影响；如果 ${\beta }_j$ 为负数，则说明 $X_j$ 对 $Y$ 为 1 的概率有负向影响。

此外，还可以通过计算预测概率 $P(Y=1|X)$ 来进行分类预测。通常我们会设置一个阈值（例如 0.5），当 $P(Y=1|X)>0.5$ 时预测 $Y=1$，否则预测 $Y=0$。

通过这种方式，逻辑回归模型在许多分类问题中得到了广泛应用，如信用风险评估、疾病预测、市场营销中的客户分类等。

使用逻辑回归模型解决泰坦尼克号上的幸存者

数据集：链接：https://pan.baidu.com/s/1sP8CYKZBqmmi04uR8D8Qgg 提取码：1111

问题描述

在1912年4月15日，泰坦尼克号在其处女航中撞上冰山并沉没。船上2224名乘客和船员中，约1500人遇难。这一灾难引发了很多关于不同乘客群体生还率差异的研究。通过分析幸存者的数据，我们可以尝试建立一个模型来预测某个乘客在这场灾难中的生还概率。逻辑回归是一种常用的二分类模型，适用于这种问题。

数据准备

变量名	描述
PassengerId	乘客唯一识别ID
Survived	是否存活，0为否，1为是
Pclass	船舱等级，1、2、3等
Name	姓名
Sex	性别
Age	年龄
SibSp	和该乘客一起旅行的兄弟姐妹和配偶的数量
Parch	和该乘客一起旅行的父母和孩子的数量
Ticket	船票号
Fare	船票价格
Cabin	船舱号
Embarked	登船港口，S=英国南安普顿（Southampton, 起航点），C=法国瑟堡市（Cherbourg, 途经点），Q=爱尔兰昆士镇（Queenstown, 途经点）

这些变量中，Survived是我们要预测的目标变量，其余变量为模型的自变量。在数据准备过程中，我们需要处理数据中的缺失值，将分类变量转换为数值格式，并确保数据可以被模型正确解析。

数学模型的建立

假设我们有一个响应变量$Y$（乘客是否幸存）和一组自变量$X_1,X_2,...,X_n$（乘客的特征，如年龄、性别等）。逻辑回归模型的基本形式为：
$$P(Y=1|X)=\frac{1}{1+e^{-({\beta }_0+{\beta }_1{X}_1+{\beta }_2{X}_2+\cdots +{\beta }_n{X}_n)}}$$
其中：

• $P(Y=1|X)$ 是在给定自变量$X$的情况下，乘客生还的概率。
• ${\beta }_0$ 是截距项，${\beta }_1,{\beta }_2,...,{\beta }_n$ 是回归系数，表示每个自变量对生还概率的影响大小。
• $e$ 是自然对数的底数。

为了简化计算，我们可以将模型转换为对数几率（Logit）形式：
$$\mathrm{Logit}(P)=\log \left(\frac{P}{1-P}\right)={\beta }_0+{\beta }_1{X}_1+{\beta }_2{X}_2+\cdots +{\beta }_n{X}_n$$
在这个模型中，目标是通过最大似然估计（MLE）来求解回归系数$\beta$。

模型的求解

逻辑回归模型的回归系数$\beta$是通过最大似然估计（MLE）来求解的。最大似然估计的目标是找到使观察到的数据在模型下的似然函数最大的回归系数。

对于给定的样本 $(x_i,y_i)$ $(i=1,2,...,m)$，似然函数定义为：
$$L(\beta )=\prod _{i=1}^{m}P(Y_i=1|X_i{)}^{y_i}\cdot {\left[1-P(Y_i=1|X_i)\right]}^{1-y_i}$$
其中，$P(Y_i=1|X_i)$ 是第 $i$ 个样本属于类别 1 的概率。

为了简化计算，通常使用对数似然函数（Log-Likelihood Function）：
$$\ell (\beta )=\sum _{i=1}^m\left[y_i\log (P(Y_i=1|X_i))+(1-y_i)\log \left[1-P(Y_i=1|X_i)\right]\right]$$
通过最大化对数似然函数，我们可以获得回归系数 $\beta$ 的估计值。

使用Python语言进行模型的求解

数据读取

import pandas as pd # 导入pandas库
import numpy as np  # 导入numpy库


df_titanic = pd.read_csv('../dataset/train.csv') # 读取文件
print("训练集数据：", df_titanic.shape)
df_titanic.head() # 显示前5行数据

数据预处理

查看数据缺失值信息

# 缺失值查询
df_titanic.isnull().sum()

从上面结果可以得出，乘客的年龄、船舱号和登船港口存在缺失值，由于船舱号和登船港口对乘客的存活率没有影响，所以可以忽略不计，这里使用平均值来填充年龄特征。

df_titanic['Age'] = df_titanic['Age'].fillna(df_titanic['Age'].mean())  # 使用均值数据填充年龄特征
df_titanic.isnull().sum()   # 缺失值查询

数据整理

由于和该乘客一起旅行的兄弟姐妹和配偶的数量以及和该乘客一起旅行的父母和孩子的数量属于同一类，可能存在多重共线性问题，这里我们使用家庭成员（兄弟姐妹和配偶的数量+父母和孩子的数量）来表示。

df_titanic['FamilyNum']=df_titanic['SibSp']+df_titanic['Parch']

乘客的ID、姓名、船票价格、船舱号、登船港口等信息对乘客存活概率影响不大，所以可以不做考虑。

lr_titanic_train=df_titanic.copy()
lr_titanic_train=lr_titanic_train.drop(['PassengerId','Name','Ticket','Cabin','Embarked'],axis=1)
lr_titanic_train.head()

性别和船舱等级信息属于分类变量，$n$项分类信息生成$n-1$个哑元变量。

lr_titanic_train=pd.get_dummies(lr_titanic_train,drop_first=True,columns=['Pclass','Sex'],dtype=int)
lr_titanic_train=lr_titanic_train.drop(['SibSp','Parch'],axis=1)

• pd.get_dummies() 是 Pandas 库中的一个函数，用于将分类变量转换为哑变量（也叫虚拟变量或指示变量）。哑变量是用于分类数据的数值表示，一般是 0 或 1，表示某个类别的存在或不存在。

• drop_first=True 参数用于在生成哑变量时丢弃第一个类别的哑变量列，以避免虚拟变量陷阱（dummy variable trap）。虚拟变量陷阱是指在回归分析中，由于哑变量的完全共线性（即哑变量相互之间完全相关），导致模型无法准确估计参数。通过丢弃一个哑变量列，可以避免这种共线性问题。

整合后的数据如下图如所示：

模型的构建

以是否存活为因变量，其余信息为自变量，建立逻辑回归模型。

y=lr_titanic_train['Survived']
X=lr_titanic_train.drop(['Survived'],axis=1)

计算相关系数，检查是否有相关性较高的变量。

X.corr().abs()>0.8

使用statsmodels 库进行逻辑回归建模。

X=sm.add_constant(X)
model=sm.Logit(y,X).fit()
model.summary()

• sm.add_constant(X) 会在特征矩阵 X 中添加一列值为 1 的列，这一列对应于模型中的截距项 ${\beta }_0$。

• 使用逻辑回归（Logit）模型来拟合数据。sm.Logit 用于创建逻辑回归模型对象，然后通过 .fit() 方法对模型进行拟合。

• summary() 方法生成并打印模型的摘要信息。这包括模型的拟合结果、各个系数的估计值及其统计显著性等。

• 具体摘要信息如下表所示

名称	作用
Dep. Variable	目标变量，即模型中要预测的变量名称。
No. Observations	观测值的数量，即样本数。
Df Residuals	剩余自由度，等于样本数减去估计参数的数量。
Df Model	模型的自由度，等于模型中估计的参数数量。
Pseudo R-squared	类似于线性回归中的 R-squared，但适用于逻辑回归，表示模型解释数据的能力。
Log-Likelihood	对数似然值，用于衡量模型的拟合好坏。
coefficients (Coef.)	每个自变量的回归系数估计值，表示每个自变量对响应变量的影响。
Std Err	回归系数的标准误差，表示回归系数的估计值的精确度。
z	z 统计量，用于检验回归系数是否显著不为零。
P>|z|	p 值，表示回归系数显著性水平的小概率值，p 值越小，回归系数越显著。
[0.025 0.975]	置信区间（通常为 95% 置信水平），表示在这个范围内回归系数的真实值有 95% 的可能性。

通过观察**P>|z|**这列信息，可以发现Fare的p值较高，该自变量可能对目标变量的预测能力不强，或者与目标变量的关系不显著。所以舍弃这列数据，再进行回归。

X=X.drop(['Fare'],axis=1)
model=sm.Logit(y,X).fit()
model.summary()

根据计算出的截距可以得出下面信息，由于模型计算出来的是难以解释的对数几率，所以通过np.exp()方法把它转化为胜算比。

以上数据表明：

• 年龄每增加1岁，生还概率降低4%左右。
• 每多一名家庭成员，生还概率降低20%左右。
• 二等舱乘客的生还概率比一等舱低71%左右。
• 三等舱乘客的生还概率比一等舱低91%左右。
• 男性乘客的生还概率比女性低94%左右。

模型的测试

# 复制数据以便进行测试
lr_titanic_test = df_titanic.copy()
# 将分类变量转换为类别类型，并设置类别的顺序
lr_titanic_test['Pclass'] = pd.Categorical(lr_titanic_test['Pclass'], categories=['1', '2', '3'])
lr_titanic_test['Sex'] = pd.Categorical(lr_titanic_test['Sex'], categories=['female', 'male'])
lr_titanic_test['Embarked'] = pd.Categorical(lr_titanic_test['Embarked'], categories=['C', 'Q', 'S'])
# 使用独热编码将分类变量转换为虚拟变量，同时丢弃第一个类别以避免多重共线性
lr_titanic_test = pd.get_dummies(lr_titanic_test, drop_first=True, columns=['Pclass', 'Sex'], dtype=int)
# 计算家庭成员数量，家庭成员包括兄弟姐妹/配偶以及父母/子女
lr_titanic_test['FamilyNum'] = lr_titanic_test['SibSp'] + lr_titanic_test['Parch']
# 提取特征变量进行预测
X_test = lr_titanic_test[['Age', 'FamilyNum', 'Pclass_2', 'Pclass_3', 'Sex_male']]
# 为特征矩阵添加常数项（截距）
X_test = sm.add_constant(X_test)
# 使用逻辑回归模型进行预测，得到每个样本的生存概率
predicted_value = model.predict(X_test)
# 将预测的概率转换为布尔值（>0.5为True，表示生存；否则为False）
boolean_result = predicted_value > 0.5
# 将布尔值转换为0和1，其中True为1，False为0
converted_result = boolean_result.astype(int)
# 获取实际的Survived列的数据，用于比较
survived_data = lr_titanic_test['Survived']
# 计算预测结果与实际结果相同的百分比
matching_percentage = (converted_result == survived_data).mean() * 100
# 打印相同的百分比
print(f"相同的百分比：{matching_percentage:.2f}%")

经过计算，预测解结果与正确结果相同的百分比为79.46%，证明模型的效果良好。

预测数据

# 读入预测数据
df_titanic_test = pd.read_csv('../dataset/test.csv') # 读取文件
print("预测集数据：", df_titanic_test.shape)
df_titanic_test.head() # 显示前5行数据
# 缺失值查询
df_titanic_test.isnull().sum()
# 用均值填充Age和Fare
df_titanic_test['Age']=df_titanic_test['Age'].fillna(df_titanic_test['Age'].mean())
df_titanic_test['Fare']=df_titanic_test['Fare'].fillna(df_titanic_test['Fare'].mean())
df_titanic_test.isnull().sum()
# 将分类变量转换为类别类型，并设置类别的顺序
df_titanic_test['Pclass']=pd.Categorical(df_titanic_test['Pclass'],categories=['1','2','3'])
df_titanic_test['Sex']=pd.Categorical(df_titanic_test['Sex'],categories=['female','male'])
df_titanic_test['Embarked']=pd.Categorical(df_titanic_test['Embarked'],categories=['C','Q','S'])
# 使用独热编码将分类变量转换为虚拟变量，同时丢弃第一个类别以避免多重共线性
df_titanic_test=pd.get_dummies(df_titanic_test,drop_first=True,columns=['Pclass','Sex'],dtype=int)
# 计算家庭成员数量，家庭成员包括兄弟姐妹/配偶以及父母/子女
df_titanic_test['FamilyNum']=df_titanic_test['SibSp']+df_titanic_test['Parch']
X_test=df_titanic_test[['Age','FamilyNum','Pclass_2','Pclass_3','Sex_male']]X_test=df_titanic_test[['Age','FamilyNum','Pclass_2','Pclass_3','Sex_male']]
X_test=sm.add_constant(X_test)
predicted_value=model.predict(X_test)
predicted_value>0.5
(predicted_value>0.5).sum()

最终得出结论，测试集中有316名幸存者。