这是一道 2019 年的某公司笔试题，面试官要求用不少于 5 种统计模型来建模，在三天内完成并给一个展示。不过，本文仅给出探索分析和广义线性模型的结果，算是前一篇广义线性模型和指数族的应用。
这个数据的背景来自某医院，问题是要分析影响病人入院的因素。这里的响应变量是病人在多久后办理入院手续，医院的床位是非常重要的资源，此问题有实际意义。
- 借助图形来探索数据，并给出分析。
- 使用不同的统计模型来建模，比较模型效果。
- 给出最终的模型结果及解释。
这个数据集的说明很少，各个列除了列名能看出些信息，相关介绍都没有，只能根据脱敏后的数据探索，反推。
1 探索性数据分析
首先导入数据
hospital_waiting_time <- readRDS(file = "data/hospital_waiting_time.rds")
# 查看数据
str(hospital_waiting_time)## 'data.frame': 2625 obs. of 11 variables:
## $ 等待时间 : num 1 1.2 20 6 8.9 2.9 7.9 2.8 2.7 5 ...
## $ 门诊次 : int 2 7 43 1 3 1 10 3 6 2 ...
## $ 住院次 : int 1 1 1 1 1 1 1 1 1 1 ...
## $ 开住院条日期: int 3 3 3 3 3 3 3 3 3 3 ...
## $ 性别 : int 0 0 1 1 1 1 0 1 1 1 ...
## $ 年龄 : int 42 32 59 9 45 73 50 25 14 20 ...
## $ 入院疾病分类: int 3 1 1 3 3 3 4 1 2 3 ...
## $ 入院目的 : int 1 1 1 1 1 1 1 1 1 1 ...
## $ 住院类别 : int 2 2 2 2 2 2 2 2 2 2 ...
## $ 入院病情 : int 1 1 1 1 1 1 1 1 1 1 ...
## $ 医生 : int 2 2 2 2 2 4 2 2 4 4 ...一共 11 个变量 2625 条记录，各个变量已经脱敏，而根据变量名，可以初步推断出变量的类别和信息。
# 变量类型转化
hospital_waiting_time <- transform(
hospital_waiting_time,
性别 = as.factor(性别),
入院疾病分类 = as.factor(入院疾病分类),
入院目的 = as.factor(入院目的),
住院类别 = as.factor(住院类别),
入院病情 = as.factor(入院病情),
医生 = as.factor(医生)
)
# 汇总数据信息
summary(hospital_waiting_time)## 等待时间 门诊次 住院次 开住院条日期 性别
## Min. : 0.0 Min. : 0.0 Min. : 1.00 Min. :1.00 0:1103
## 1st Qu.: 3.2 1st Qu.: 2.0 1st Qu.: 1.00 1st Qu.:1.00 1:1522
## Median : 9.6 Median : 4.0 Median : 1.00 Median :2.00
## Mean : 19.4 Mean : 14.9 Mean : 1.38 Mean :2.46
## 3rd Qu.: 23.8 3rd Qu.: 12.0 3rd Qu.: 1.00 3rd Qu.:4.00
## Max. :359.7 Max. :670.0 Max. :150.00 Max. :4.00
## NAs :380
## 年龄 入院疾病分类 入院目的 住院类别 入院病情 医生
## Min. : 2.0 1:750 1:2512 1: 318 1:2602 2:1908
## 1st Qu.:23.0 2:681 2: 99 2:2307 2: 23 3: 380
## Median :36.0 3:997 3: 14 4: 337
## Mean :37.3 4:197
## 3rd Qu.:54.0
## Max. :89.0
## - 等待时间为 0 意味着立即办理住院手续。
- 门诊次为 0 意味着就诊次数 0 次。
- 住院次为 1 意味着住院次数为 1 次。
- 开住院条日期存在 380 个缺失记录，开住院条到办理住院手续有一个等待时间。这个变量对等待时间应该没有意义。
- 开住院条日期：日期数据。
- 性别：分类变量，男或女。
- 入院疾病分类：分类变量，比如肛肠疾病、皮肤病等。
- 入院目的：分类变量，比如休息疗养，等待手术等。
- 住院类别：分类变量，比如普外科、眼科等。
- 入院病情：有序分类变量，比如轻度、中度、重度。
- 医生：分类变量，比如医生 A、医生 B。
1.1 等待时间的分布
先加载 gglite 包来探索这个数据集。
library(gglite)绘制等待时间的直方图，如下图所示。
g2(hospital_waiting_time, ~等待时间)等待时间呈现指数衰减，从量级来看，时间单位很可能是小时。
对于等待时间为 0 的情况，一共有 260 条记录。怎么理解这个现象？是当作直接手术、住院呢？还是医生疏忽没有记录？类似这种涉及数据质量的问题，只有结合原始真实的数据才能判断。这里不把等待时间为 0 视为异常数据。
subdata = subset(hospital_waiting_time, subset = 等待时间 == 0)
nrow(subdata)## [1] 260这些等待时间为 0 的记录中，门诊次的分布如下：
table(subdata$门诊次)##
## 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
## 1 50 47 21 22 6 13 13 10 2 2 5 1 3 4 2 2 4 5 4
## 20 21 22 23 24 29 30 31 34 35 37 42 44 45 46 52 53 64 77 78
## 3 2 1 2 2 1 1 1 1 1 3 1 1 1 1 2 1 1 1 1
## 81 84 89 97 104 110 114 156 164 168 169 189 270 319 340
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 11.2 等待时间与门诊次的关系
g2(hospital_waiting_time, 等待时间~门诊次)从图上来看，门诊次数越多的人，等待时间明显越短一些。
1.3 等待时间与住院次的关系
g2(hospital_waiting_time, 等待时间~住院次)将横轴做对数变换，见下图。
g2(hospital_waiting_time, 等待时间~住院次) |>
scale_x(type = 'log')从图上来看，住院次数越多的人，等待时间明显越短一些。
1.4 等待时间随性别的分布
g2(hospital_waiting_time, 等待时间~性别) |>
mark_boxplot()g2(hospital_waiting_time, ~等待时间, color = ~性别) |>
mark_density()如图所示，性别对等待时间的影响不太显著。如果做检验，差别也达不到统计上的显著性。
wilcox.test(等待时间 ~ 性别, data = hospital_waiting_time)##
## Wilcoxon rank sum test with continuity correction
##
## data: 等待时间 by 性别
## W = 855248, p-value = 0.4
## alternative hypothesis: true location shift is not equal to 01.5 等待时间与患者年龄的关系
g2(hospital_waiting_time, 等待时间~年龄)从图上来看，[10-20] 和 [70+] 两个年龄段的等待时间与其它年龄段相比不一样。统计局有一个对人口分年龄段统计的数据，这里采纳统计局的年龄分段方式，分成 (0-17]、[18-34]、[35-59]、[60+] 四个年龄组，且年龄组是有序的。
hospital_waiting_time <- transform(
hospital_waiting_time,
年龄段 = cut(年龄, breaks = c(0, 17, 34, 59, 100),
# 分成左开右闭的四个年龄区间
labels = c("未成年", "青年人", "壮年人", "老年人"),
ordered_result = TRUE)
)
table(hospital_waiting_time$年龄段)
##
## 未成年 青年人 壮年人 老年人
## 462 757 979 427
levels(hospital_waiting_time$年龄段)
## [1] "未成年" "青年人" "壮年人" "老年人"
class(hospital_waiting_time$年龄段)
## [1] "ordered" "factor"g2(hospital_waiting_time, ~等待时间, color = ~年龄段) |>
mark_density()从图中来看，年龄对等待时间有一定的影响，特别是未成年人。
table(hospital_waiting_time$年龄段)##
## 未成年 青年人 壮年人 老年人
## 462 757 979 427head(hospital_waiting_time)## 等待时间 门诊次 住院次 开住院条日期 性别 年龄 入院疾病分类 入院目的 住院类别
## 1 1.0 2 1 3 0 42 3 1 2
## 2 1.2 7 1 3 0 32 1 1 2
## 3 20.0 43 1 3 1 59 1 1 2
## 4 6.0 1 1 3 1 9 3 1 2
## 5 8.9 3 1 3 1 45 3 1 2
## 6 2.9 1 1 3 1 73 3 1 2
## 入院病情 医生 年龄段
## 1 1 2 壮年人
## 2 1 2 青年人
## 3 1 2 壮年人
## 4 1 2 未成年
## 5 1 2 壮年人
## 6 1 4 老年人年龄与医生的交叉作用
g2(hospital_waiting_time, 等待时间~年龄段, color = ~医生) |>
mark_boxplot() |>
transform("dodgeX")年龄与入院疾病分类的交叉作用
g2(hospital_waiting_time, 等待时间~年龄, color = ~入院疾病分类)g2(hospital_waiting_time, 等待时间~年龄, color = ~入院病情)分箱图
g2(hospital_waiting_time, 等待时间~年龄) |>
mark_rect() |>
transform(type = "bin", color = "count",
thresholdsX= 50,
thresholdsY= 30)tooltips 显示不正常。
g2(hospital_waiting_time, 等待时间~年龄) |>
mark_point() |>
transform(type = "bin", color = "count",
thresholdsX= 50,
thresholdsY= 30)tooltips 显示正常。但是可以更加紧凑易读，示例。
1.6 等待时间随入院疾病分类的分布
g2(hospital_waiting_time, ~等待时间, color = ~入院疾病分类) |>
mark_density()1.7 等待时间随入院目的的分布
g2(hospital_waiting_time, ~等待时间, color = ~入院目的) |>
mark_density()1.8 等待时间随住院类别的分布
g2(hospital_waiting_time, ~等待时间, color = ~住院类别) |>
mark_density()1.9 等待时间随入院病情的分布
g2(hospital_waiting_time, ~等待时间, color = ~入院病情) |>
mark_density()g2(hospital_waiting_time, 等待时间 ~ 入院病情, color = ~入院病情) |>
mark_beeswarm()1.10 等待时间随医生的分布
g2(hospital_waiting_time, ~等待时间, color = ~医生) |>
mark_density()g2(hospital_waiting_time, 等待时间 ~ 医生) |>
mark_boxplot()2 广义线性模型
对于这种等待时间的建模，可以让人想到生存分析，以及指数分布、伽马分布、威布尔分布等。
# 去掉开住院条日期
hospital_waiting_time$开住院条日期 <- NULL
# 查看各个变量的类型
str(hospital_waiting_time)## 'data.frame': 2625 obs. of 11 variables:
## $ 等待时间 : num 1 1.2 20 6 8.9 2.9 7.9 2.8 2.7 5 ...
## $ 门诊次 : int 2 7 43 1 3 1 10 3 6 2 ...
## $ 住院次 : int 1 1 1 1 1 1 1 1 1 1 ...
## $ 性别 : Factor w/ 2 levels "0","1": 1 1 2 2 2 2 1 2 2 2 ...
## $ 年龄 : int 42 32 59 9 45 73 50 25 14 20 ...
## $ 入院疾病分类: Factor w/ 4 levels "1","2","3","4": 3 1 1 3 3 3 4 1 2 3 ...
## $ 入院目的 : Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
## $ 住院类别 : Factor w/ 2 levels "1","2": 2 2 2 2 2 2 2 2 2 2 ...
## $ 入院病情 : Factor w/ 2 levels "1","2": 1 1 1 1 1 1 1 1 1 1 ...
## $ 医生 : Factor w/ 3 levels "2","3","4": 1 1 1 1 1 3 1 1 3 3 ...
## $ 年龄段 : Ord.factor w/ 4 levels "未成年"<"青年人"<..: 3 2 3 1 3 4 3 2 1 2 ...根据前面的探索，假定响应变量服从指数分布，建立广义线性模型。
mod_glm <- glm(等待时间 ~ ., data = hospital_waiting_time,
family = quasi(variance = "mu^2", link = "log"))
# 固定离散参数 dispersion = 1
summary(mod_glm, dispersion = 1)##
## Call:
## glm(formula = 等待时间 ~ ., family = quasi(variance = "mu^2",
## link = "log"), data = hospital_waiting_time)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.774084 0.153443 11.56 < 2e-16 ***
## 门诊次 0.000486 0.000599 0.81 0.41758
## 住院次 -0.000186 0.006338 -0.03 0.97658
## 性别1 0.064923 0.040364 1.61 0.10774
## 年龄 -0.009371 0.003368 -2.78 0.00540 **
## 入院疾病分类2 0.556442 0.055381 10.05 < 2e-16 ***
## 入院疾病分类3 -0.110751 0.051975 -2.13 0.03310 *
## 入院疾病分类4 0.004240 0.082144 0.05 0.95884
## 入院目的2 -2.831399 0.114355 -24.76 < 2e-16 ***
## 入院目的3 -1.203519 0.274640 -4.38 1.2e-05 ***
## 住院类别2 1.507240 0.070593 21.35 < 2e-16 ***
## 入院病情2 -0.682598 0.219636 -3.11 0.00188 **
## 医生3 0.058337 0.057876 1.01 0.31347
## 医生4 -0.599757 0.060942 -9.84 < 2e-16 ***
## 年龄段.L 0.498686 0.156148 3.19 0.00140 **
## 年龄段.Q -0.158580 0.042974 -3.69 0.00022 ***
## 年龄段.C 0.074727 0.036178 2.07 0.03887 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for quasi family taken to be 1)
##
## Null deviance: 4623.3 on 2624 degrees of freedom
## Residual deviance: 3807.5 on 2608 degrees of freedom
## AIC: NA
##
## Number of Fisher Scoring iterations: 24离散参数的估计值大于 1 意味着指数分布的方差比均值大。
从结果来看，门诊次、住院次、性别、年龄等变量不太显著，也符合前面探索分析的结果。年龄不显著，但是年龄段是显著的，将连续变量切分成有一定实际意义的分类变量是合理的。
去掉不显著的变量后，看模型的结果。
mod_glm <- glm(等待时间 ~ 入院疾病分类 + 入院目的 + 住院类别 + 入院病情 + 医生 + 年龄段,
data = hospital_waiting_time,
family = quasi(variance = "mu^2", link = "log"))
summary(mod_glm, dispersion = 1)##
## Call:
## glm(formula = 等待时间 ~ 入院疾病分类 + 入院目的 +
## 住院类别 + 入院病情 + 医生 + 年龄段, family = quasi(variance = "mu^2",
## link = "log"), data = hospital_waiting_time)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.49889 0.08041 18.64 < 2e-16 ***
## 入院疾病分类2 0.53634 0.05451 9.84 < 2e-16 ***
## 入院疾病分类3 -0.12935 0.05143 -2.52 0.01190 *
## 入院疾病分类4 0.00188 0.08202 0.02 0.98169
## 入院目的2 -2.85724 0.11431 -24.99 < 2e-16 ***
## 入院目的3 -1.18269 0.27460 -4.31 1.7e-05 ***
## 住院类别2 1.49472 0.07050 21.20 < 2e-16 ***
## 入院病情2 -0.72631 0.21887 -3.32 0.00091 ***
## 医生3 0.07000 0.05777 1.21 0.22562
## 医生4 -0.58747 0.06091 -9.64 < 2e-16 ***
## 年龄段.L 0.08756 0.04782 1.83 0.06710 .
## 年龄段.Q -0.14146 0.04265 -3.32 0.00091 ***
## 年龄段.C 0.07303 0.03616 2.02 0.04344 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for quasi family taken to be 1)
##
## Null deviance: 4623.3 on 2624 degrees of freedom
## Residual deviance: 3813.5 on 2612 degrees of freedom
## AIC: NA
##
## Number of Fisher Scoring iterations: 24去掉一些变量后，Residual deviance 变化很小，但是变量达到统计显著性的明显多了。可见不同的变量对入院等待时间的影响程度是不同的，有的变量是减少作用，有的变量是增加作用。
- 门诊次对等待时间无影响。可以理解，诊断一个病情可能需要多次诊断。
- 住院次比较多的等待时间越长，可能病人并不想住院，住院多的不愿意住院。
- 医生越好，等待时间越短。年龄越大，等待时间越长，且年龄因素比医生因素对等待时间的影响更大。
- 入院病情越严重，等待时间越短。
入院疾病分类 3 相比于入院疾病分类 2、4 对应的疾病更加严重，病人才会更短时间决定住院。 入院目的 2、3 都是促使病人决定住院的因素。 入院病情 2 类似。 医生 4 相比于医生 3、2 ，可能是个专家号，病人住院的决策时间更短。
总之，所有回归系数的估计结果皆得到前面探索分析的验证。
3 运行环境
sessionInfo()## R version 4.6.0 (2026-04-24)
## Platform: aarch64-apple-darwin23
## Running under: macOS Tahoe 26.4.1
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.6/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.6/Resources/lib/libRlapack.dylib; LAPACK version 3.12.1
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## time zone: Asia/Shanghai
## tzcode source: internal
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] gglite_0.0.31
##
## loaded via a namespace (and not attached):
## [1] digest_0.6.39 R6_2.6.1 bookdown_0.46 fastmap_1.2.0
## [5] xfun_0.57 blogdown_1.23 cachem_1.1.0 knitr_1.51
## [9] htmltools_0.5.9 rmarkdown_2.31 lifecycle_1.0.5 cli_3.6.6
## [13] sass_0.4.10 jquerylib_0.1.4 compiler_4.6.0 rstudioapi_0.18.0
## [17] tools_4.6.0 evaluate_1.0.5 bslib_0.10.0 yaml_2.3.12
## [21] otel_0.2.0 jsonlite_2.0.0 rlang_1.2.0