全文链接：http://tecdat.cn/?p=22596

本报告是对心脏钻研的机器学习/数据迷信考察剖析。更具体地说，咱们的指标是在心脏钻研的数据集上建设一些预测模型，并建设探索性和建模办法。但什么是心脏钻研？ （ 点击文末“浏览原文”获取残缺代码数据 ）

[]()钻研纲要

  • 介绍数据集和钻研的指标
  • 摸索数据集
  • 可视化
  • 应用Chi-Square独立测验、Cramer's V测验和GoodmanKruskal tau值对数据集进行摸索
  • 预测模型，Logisitic回归和RandomForest
  • step()
  • bestglm()
  • 两个逻辑回归的实例
  • 应用5折穿插验证对模型实例进行评估
  • 变量抉择改良
  • 随机森林模型
  • 用RandomForest和Logisitc回归进行预测
  • 应用可视化进行最终的模型摸索
  • 论断和下一步改良

相干视频

拓端

，赞25

拓端

，赞15

拓端

，赞20

[]()1.简介

咱们浏览了对于FHS的材料：

心脏钻研是对社区自在生存的人群中心血管疾病病因的长期前瞻性钻研。心脏钻研是流行病学的一个里程碑式的钻研，因为它是第一个对于心血管疾病的前瞻性钻研，并确定了危险因素的概念。

该数据集是FHS数据集的一个相当小的子集，有4240个观测值和16个变量。这些变量如下：

  1. 观测值的性别。该变量在数据集中是一个名为 "男性 "的二值。
  2. 年龄：体检时的年龄，单位为岁。
  3. 教育 : 参与者教育水平的分类变量，有不同的级别。一些高中（1），高中/GED（2），一些大学/职业学校（3），大学（4)
  4. 目前吸烟者。
  5. 每天抽的烟的数量
  6. 查看时应用抗高血压药物的状况
  7. 流行性中风。流行性中风（0 = 无病）。
  8. 流行性高血压（prevalentHyp）。流行性高血压。如果承受医治，受试者被定义为高血压
  9. 糖尿病。依据第一次查看的规范医治的糖尿病患者
  10. 总胆固醇(mg/dL)
  11. 收缩压(mmHg)
  12. 舒张压(mmHg)
  13. BMI: 身材品质指数，体重（公斤）/身高（米）^2
  14. 心率（次/分钟）
  15. 葡萄糖。血糖程度(mg/dL)

最初是因变量：冠心病（CHD）的10年危险。

这4240条记录中有3658条是残缺的病例，其余的有一些缺失值。

[]()2.理解数据的意义

在每一步之前，要加载所需的库。

require(knitr)require(dplyr)require(ggplot2)require(readr)require(gridExtra) #出现多幅图

而后，加载心脏钻研的数据集。

[]()2.1 变量和数据集构造的查看

咱们对数据集进行一次查看。

dim(dataset)
kable(head(dataset))

str(dataset)

##查看变量的摘要summary(dataset)

[]()2.2 数据集的单变量图

生成一个数据集的所有单变量图。

# 须要删除字符、工夫和日期等变量                                                     geom_bar(data = dataset,                                                     theme_linedraw()+#colnames(dataset)        marrangeGrob(grobs=all_plots, nrow=2, ncol=2)

转存失败从新上传勾销

这是为了取得对变量，对整个问题和数据集的了解，将通过多变量或至多双变量的可视化来实现。

点击题目查阅往期内容

数据分享|R语言逻辑回归、Naive Bayes贝叶斯、决策树、随机森林算法预测心脏病

转存失败从新上传勾销

左右滑动查看更多

转存失败从新上传勾销

01

02

03

转存失败从新上传勾销

04

[]()2.3 数据集的双变量图：因变量和预测因素之间的关系

当初咱们能够进行一些双变量的可视化，特地是为了看到因变量（TenYearCHD）和预测因素之间的关系。因为图的数量太多，不是所有的一对变量都能被考察到！咱们能够在前面的步骤中持续考察。咱们能够稍后再回到这一步，深刻理解。

上面的代码能够生成因变量的所有双变量图。因为因变量是一个二元变量，所以当预测变量是定量的时候，咱们会有boxplots，或者当预测变量是定性的时候，咱们会有分段的bar图。

for (var in colnames(dataset) ){        if (class(dataset[,var]) %in% c("factor","logical") ) {ggplot(data = dataset)  +                         geom_bar( aes_string(x = var,        } else if (class(dataset[,var]) %in% c("numeric","double","integer") ) { ggplot(data = dataset) +                         geom_boxplot()

依据咱们把握的状况，男性与TenYearCHD间接相干，因而男性这个变量仿佛是一个绝对较好的预测因素。同样，年龄仿佛也是一个很好的预测因素，因为TenYearCHD == TRUE的病人有较高的年龄中位数，其散布简直类似。相同，不同类别的教育和因变量之间仿佛没有关系。目前的吸烟者变量与因变量有轻微的关系，因为目前的吸烟者患TenYearCHD的危险略高。

[]()2.4 应用Goodman&Kruskal tau测验定性变量之间的关系

然而，除了这些实质上是定性办法的图表外，人们可能心愿对这种关联有一个数字值。为了有这样的数字测量，我想应用Goodman&Kruskal的tau测量，这是两个无序因子，即两个分类/名义变量之间的关联测量。在咱们这个数据集中的因子变量中，只有教育是_序数变量_，即它的类别有意义。这种测量方法比Cramer's V或chi-square测量方法更具信息量。

GKtauData(cat_variables)plot(dataset)

能够看出，对于因变量的变异性，预测因素的解释力十分小。换句话说，依据Goodman和Kruskal's tau度量，咱们的预测因素和因变量之间简直没有关联。这能够从TenYearCHD一栏的数值中看出。

假如我的G&Ktau测验正确的话，这对模型来说并不是一个好消息。

为了测验这些发现，咱们能够用Chi-square测验来测验分类变量与因变量的关联的显著性，而后用Phi相关系数来评估可能的关联的强度。Phi用于2x2等值表。对于更大的表格，即有更多层次的变量，能够利用Cramer's V。

chisq.test(table(dataset_cat$p.value ))

        phi(matrix(table(dataset_cat_variables[,7],

奇怪的是，当Chi-square的P值如此之低时，可能的关联的显著性为零。这两个测试（Chi-square和Phi相干）在大量的察看中基本上得出雷同的后果，因为一个是基于正态分布的，另一个是基于t散布的。

[]()2.5 多重共线性的双变量剖析

该模型的真正问题在于共线性现象。共线性关系产生在两个预测因子高度相干的状况下。咱们须要查看这种个性，而后持续建设对数回归模型。

依据Goodman和Kruskal's tau图，咱们不应该放心共线性。然而，有序变量的教育变量呢？Cramer's V测验显示，其强度不大。

#  教育与其余分类变量的Chi square独立性测试  chisq.test(table(education,variables[,x]))$p.value )
#将教育变量从新定位到数据集的第一个变量上         assocstats(x = table(dataset_cat_variables[,1], dataset_$cramer ) )

没有一个变量显示与教育有很强的关联。Cramer's V的最高值是0.145，这在教育和性别之间是相当弱的。

然而诸如currentSmoker和cigsPerDay这样的变量呢？很显著，其中一个是能够预测的。有一个数字变量和一个分类变量，咱们能够把数字变量分成几个类别，而后应用Goodman和Kruskal's tau。GroupNumeric()函数能够帮忙将定量变量转换成定性变量，然而，基于对数据的主观了解，以及之前看到的cigsPerDay的多模态散布，在这里应用cut()函数很容易。
当初让咱们检查一下GKtau的数值

class_list <- lapply(X = 1:ncol(dataset_2), function(x) class(dataset_2[,x]))t <- sapply(X = names(class_list) , FUN = function(x) TRUE %in% ( class_list[x] %in% c("factor","logical")) )dataset_cat_variables_2 <- subset(x = dataset_2, select = t )plot(dataset_2)

从矩阵图上的tau值及其背景形态，咱们能够看到cigsPerDay能够齐全解释currentSmoker的变异性。这并不奇怪，因为如果咱们晓得一个人每天抽多少支烟就能够断言咱们晓得一个人是否是吸烟者!

第二个关联是cigsPerDay与男性的关系，但它并不强烈。因而，前者能够解释后者的较小的变动性。

在下一个数据集中，我把所有定量变量转换成定性/分类变量。当初咱们能够有一个全面的矩阵，只管因为转换，一些信息会失落。

dataset_3$totChol <- GroupNumeric(x = dataset$totChol , n = 5 )

咱们能够看到，sysBP和diaBP能够预测prevalentHyp，但不是很强。(0.5左右)。因而咱们能够在模型中保留prevalentHyp。第二点是对于GK tau的输入。

[]()3.预测模型:Logistic回归和RandomForest

当初是评估模型实例的时候了。在这里，咱们把逻辑回归称为模型。

咱们有两个实例。

  1. 一个包含所有原始变量的模型实例，特地是cigsPerday和currentSmoker变量
  2. 一个包含所有原始变量的模型实例，除了currentSmoker，cigsPerday被转换为一个因子变量
    为了评估模型实例，咱们能够应用数学调整训练误差率的办法，如AIC。另一种办法是应用验证数据集，依据模型在这个数据集上的体现来评估模型。在后一种办法中，我抉择应用K-fold Cross-Validation(CV)技术，更具体地说是5-fold CV。在这里，还有其余一些技术，如留一法穿插验证。

[]()3.1 两个Logistic回归模型实例

# 因为下一步的cv.glm()不能解决缺失值。# 我只保留模型中的残缺案例。dataset_1 <- dataset[complete.cases(dataset),] glm(TenYearCHD ~ . , family = "binomial")

这个模型是基于原始数据集的。有缺失值的记录被从数据集中省略，模型显示变量男性、年龄、cigsPerDay、totChol、sysBP和葡萄糖是显著的，而prevalentHyp在某种程度上是显著的。

 glm(formula = TenYearCHD ~ . , family = "binomial")

在第二个模型实例中，重要变量与前一个模型实例雷同。

一个十分重要的问题是，如何掂量这两个模型实例的性能以及如何比拟它们？有各种办法来掂量性能，但我在这里抉择了5折穿插验证法。

为了进行穿插验证和评估模型实例，咱们须要一个老本函数。boot软件包举荐的一个函数，是一个简略的函数，它能够依据一个阈值返回谬误分类的平均数。阈值默认设置为0.5，这意味着任何察看到的超过50%的CHD机会都被标记为有继续疾病的TRUE病例。从医学的角度来看，我把阈值升高到0.4，这样即便是有40%机会得心脏病的病例，也会被标记为承受进一步的医疗关注。升高阈值，减少了假阳性率，从而减少了医疗费用，但缩小了假阴性率，解救了生命。咱们能够应用敏感度或特异性作为老本函数。此外，也能够应用cvAUC软件包将曲线下面积（AUC）与CV联合起来。

[]()3.2 模型实例的穿插验证评估

model1_cv_delta <- cv.glm(  model1, cost = cost, K = 5)$delta[1]kable(data.frame("model1" = model1_cv_delta ,

kable(                caption = "CV-Accuracy", digits = 4)

咱们能够看到，两个模型十分类似，然而，模型2显示出轻微的劣势。准确率的确相当高。然而，让咱们看看咱们是否能够通过删除一些变量来改良model1。

[]()3.3 通过变量抉择改良模型

咱们看一下model1的总结。

summary(model1)

到当初为止，咱们始终假如所有的变量都必须蕴含在模型中，除非是共线性的状况。当初，咱们被容许通过删除不重要的变量。这里有几种办法，如前向抉择和后向抉择。

例如，后向抉择法是基于不显著变量的P值。淘汰持续进行，直到AIC显示没有进一步改善。还有stats::step()和bestglm::bestglm()函数来主动进行变量抉择过程。后者的软件包及其次要函数有许多抉择信息规范的选项，如AIC、BIC、LOOCV和CV，而前者的逐渐算法是基于AIC的。

bestglm(Xy = dataset_1 , family = binomial , IC = "BIC")step(object = model1 )

当初让咱们来看看这两个模型和它们的穿插验证误差。

bestglm_bic_model

基于BIC的bestglm::bestglm()将模型变量缩小到5个：男性、年龄、cigsPerDay、sysBP和葡萄糖。所有的变量都是十分显著的，正如预期的那样。

summary(step_aic_model)

基于AIC的step()函数将模型变量缩小到8个：男性、年龄、cigsPerDay,prevalentStroke、prevalentHyp、totChol、sysBP和glucose。值得注意的是，通过step()找到的最佳模型实例具备不显著的变量。

glm_cv_error <- cv.glm(                                 glmfit = glm(formula                                               family = binomial, data = dataset_1),step_cv_error <- cv.glm(glmfit = step_aic_model, cost = cost, K = 5)$delta[1]kable(bestglm_model_cv_error ,                step_model_cv_error )        )

穿插验证误分类误差

转存失败从新上传勾销

kable(data.frame("bestglm() bic model"                  "step() aic model"

穿插验证-准确度

AIC办法和BIC办法都能产生雷同的准确性。该抉择哪种办法呢？我宁愿抉择AIC，因为该模型实例有更多的预测因素，因而更有洞察力。然而，抉择BIC模型实例也是正当的，因为它更扼要。与model1的准确度相比，咱们通过变量抉择在准确度上有0.8475-0.842=0.00550.8475-0.842=0.0055的进步。然而，咱们失去了对于其余预测因子和因变量关系的信息。

[]()3.4 RandomForest模型

到目前为止，我只做了逻辑回归模型。有更多的模型能够用来为以后的问题建模，而RandomForest是一个受欢迎的模型。让咱们试一试，并将后果与之前的模型进行比拟。

#---- 差是每个RF模型实例的CV输入的谬误分类率#---- 每个选定的树的CV谬误分类率的最终后果被绘制进去 # 对于不同数量的树，咱们计算CV误差。for (n in seq(50,1000,50))             for (k in 1:5)                      rf_dataset_train <- dataset_1[fold_seq != k ,]                     rf_dataset_test <- dataset_1[fold_seq == k , ]                     rf_model <- randomForest( formula,kable(rf_df[sort(x = rf_df[,2])

#----- 误差基于RandomForest OOB，即RandomForest输入的混同矩阵for (n in seq(50,1000,50)) {        counter <- counter + 1                 rf_model <- randomForest( formula ntree = n, x =}ggplot() +         geom_point(data = rf_df , aes(x = ntree , y = accuracy)

在这里，我同时应用了CV和out-of-bag(OOB)来评估随机森林性能。

咱们能够看到，在50到1000棵树的范畴内，RandomForest模型的最高精度能够通过设置CV办法的树数等于400来取得。图中的红线显示了咱们从逻辑回归模型实例中失去的最佳CV精度。因为OOB的最高准确率高于CV的最高准确率，所以我抉择了CV的准确率，使其更加审慎。ntree=400的CVaccuracy=0.8486CVaccuracy=0.8486，比最好的逻辑回归模型差0.00020.0002! 然而，如果咱们思考OOB的准确性，那么RandomForest模型比最佳逻辑回归模型好0.00120.0012。

在RF中，模型的准确性有所提高，但代价是失去了可解释性。RF是一个黑箱，咱们无法解释预测因子和因变量之间的关系。

[]()3.5 模型对集体数据如何预测？

这里为了实现这个报告，我想在一个新的数据集上减少一个预测局部。该数据集只有一条记录，其中包含我本人的集体数据。换句话说，我曾经创立了一个模型，我想晓得它是否预测了我的CHD。

> pred_data$年龄 <- 31> pred_data$教育 <- factor(4, levels = c(1,2,3,4)) > pred_data$以后吸烟者 <- FALSE > pred_data$每日吸烟量 <- 0 > pred_data$抗高血压药物 <- FALSE> pred_data$流行性中风 <- FALSE > pred_data$流行性高血压 <- FALSE

逻辑回归模型的预测输入。

glm_BIC_opt <- glm(data = dataset_1 , formula ,family = binomial )predict(glm_BIC_opt, newdata = pred_data)

转存失败从新上传勾销

随机森林预测。

rf_model <- randomForest( formula = . ,predict(rf_model, pred_data)

因而，当初看来，我没有危险! 然而，正如我之前提到的，这些模型是为了教育和机器学习的实际，而不是为了医学预测！所以，我认为这些模型是有价值的。

[]()4.最终模型摸索

让咱们最初看一下这个模型

dataset_3 <- dataset_2[complete.cases(dataset_2),]dataset_3_GK <- plot(dataset_3_GK)

ggpplot(data = dataset , text.angle = 0,label.size =2 , order = 0  ) +        scale_colour_manual(values = color)+        scale_fill_manual(values = color)

左右滑动查看更多

01

02

03

04

后果大多合乎预期。依据GKtau值，预测因子之间的关联最小。这正是咱们想要的，以防止共线性现象。
然而，平行坐标依然显示了一些乏味的点。例如，年龄组与 "十年衰弱倒退 "后果之间的关联很有意思。较低的年龄组在TenYearCHD==TRUE中的参与度很低，这意味着年龄与该疾病有正相干。另一方面，与男性相比，女性（男性==FALSE）在0支烟和[1,20]支烟组的奉献更大。换句话说，男性偏向于抽更多的烟，而且是重度吸烟者。

桑吉图能够产生更好的洞察力，因为咱们能够沿着坐标轴察看样本。

[]()5.论断

在这项钻研中，为了建设预测模型，应用了包含4240个观测值和16个变量的心脏钻研的数据集。这些模型旨在预测十年后的冠心病（CHD）。
在对数据集进行摸索后，利用逻辑回归和随机森林模型来建设模型。应用K-Fold Cross-Validation对模型进行了评估。

为了扩大这项钻研，能够应用进一步的分类办法，如反对向量机（SVM）、梯度晋升（GB）、神经网络模型、K-近邻算法，甚至决策树。

点击文末 “浏览原文”

获取全文残缺代码材料。

本文选自《R语言随机森林RandomForest、逻辑回归Logisitc预测心脏病数据和可视化剖析》。

本文中的数据分享到会员群，扫描上面二维码即可加群！

点击题目查阅往期内容

数据分享|R语言逻辑回归、线性判别分析LDA、GAM、MARS、KNN、QDA、决策树、随机森林、SVM分类葡萄酒穿插验证ROC

MATLAB随机森林优化贝叶斯预测剖析汽车燃油经济性

R语言用Rcpp减速Metropolis-Hastings抽样预计贝叶斯逻辑回归模型的参数

R语言逻辑回归、Naive Bayes贝叶斯、决策树、随机森林算法预测心脏病

R语言中贝叶斯网络（BN）、动静贝叶斯网络、线性模型剖析错颌畸形数据

R语言中的block Gibbs吉布斯采样贝叶斯多元线性回归

Python贝叶斯回归剖析住房累赘能力数据集

R语言实现贝叶斯分位数回归、lasso和自适应lasso贝叶斯分位数回归剖析

Python用PyMC3实现贝叶斯线性回归模型

R语言用WinBUGS 软件对学术能力测验建设档次（分层）贝叶斯模型

R语言Gibbs抽样的贝叶斯简略线性回归仿真剖析

R语言和STAN,JAGS：用RSTAN,RJAG建设贝叶斯多元线性回归预测选举数据

R语言基于copula的贝叶斯分层混合模型的诊断准确性钻研

R语言贝叶斯线性回归和多元线性回归构建工资预测模型

R语言贝叶斯推断与MCMC：实现Metropolis-Hastings 采样算法示例

R语言stan进行基于贝叶斯推断的回归模型

R语言中RStan贝叶斯层次模型剖析示例

R语言应用Metropolis-Hastings采样算法自适应贝叶斯预计与可视化

R语言随机搜寻变量抉择SSVS预计贝叶斯向量自回归（BVAR）模型

WinBUGS对多元随机稳定率模型：贝叶斯预计与模型比拟

R语言实现MCMC中的Metropolis–Hastings算法与吉布斯采样

R语言贝叶斯推断与MCMC：实现Metropolis-Hastings 采样算法示例

R语言应用Metropolis-Hastings采样算法自适应贝叶斯预计与可视化

视频：R语言中的Stan概率编程MCMC采样的贝叶斯模型

R语言MCMC:Metropolis-Hastings采样用于回归的贝叶斯预计

R语言用lme4多层次（混合效应）狭义线性模型（GLM），逻辑回归剖析教育留级考察数据

R语言随机森林RandomForest、逻辑回归Logisitc预测心脏病数据和可视化剖析

R语言基于Bagging分类的逻辑回归(Logistic Regression)、决策树、森林剖析心脏病患者、决策树、森林剖析心脏病患者")

R语言用主成分PCA、 逻辑回归、决策树、随机森林剖析心脏病数据并高维可视化