心血管疾病(cardiovascular disease，简称CVD)是指心脏和血管疾病的总称，常见的心血管疾病包括：高血压(血压升高)、冠心病(心脏病发作)、脑血管疾病(中风)、周围血管疾病、心力衰竭、风湿性心脏病、先天性心脏病和心肌病等。据世界卫生组织报道，在2012年约有1750万人死于心血管疾病，占全球死亡总数的31％。这些死者中，估计740万人死于冠心病，670万人死于中风；由于心血管疾病患者不断增加，心血管疾病的诊断与治疗成为医疗行业的一大重要问题。
本案例基于Kaggle公开的心血管患者诊断数据，从患者生理指标、医疗检测指标和患者提供的主观信息出发，使用机器学习中的集成方法对患者是否患心血管疾病进行预测。首先我们对数据进行预处理和探索性分析；然后借助sklearn中的分类模型进行预测；最后比较不同分类器下的预测效果。
1.数据说明与预处理¶
1.1 数据说明
本案例数据基于患者的生理指标(性别、年龄、体重、身高等)、医疗检测指标(血压、血糖、胆固醇水平等)和患者提供的主观信息(吸烟、饮酒、运动等)共计12个特征，对患病情况进行分析。数据共计70000条，其中心血管疾病患者人数为34979，未患病人数为35021，各数据指标含义如下表所示：
| 特征来源 | 列名 | 含义说明 |
|---|---|---|
| 患者编号 | id | |
| 生理指标 | age | 患者年龄，单位天(day) |
| 生理指标 | height | 患者身高，单位cm |
| 生理指标 | weight | 患者体重，单位kg |
| 生理指标 | gender | 患者性别，1为女性，2为男性 |
| 医疗指标 | ap_hi | 收缩压(心脏收缩时的动脉血压最高值) |
| 医疗指标 | ap_lo | 舒张压(心脏舒张时的动脉血压最低值) |
| 医疗指标 | cholesterol | 胆固醇水平，1为正常, 2为超出正常水平, 3为大量超出正常水平 |
| 医疗指标 | gluc | 血糖浓度(血液中的葡萄糖含量)，1为正常, 2为超出正常水平, 3为大量超出正常水平 |
| 主观信息 | smoke | 患者是否经常吸烟，1代表经常吸烟，0代表不经常吸烟 |
| 主观信息 | alco | 患者是否经常饮酒，1代表经常饮酒，0代表不经常饮酒 |
| 主观信息 | active | 患者是否经常运动，1代表经常运动，0代表不经常运动 |
| 患病情况 | cardio | 当患有心血管疾病时，目标类别cardio为1，否则cardio为0 |
首先，我们载入需要使用的模块和读入数据。
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
data = pd.read_csv("./input/cardio_train.csv")
data.head()
data.info()
可以看到：所有特征字段均为数字类型，12个整数型和1个浮点小数型，同时每个特征字段均不存在缺失值。接下来我们使用描述性统计函数观察一下各特征的数值分布。
data.describe(include='all')
可以看到：
- 由于数据中的年龄
age用天为单位进行表示，不太方便我们的分析与观察，之后我们需要将年龄单位转换为年。 - 身高
height、体重weight数值分布可能存在异常，如最小身高为55cm，最低体重为10kg。 - 收缩压
ap_hi和舒张压ap_lo存在明显异常，最小值均为负数且最大值分别为16020和11000。
data['age'] = round(data['age']/365).astype(int)
fig = plt.figure(figsize = (6,4))
plt.hist(data['age'],bins=15)
plt.title("年龄分布")
plt.show()
可以看到，我们已经将数据中的年龄字段由天转换为年，该数据集中患者年龄的取值区间为30-65岁，其中50-60岁的人数最多。
接下来我们处理身高和体重的异常值，数据中患者均为成年人，根据生活常识剔除异常值,身高低于140cm且体重低于40kg则认为是异常值，按照这个方法从数据中剔除异常数据。
data = data[(data['height']>=140)&(data['weight']>=40)]
data.shape
可以看到，通过身高和体重的限定，我们剔除了198个异常数据，接下来我们处理血压的异常数据。
血压的单位为千帕，1千帕=7.6mmHg，成人正常的收缩压为90~130mmHg，正常的舒张压为60~90mmHg。根据世界卫生组织规定，成人收缩压小于90mmHg则可确诊为低血压;90到130之间称为正常血压；130到140之间称为临界高血压；大于等于140mmHg称为高血压。如下所示：
- 低血压：收缩压<90
- 正常血压：90<=收缩压<130
- 临界高血压：130<=收缩压<140
- 高血压：140<=收缩压
考虑到高低血压的变化范围，我们将收缩压ap_hi变化区间限定为[60,250]，将舒张压ap_lo的变化区间限定为[30,180]。
data = data[(60<=data['ap_hi'])&(data['ap_hi']<=250)&(30<=data['ap_lo'])&(data['ap_lo']<=180)]
data.shape
可以看到，通过血压变化的限定，我们又剔除了1215个异常数据。
2.探索性分析¶
由于数据特征较多，我们需要探究一下各特征之间的关系，这将有助于后面建立合理有效的模型。从以下几个方面进行入手。
- 采用直方图和核密度曲线观察各个特征的分布情况。
- 使用分组箱线图和分组条形图分析患病与各个特征之间的关系。
- 使用特征的相关系数图分析各个特征之间的相关性。
col = data.columns[[1,2,7,8]].tolist()
fig,ax = plt.subplots(2,2,figsize = (15,8))
for i in range(0,len(col)):
plt.subplot("22"+str(i+1))
plt.hist(data[col[i]])
plt.xlabel(col[i],fontsize = 15)
plt.ylabel("Count",fontsize = 13)
plt.tight_layout()
plt.show()
可以看到在整体数据集中，女性人数约为男性的两倍，大部分人的胆固醇水平和血糖浓度均处于正常水平。接下来我们使用核密度曲线观察连续性特征的数值分布情况。
col = data.columns[[3,4,5,6]].tolist()
fig,ax = plt.subplots(2,2,figsize = (15,8))
for i in range(0,len(col)):
plt.subplot('22'+str(i+1))
x = data[col[i]]
sns.distplot(x,color = 'lightpink')
plt.tight_layout()
plt.show()
可以看到身高大部分集中在150-180cm，体重大部分集中在50-100kg，收缩压和舒张压也在正常的变化范围内。接下来我们将年龄和性别进行分组，探索不同分组的患病情况。
fig,[ax1,ax2] = plt.subplots(1,2,figsize=(15,5))
sns.countplot(x='age', hue='cardio', data = data, palette="Set2",ax=ax1)
sns.countplot(x='gender', hue='cardio', data = data, palette="Set2",ax=ax2)
ax1.set_xlabel('年龄',fontsize = 13)
ax2.set_xlabel('性别',fontsize = 13)
ax1.set_title('年龄与患病情况的分布',fontsize = 13)
ax2.set_title('性别与患病情况的分布',fontsize = 13)
ax1.legend(['不患病','患病'],fontsize = 13)
ax2.legend(['不患病','患病'],fontsize = 13)
plt.show()
可以看到，随着年龄的增加，患病的比例逐渐上升；性别特征对患病的影响不大。
data['BMI'] = data['weight'] / np.square(data['height']/100)
print(data.head())
可以看到，我们已经计算出每个人的BMI指标，接下来我们使用分组箱线图查看一下性别、BMI与患病的取值情况。
fig = plt.figure(figsize=(4,8))
sns.catplot(x="gender", y="BMI", hue="cardio", data=data, color = "yellow",kind="box",height=10, aspect=.7)
plt.show()
箱线图是描述数据的一种方法，每一个"箱子"从上往下的五条线分别代表数据的五个统计量：最大非异常值、上四分位数、中位数、下四分位数和最小非异常值。 从整体来看，大部分人的BMI主要分布在20-35之间，女性(gender=1)比男性(gender=0)的BMI值略高，且患病群体(cardio=1，箱体为黄色)相较于未患病群体(cardio=0，箱体为白色)的BMI值更高。
最理想的体重指数是22，通过BMI判断体重的方法如下：
- 轻体重：BMI<18.5
- 健康体重：18.5<=BMI<24
- 超重：24<=BMI<28
- 肥胖：28<=BMI
我们新建一列obesity，根据BMI数值指标将每个人进行标记，轻体重标记为1，健康体重标记为2，超重标记为3，肥胖标记为4。
# 定义一个标记函数
def fill_obesity(BMI):
if BMI<18.5:return 1
elif BMI<24:return 2
elif BMI<28:return 3
else:return 4
data['obesity'] = data['BMI'].apply(fill_obesity)
print(data.head())
# 定义一个标记函数
def fill_pressure(ap_hi):
if ap_hi<90:return 1
elif ap_hi<130:return 2
elif ap_hi<140:return 3
else:return 4
data['pressure'] = data['ap_hi'].apply(fill_pressure)
print(data.head())
分析是否患病、BMI和血压的分布情况
# 计算交叉表
counts1 = pd.crosstab(data["obesity"],data["cardio"])
counts1.index = ["轻体重","健康体重","超重","肥胖"]
counts1.columns= ["不患病","患病"]
fig,[ax1,ax2] = plt.subplots(1,2,figsize=(12,4))
counts1.plot(kind ="bar",rot = 360,width = 0.5,color = ['lightskyblue','lightsalmon'],ax=ax1)
# 计算交叉表
counts2 = pd.crosstab(data["pressure"],data["cardio"])
counts2.index = ["低血压","正常血压","临界高血压","高血压"]
counts2.columns= ["不患病","患病"]
counts2.plot(kind ="bar",rot = 360,width = 0.5,color = ['lightskyblue','lightsalmon'],ax=ax2)
ax1.set_ylabel("人数",fontsize = 13)
ax1.set_xlabel("体重指数情况",fontsize = 13)
ax1.set_title("是否患病与体重指数的分布情况",fontsize = 13)
ax2.set_ylabel("人数",fontsize = 13)
ax2.set_xlabel("血压情况",fontsize = 13)
ax2.set_title("是否患病与血压的分布情况",fontsize = 13)
plt.show()
可以看到，肥胖和超重的人群中患病的占比更大，患病的人数也更多；临界高血压和高血压的人群患病比率也非常高。
import seaborn as sns
import matplotlib.pyplot as plt
import matplotlib
matplotlib.rcParams['axes.unicode_minus']=False
selected_features = data.drop(['id'],axis = 1)
corr = selected_features.corr()
fig = plt.figure(figsize = (16,10))
cmap = sns.diverging_palette(220, 10, as_cmap=True)
mask = np.zeros_like(corr, dtype=np.bool)
mask[np.triu_indices_from(mask)] = True
# 绘制热力图
sns.heatmap(corr, mask=mask, cmap=cmap,vmin=-0.5,vmax=0.7, center=0,annot = True,
square=True, linewidths=.5,cbar_kws={"shrink": .5});
plt.title("特征间的相关性",fontsize = 13)
plt.show()
颜色越深，特征之间的相关性越强。从特征之间的相关系数图可以看出，是否患病与血压情况、年龄、体重指数等特征相关性较强。
3.分类建模¶
判断患者是否患病是一个二分类问题，我们将使用逻辑回归、随机森林和GBDT这三种分类模型进行建模。首先，将数据集按照4：1划分成训练集和测试集，使用的方法为train_test_split()；之后在评价模型好坏时，我们分别使用函数classification_report()、confusion_matrix()和accuracy_score,用于输出模型的预测报告、混淆矩阵和分类正确率。
from sklearn.model_selection import train_test_split
from sklearn.metrics import classification_report
from sklearn.metrics import confusion_matrix
from sklearn.metrics import accuracy_score
##训练集测试集的划分
x = data.drop(['id',"cardio"],axis = 1)
y = data["cardio"]
x_train, x_test, y_train, y_test = train_test_split(x, y,test_size = 0.2, random_state = 2)
## 逻辑回归
from sklearn.linear_model import LogisticRegression
lr = LogisticRegression(penalty='l1',C=0.6,random_state = 8).fit(x_train,y_train)
y_pred = lr.predict(x_test)
print(classification_report(y_test,y_pred))
print(confusion_matrix(y_test, y_pred))
print("The acc of logistic regression is %.5f"%(accuracy_score(y_test,y_pred)))
## 随机森林
from sklearn.ensemble import RandomForestClassifier
rf = RandomForestClassifier(n_estimators = 100, max_depth = 20,max_features = 10, random_state = 20).fit(x_train,y_train)
y_pred = rf.predict(x_test)
print(classification_report(y_test,y_pred))
print(confusion_matrix(y_test, y_pred))
print("The acc of randomforest is %.5f"%(accuracy_score(y_test,y_pred)))
## 梯度提升树
from sklearn.ensemble import GradientBoostingClassifier
gbdt = GradientBoostingClassifier(n_estimators=100, learning_rate=0.1,max_depth=4, random_state=10).fit(x_train, y_train)
y_pred = gbdt.predict(x_test)
print(classification_report(y_test,y_pred))
print(confusion_matrix(y_test, y_pred))
print("The acc of gbdt is %.5f"%(accuracy_score(y_test,y_pred)))
从F1值和准确率对三个模型进行对比，逻辑回归、随机森林和GBDT的F1值分别为0.72、0.72和0.73，准确率分别为72.61%、72.10%和73.46%，GBDT的预测效果最好。最后，我们通过输出随机森林和GBDT模型的特征重要性对关键特征进行筛选。
figure,[ax1,ax2] = plt.subplots(1,2,figsize=(12,6))
rf_importance = rf.feature_importances_
index = data.drop(['id',"cardio"], axis=1).columns
rf_feature_importance = pd.DataFrame(rf_importance.T, index=index,columns=['score']).sort_values(by='score', ascending=True)
# 水平条形图绘制
rf_feature_importance.plot(kind='barh', title='随机森林特征重要性',legend=False,ax=ax1)
gbdt_importance = gbdt.feature_importances_
index = data.drop(['id',"cardio"], axis=1).columns
gbdt_feature_importance = pd.DataFrame(gbdt_importance.T, index=index,columns=['score']).sort_values(by='score', ascending=True)
# 水平条形图绘制
gbdt_feature_importance.plot(kind='barh', title='GBDT特征重要性',legend=False,ax=ax2)
plt.show()
从两种模型的特征重要性排名中可以看出，收缩压ap_hi、年龄age、体重指数BMI、体重weight特征的重要性较高。
4.总结¶
本案例基机器学习模型对心血管疾病诊断数据进行了分析和预测，通过数据预处理、探索性分析和分类建模三个阶段进行了深入挖掘。在数据预处理中，根据身高和体重的变化范围、血压的变化范围剔除了部分异常值；
在数据探索性分析中，新增了BMI指标、肥胖情况和血压情况三个特征，并通过分组对比了不同种类的人群中的患病占比；
在分类建模过程中，分别使用了不同方法进行预测，通过对比预测准确率和F1值对模型进行评估，结果发现GBDT的预测效果最好。另外，从模型的预测结果来看，收缩压ap_hi、年龄age、体重指数BMI、体重weight等特征对心血管疾病的诊断起着重要作用，因而在实际心血管疾病诊断时，可重点参考这些指标。