R语言随机森林RandomForest、逻辑回归Logisitc预测心脏病数据和可视化分析（上）：https://developer.aliyun.com/article/1491743
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-近邻算法，甚至决策树。
