心脏病预测-----集成学习xgboost

    # 导入相关的宝
    import numpy as np 
    import pandas as pd
    import seaborn as sns
    from sklearn.neighbors import KNeighborsClassifier
    from sklearn.model_selection import GridSearchCV
    from sklearn.metrics import accuracy_score, confusion_matrix, classification_report
    from matplotlib import pyplot as plt
    from sklearn.linear_model import LogisticRegression
    from sklearn.metrics import roc_curve, roc_auc_score
    from sklearn.model_selection import train_test_split
    from sklearn.preprocessing import StandardScaler, RobustScaler, MinMaxScaler, LabelEncoder
    from sklearn.tree import DecisionTreeClassifier
    from sklearn.decomposition import PCA
    from sklearn.svm import SVC
    from sklearn.ensemble import RandomForestClassifier
    from xgboost import XGBClassifier
    import warnings
    
    # 导入数据
    warnings.simplefilter(action='ignore')
    df = pd.read_csv("C:/Users/David/Desktop/data/data/heart/heart.csv")
    df.head()
    
    # 基本的EDA和可视化
    df.describe()
    
    df.info()
    
    sns.set_context("talk")
    plt.figure(figsize=(15, 12))
    sns.heatmap(df.corr(), annot=True, fmt=".2f")

    df.groupby("target")[["thalach","chol","age","trestbps"]].mean()

    # 年龄分布直方图
    sns.set_style("darkgrid")
    sns.set_context("poster")
    plt.figure(figsize=(12, 9))
    sns.histplot(data=df["age"], alpha=0.6)
    plt.axvline(df["age"].mean(), color="red")
    plt.title("Distribution of age in the dataset")
    plt.text(33, 42, "Mean age: " + str(round(df["age"].mean(), 1)), color="indianred", weight="heavy", size="larger", bbox=dict(boxstyle="round", color="rosybrown"))
    plt.xlabel("Age")
    plt.ylabel("Count")
    plt.show()
    # 数据集中的年龄分布比较正常，平均年龄为54.5岁。
    
    # ST降压分布直方图
    sns.set_style("darkgrid")
    sns.set_context("poster")
    plt.figure(figsize=(12, 9))
    sns.histplot(data=df["oldpeak"], alpha=0.6)
    plt.title("Relative ST depression (Exercise-Rest, mm)")
    plt.xlabel("Difference (Exercise-Rest, mm)")
    plt.ylabel("Count")
    plt.show()

    #创建一个新的列
    df["Heart condition"] = df["target"].replace({1:"Present", 0:"Not present"})
    
    # 关于实际标绘的部分
    sns.set_style("darkgrid")
    sns.set_context("poster")
    plt.figure(figsize=(12, 9))
    sns.scatterplot("chol", "age", hue=df["Heart condition"], data=df)
    plt.title("Cases presenting heart disease compared by their age and blood cholesterol")
    plt.xlabel("Blood cholesterol (mg/dL)")
    plt.ylabel("Age")
    plt.show()
    
    sns.set_style("darkgrid")
    sns.set_context("poster")
    plt.figure(figsize=(12, 9))
    sns.scatterplot("thalach", "age", hue=df["Heart condition"], data=df)
    plt.title("Maximum heart rate achieved by presence of heart condition and distribution of age")
    plt.xlabel("Maximum heart rate (beats per minute)")
    plt.ylabel("Age")
    plt.show()

    sns.set_style("darkgrid")
    sns.set_context("poster")
    plt.figure(figsize=(12, 9))
    sns.scatterplot("trestbps", "age", hue=df["Heart condition"], data=df)
    plt.title("Change of rest heart bpm by age")
    plt.xlabel("Resting heart rate (beats per minute)")
    plt.ylabel("Age")
    plt.show()
    
    plt.figure(figsize=(12, 9))
    g = sns.countplot("exang", hue="slope", data=df)
    g.set_xticklabels(["Not present", "Present"])
    plt.title("Presence of exercise induced angina by the slope's direction on peak of the ST segment")
    plt.xlabel("Presence of exercise induced angina")
    plt.ylabel("Counts")
    plt.legend(["Upsloping", "Flat", "Downsloping"])
    plt.show()

    plt.figure(figsize=(12, 9))
    g = sns.countplot("cp", hue="sex", data=df)
    g.set_xticklabels(["Typical", "Atypical", "Non-Anginal", "Asymptomatic"])
    plt.title("Distribution of anginal pain type by gender")
    plt.xlabel("Anginal pain type")
    plt.ylabel("Count")
    plt.legend(["Female", "Male"])
    plt.show()
    
    plt.figure(figsize=(20, 9))
    g = sns.countplot(pd.cut(df["oldpeak"], 13), hue="Heart condition", data=df) #I had to bin the x values since there were too many distinct values to interprete
    plt.title("Difference ST segment depression (Exercise-Rest) by presence of heart disease")
    plt.xlabel("Difference")
    plt.ylabel("Count")
    plt.xticks(rotation=75)
    plt.legend(["Not present", "Present"])
    plt.show()

    #将连续的特征集合在一起，从而创建离散的分类列，可以帮助模型一般化数据和减少过拟合
    df["thalach"] = pd.cut(df["thalach"], 8, labels=range(1, 9))
    df["trestbps"] = pd.cut(df["trestbps"], 5, labels=range(8, 13))
    df["age"] = pd.cut(df["age"], 12, labels=range(12, 24))
    df["chol"] = pd.cut(df["chol"], 10, labels=range(24, 34))
    df["oldpeak"] = pd.cut(df["oldpeak"], 5, labels=range(34, 39))

    # 对类别变量独热编码
    a = pd.get_dummies(df, columns=["cp", "restecg", "slope", "thalach", "trestbps", "age", "chol", "thal", "oldpeak"], 
                   prefix=["cp", "restecg", "slope", "thalach", "trestbps", "age", "chol", "thal", "oldpeak"], drop_first=True)

    a = a.drop("Heart condition", axis=1)

    a.head()
    
    ## 确定特征重要性
    rf = RandomForestClassifier()
    rf.fit(a.drop("target", axis=1), a.target)
    importances = rf.feature_importances_
    features = pd.Series(importances, index=a.drop("target", axis=1).columns)
    plt.figure(figsize=(12, 16))
    features.plot(kind="barh")
    plt.show()

    # 构建模型和超参数调优
    X = a.drop(["target", "restecg_1", "thalach_2", "thalach_8", "trestbps_12", "age_13", "age_22", "age_23", 
            "chol_29", "chol_30", "chol_31", "chol_32", "chol_33", "oldpeak_37", "oldpeak_38"], axis=1)
    y = a.target
    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, stratify=y, random_state=42)
    
    #用xGBoosting运行GridSearchCV要花费很长时间，后面代码是添加了超参数的结果。
    # 查找最佳超参数
    # xgb = xgb = XGBClassifier(objective="binary:logistic", eval_metric="auc")
    # params = {"max_depth":np.arange(1, 8), "n_estimators":[100, 300, 500, 600, 900, 1000], 
    #           "learning_rate":[0.001, 0.01, 0.03, 0.05, 0.07], "colsample_bytree":[0.35, 0.5, 0.65],
    #           "subsample":[0.4, 0.5, 0.55, 0.7]
    #           }
    # grid_xgb = GridSearchCV(estimator = xgb, param_grid = params, cv=10, n_jobs=-1)
    # grid_xgb.fit(X_train, y_train)
    # xgb_pred = grid_xgb.predict(X_test)
    # print(classification_report(y_test, xgb_pred))
    
    xgb = XGBClassifier(objective="binary:logistic", eval_metric="auc", max_depth=1,
                    n_estimators=500,learning_rate=0.05,colsample_bytree=0.35,subsample=0.5)
    xgb.fit(X_train, y_train)
    xgb_pred = xgb.predict(X_test)

    # 模型评估
    print(grid_xgb.best_params_)
    
    print(grid_xgb.cv_results_["mean_test_score"].mean())
    
    print(classification_report(y_test, xgb_pred))

    sns.set_context("talk")
    plt.figure(figsize=(12, 9))
    sns.heatmap(confusion_matrix(y_test, xgb_pred), annot=True, xticklabels=["Healthy", "Sick"], yticklabels=["Healthy", "Sick"], fmt="g", cmap="icefire_r")
    plt.ylabel("Predicted")
    plt.xlabel("Actual")
    plt.title("Confusion matrix of xGBoosting")
    plt.show()

    xgb_prob = xgb.predict_proba(X_test)[:, 1]
    fpr, tpr, thresholds = roc_curve(y_test, xgb_prob)
    sns.set_style("darkgrid")
    sns.set_context("poster")
    plt.figure(figsize=(15, 10))
    plt.plot([0, 1], [0, 1], 'k--')
    sns.lineplot(fpr, tpr, alpha=0.6, ci=None)
    plt.xlabel('False positive rate')
    plt.ylabel('True positive rate')
    plt.title('ROC curve')
    plt.legend(["Baseline", "xGBoosting"])
    plt.show()
    print(roc_auc_score(y_test, xgb_prob))