数据科学导论+学习过程
大概把这次作业 6 的实现过程记录一下吧，因为在过程中发现了很多数据分析很有用的一些知识，这里记录一下，以后可能会用到。
R 语言实战 - 基于 Logistic 回归模型的构建及验证 - 文章依据
R 语言 Logistic 回归～变量筛选
我该咋筛选变量进入多因素回归
logistics 回归建立预后模型过程中自变量太多，如何筛选？
上面是学习过程中学习的文章。
具体流程如下：
使用的 R 语言包如下：
| library(mice) library(autoReg) library(rms) library(caret) library(pROC) library(rmda) library(dplyr) library(rrtable) |
根据 PPT 提供的网址获取预处理过后的数据集 NatalRiskData.rData，进入 Rstudio 并使用 load 函数加载数据集。
由于数据已经经过预处理，我们通过 skimr::skim (sdata) 可以对数据进行分析，分析结果如下：
|  |
之后我们可以用 autoReg 包制作表格来展示各个分类变量的数据特征。
ft <- gaze(atRisk~.,data=sdata) %>%myft()
print(ft)
|  PWGT: 可能代表体重（Pregnant Weight），可能是在特定的医疗背景下测量的体重，比如妊娠期体重。 UPREVIS: 可能代表了某种类型的修订或更新次数（Update Revisions），或者是预定的访问次数（Prenatal Visits）。 CIG_REC: 可能指的是吸烟记录（Cigarette Record），表示个体是否有吸烟的习惯。 GESTREC3: 这个变量可能与妊娠期（Gestation Record）有关，特别是与妊娠周期的长度有关。 DPLURAL: 可能代表的是出生时胎儿的数量（Delivery Plurality），如单胞胎、双胞胎或三胞胎及以上。 ULD_MECO: 可能代表的是与胎儿的粪便（Meconium）有关的紧急劳动解除（Urgent Labor Decisions Meconium）。 ULD_PRECIP: 可能指的是急速分娩（Urgent Labor Decisions Precipitate），即分娩过程异常迅速的情况。 ULD_BREECH: 可能指的是臀位分娩（Urgent Labor Decisions Breech），即胎儿在分娩时臀部或脚部先行的情况。 URF_DIAB: 可能代表妊娠期糖尿病（Urgent Risk Factor Diabetes）的记录。 URF_CHYPER: 可能代表妊娠期慢性高血压（Urgent Risk Factor Chronic Hypertension）的记录。 URF_PHYPER: 可能代表妊娠期短暂高血压（Urgent Risk Factor Preeclampsia Hypertension）的记录。 URF_ECLAM: 可能代表子痫前期（Urgent Risk Factor Eclampsia），这是一种怀孕晚期的高血压疾病。 DBWT: 可能代表出生时的体重（Delivery Birth Weight）。 ORIGRANDGROUP: 这个变量不太明确，但可能代表某个原始的随机分组（Original Random Group）。 |
在建立测试集和训练集时，我们可以先对数据进变量筛选。这里我选择先单后多的变量筛选方式。如果某个变量单因素分析时 P < 0.05 ，就纳入多因素模型。同样，我们可以直接利用 autoReg 包进行实现。
我们先建立逻辑回归模型。
mod <- glm(atRisk~.,data = sdata,family = “binomial”)
summary(mod)
分析结果:
|  |
通过上图，我们发现 PWGT,UPREVIS,GESTREC3,ULD_MECO, ULD_BREECH, DBWT 这几个的 p 值都低于 0.05，故我们建立逻辑回归模型时可以将其纳入。
之后，我们需要把数据分为训练集和测试集，在训练集中训练模型并在测试集中进行检验。
set.seed(1234)
ind <- createDataPartition(sdata$atRisk,
p=0.8,
list = FALSE)
train_data <- sdata [ind,] # 训练集
test_data <- sdata [-ind,] # 测试集
然后通过基于 rms 包构建 Logistic 回归模型的构建，绘制列线图和校准曲线。
fit <- lrm(atRisk~PWGT+UPREVIS+GESTREC3+ULD_MECO+ULD_BREECH+DBWT,
data = train_data,x=TRUE,y=TRUE)
fit
|  |
通过结果可以知道每个自变量对结果的影响。
Intercept（截距）: -1.0618。这是当所有自变量值为 0 时，对数几率（log-odds）的估计值。
PWGT: 0.0042。这意味着 PWGT 的每个单位增加，与被分类为 atRisk 相关的对数几率增加 0.0042。
UPREVIS: -0.0318。这意味着 UPREVIS 的每个单位增加，与被分类为 atRisk 相关的对数几率减少 0.0318。
GESTREC3=< 37 weeks: 0.4968。这表明如果 GESTREC3 小于 37 周，则与被分类为 atRisk 相关的对数几率比 GESTREC3 大于等于 37 周的情况高 0.4968。
ULD_MECO: 1.0680。这表示 ULD_MECO 为真时，与被分类为 atRisk 相关的对数几率比 ULD_MECO 为假时高 1.0680。
ULD_BREECH: 0.3340。这表示 ULD_BREECH 为真时，与被分类为 atRisk 相关的对数几率比 ULD_BREECH 为假时高 0.3340。
DBWT: -0.0012。这意味着 DBWT 的每个单位增加，与被分类为 atRisk 相关的对数几率减少 0.0012。
通过图中可以知道伪 R - 平方的值为 0.144，，AIC 的值为 3353.608。
之后进行预测并绘制 ROC 曲线。
train_pred <- predict(fit,
newdata=train_data,
type = “fitted”)
test_pred <- predict(fit,
newdata=test_data,
type=“fitted”)
train_roc <- roc(train_data$atRisk,train_pred)
Setting levels: control = FALSE, case = TRUE
Setting direction: controls < cases
auc(train_roc)
Area under the curve: 0.7034
plot(train_roc)
test_roc <- roc(test_data$atRisk,test_pred)
Setting levels: control = FALSE, case = TRUE
Setting direction: controls < cases
auc(test_roc)
Area under the curve: 0.7135
plot(test_roc)
Fit 模型训练集 ROC 曲线
Fit 模型测试集 ROC 曲线
我们通过 ROC 曲线来获取最佳阈值。
coords_test <- coords(test_roc, “best”)
coords_test
threshold specificity sensitivity
1 0.02988667 0.9086334 0.4583333
并根据此建立混淆矩阵计算精度，准确率和召回率以及 F-score。
| > predicted_labels_test <- ifelse(test_pred >= best_threshold,TRUE,FALSE) > confusion_matrix <- table(predicted_labels_test, test_data$atRisk) > confusion_matrix predicted_labels_test FALSE TRUE FALSE 4694 52 TRUE 472 44 > accuracy <- (confusion_matrix[1, 1] + confusion_matrix[2, 2]) / sum(confusion_matrix) > recall <- confusion_matrix[2, 2] / (confusion_matrix[2, 2] + confusion_matrix[2, 1]) > cat(“Accuracy:”, accuracy, “\n”) Accuracy: 0.9004181 > cat(“Recall:”, recall, “\n”) Recall: 0.08527132 > precision <- confusion_matrix[2, 2] / (confusion_matrix[2, 2] + confusion_matrix[1, 2]) > f1_score <- 2 * (precision * recall) / (precision + recall) > cat(“Precision:”, precision, “\n”) Precision: 0.4583333 > cat(“F1 Score:”, f1_score, “\n”) F1 Score: 0.1437908 |
通过上述来看，筛选了一部分非显著性因素建立的模型精度较为高，但是召回率和准确率相对较低，这也导致了其 F-score 很低。这可能会导致很多漏判误判的情况。
下面作出双密度图：
test_data$pred <- test_pred
ggplot(data = test_data, aes(x = pred , fill = factor(atRisk))) +
-
geom_density(alpha = 0.5) +
-
labs(x = “Predicted Probabilities”, y = “Density”) +
-
ggtitle(“Dual Density Plot of Predicted Probabilities”)
Fit 模型预测双密度图
可以看出图中有以下特征。
“FALSE” 类的预测概率主要集中在 0 附近，这表明模型对于这类的预测非常有信心，认为它们几乎不可能是 “TRUE”。
“TRUE” 类的预测概率分布则相对更平坦，但大多数预测值仍然靠近 0。
在预测概率较低的区域（接近 0），两个分类的曲线有所重叠，这意味着在这些概率值下，模型在区分 “TRUE” 和 “FALSE” 上存在一定的不确定性。
这里我也对先前建立的概括了全部因素的模型 mod 进行了相同的操作，分别对其在测试数据集下的性能做出了评价。
Mod 模型测试集 ROC 曲线
对其取最佳阈值并建立混淆矩阵查看
| > coords_test_sum <- coords(roc_curve, “best”) > coords_test_sum threshold specificity sensitivity 1 0.01760564 0.7667441 0.59375 > best_threshold_sum <- 0.01760564 > predicted_labels_test_sum <- ifelse(sum_testpred >= best_threshold_sum , TRUE, FALSE) > confusion_matrix_sum <- table(predicted_labels_test_sum , test_data$atRisk) > confusion_matrix_sum predicted_labels_test_sum FALSE TRUE FALSE 3961 39 TRUE 1205 57 > accuracy_sum <- (confusion_matrix_sum [1, 1] + confusion_matrix_sum [2, 2]) / sum(confusion_matrix_sum ) > accuracy_sum [1] 0.763588 > recall_sum <- confusion_matrix_sum [2, 2] / (confusion_matrix_sum [2, 2] + confusion_matrix_sum [2, 1]) > recall_sum [1] 0.0451664 |
可以发现其召回率和准确率相对于 Fit 模型都低，这里可以得出提取主要因素对模型的性能还是有一定程度上的优化的，如果能够筛选出最相关的因素进行建模，模型的性能也会一定程度的提升。