有一个非常有趣的实用解决方法
除了一些在网上找到的过度清理的数据集外，缺失值无处不在。实际上，数据集越复杂和庞大，缺失值的出现就越有可能。缺失值是一个迷人的统计研究领域，但在实践中，它们经常是一种麻烦。
如果您处理的预测问题中，您想要从 p 维协变量 X =(X_1,…,X_p) 预测变量 Y，并且您面临着 X 中的缺失值，则树形方法有一个有趣的解决方案。这种方法实际上相当古老，但在各种数据集中似乎表现得非常出色。我正在谈论的是“缺失值合并到属性准则”（MIA; [1]）。虽然有许多关于缺失值的好文章（例如这篇文章），但这种强大的方法似乎被少用。特别地，您不需要以任何方式插补、删除或预测您的缺失值，而是可以像有完全观测数据一样运行您的预测。
我将快速解释一下这种方法本身的工作原理，并提供一个使用分布式随机森林（DRF）的示例，该分布式随机森林在此处进行了解释。我选择 DRF，因为它是随机森林的一个非常通用的版本（特别是它也可以用于预测随机向量 Y），并且因为我在这里有些偏见。MIA 实际上是为广义随机森林（GRF）实现的，它涵盖了广泛的森林实现。特别是，由于 CRAN 上的 DRF 实现是基于 GRF 的，在稍作修改后，它也可以使用 MIA 方法。
当然，要注意这是一个快速修复，（据我所知）没有理论保证。根据缺失机制的不同，它可能会严重影响分析结果。另一方面，处理缺失值的大多数常用方法都没有任何理论保证，或者已知会影响分析结果，至少在经验上，MIA 似乎表现良好且
工作原理
回想一下，在随机森林中，分裂的形式是 X_j < S 或 X_j ≥ S，其中 j=1,…,p。为了找到这个分裂值 S，它会基于 Y 的某种标准进行优化，例如 CART 标准。因此，观测值会通过依赖于 X 的决策规则逐步分裂。
原始论文有点混乱，但据我所知，MIA 的工作原理如下：让我们考虑一个样本 ( Y_1 , X _1),…, (Y_n, X _n)，其中
X _i=(X_i1,…,X_ip)’。
没有缺失值的拆分只是像上面那样查找值 S，然后将所有 X_ij < S 的 Y_i 放到节点 1 中，将所有 X_ij ≥ S 的 Y_i 放到节点 2 中。为每个值 S 计算目标标准（例如 CART），我们可以选择最佳值。对于缺失值，每个候选拆分值 S 都有三个选项：
- 对于所有观测值 i，如果 X_ij 是已观测到的，则使用通常的规则并将其发送到节点 1，如果 X_ij 是缺失的，则将其发送到节点 1。
- 对于所有观测值 i，如果 X_ij 是已观测到的，则使用通常的规则并将其发送到节点 2，如果 X_ij 是缺失的，则将其发送到节点 2。
- 忽略通常的规则，如果 X_ij 是缺失的，则将 i 发送到节点 1，如果它是已观测到的，则将其发送到节点 2。
我们要遵循哪些规则取决于我们使用的 Y_i 标准。
一个小例子
需要在此指出的是，CRAN上的drf软件包尚未更新为最新的方法。未来将会有一点，其中所有这些都将在CRAN上实现一个软件包(!)。然而，目前有两个版本：
如果您想使用带有缺失值的快速drf实现（不包括置信区间），则可以使用本文末尾附加的“drfown”函数。这个代码是从
lorismichel/drf: Distributional Random Forests (Cevid et al., 2020) (github.com)中改编的
另一方面，如果您想要使用您的参数置信区间，则使用这个（更慢）代码
drfinference/drf-foo.R at main · JeffNaef/drfinference (github.com)
特别是，在后一种情况下，drf-foo.R包含您所需的所有内容。
我们将重点放在具有置信区间的较慢代码上，如本文所述，并考虑与上述文章中相同的示例：
set.seed(2)n<-2000beta1<-1beta2<--1.8# 模型模拟X<-mvrnorm(n = n, mu=c(0,0), Sigma=matrix(c(1,0.7,0.7,1), nrow=2,ncol=2))u<-rnorm(n=n, sd = sqrt(exp(X[,1])))Y<- matrix(beta1*X[,1] + beta2*X[,2] + u, ncol=1)请注意，这是一个异方差线性模型，其中p=2，误差项的方差取决于X_1的值。现在，我们还以Missing at Random (MAR)的方式向X_1添加丢失的值：
prob_na <- 0.3X[, 1] <- ifelse(X[, 2] <= -0.2 & runif(n) < prob_na, NA, X[, 1]) 这意味着当X_2的值小于-0.2时，X_1的丢失概率为0.3。因此，X_1的丢失概率取决于X_2，这就是所谓的“Missing at Random”。这已经是一个复杂的情况，通过查看缺失值的模式可以获得信息。也就是说，缺失不是“Missing Completely at Random (MCAR)”，因为X_1的缺失取决于X_2的值。这反过来意味着我们从中抽取的X_2的分布是不同的，具体取决于X_1是否丢失。这特别意味着删除具有缺失值的行可能会严重影响分析结果。
我们现在固定x并估算给定X = x的条件期望和方差，就像上一篇文章中一样。
# 选择不太远的x值x<-matrix(c(1,1),ncol=2)# 选择CIs的alpha值alpha<-0.05然后我们还拟合DRF并预测测试点x的权重（对应于预测Y | X=x的条件分布）：
## 拟合新的DRF框架drf_fit <- drfCI(X=X, Y=Y, min.node.size = 5, splitting.rule='FourierMMD', num.features=10, B=100)## 预测权重DRF = predictdrf(drf_fit, x=x)weights <- DRF$weights[1,]例1：条件期望
我们首先估计Y | X = x的条件期望。
# Estimate the conditional expectation at x:condexpest<- sum(weights*Y)# Use the distribution of weights, see belowdistofcondexpest<-unlist(lapply(DRF$weightsb, function(wb) sum(wb[1,]*Y) ))# Can either use the above directly to build confidence interval, or can use the normal approximation.# We will use the lattervarest<-var(distofcondexpest-condexpest)# build 95%-CIlower<-condexpest - qnorm(1-alpha/2)*sqrt(varest)upper<-condexpest + qnorm(1-alpha/2)*sqrt(varest)round(c(lower, condexpest, upper),2)# without NAs: (-1.00, -0.69 -0.37)# with NAs: (-1.15, -0.67, -0.19)值得注意的是，带有NAs的值非常接近于前一篇文章中没有NAs的第一次分析中的值！这对我来说真的很令人惊讶，因为这个缺失机制不容易处理。有趣的是，估计量的估计方差也翻了一番，从没有缺失值的0.025左右增加到有缺失值的大约0.06左右。
事实如下：
所以我们有一个轻微的误差，但置信区间包含真实值，就像它们应该的那样。
对于更复杂的目标，如条件方差，结果看起来类似：
# Estimate the conditional expectation at x:condvarest<- sum(weights*Y^2) - condexpest^2distofcondvarest<-unlist(lapply(DRF$weightsb, function(wb) { sum(wb[1,]*Y^2) - sum(wb[1,]*Y)^2} ))# Can either use the above directly to build confidence interval, or can use the normal approximation.# We will use the lattervarest<-var(distofcondvarest-condvarest)# build 95%-CIlower<-condvarest - qnorm(1-alpha/2)*sqrt(varest)upper<-condvarest + qnorm(1-alpha/2)*sqrt(varest)c(lower, condvarest, upper)# without NAs: (1.89, 2.65, 3.42)# with NAs: (1.79, 2.74, 3.69)这里估计值之间的差异稍大。由于事实如下：
带有NAs的估计值甚至更精确（尽管这当然可能只是随机性）。同样，（方差）估计量的方差估计随着缺失值的增加而增加，从没有缺失值的0.15增加到0.23。
结论
在本文中，我们讨论了MIA，它是Random Forest中拆分方法的一种适应缺失值的方法。由于它在GRF和DRF中实现，因此它可以广泛使用，我们所看到的小例子表明它的效果非常好。
然而，我想再次指出，即使对于非常大的数据点，也没有一致性或置信区间有意义的理论保证。缺失值的原因很多，必须非常小心，以免通过粗心处理此问题而产生偏差。MIA方法绝不是这个问题的一个被充分理解的修复程序。但是，它似乎是目前一个合理的快速修复程序，可以利用数据中的缺失模式。如果有人进行了/有更广泛的模拟分析，我会对结果感到好奇。
代码
require(drf) drfown <- function(X, Y, num.trees = 500, splitting.rule = "FourierMMD", num.features = 10, bandwidth = NULL, response.scaling = TRUE, node.scaling = FALSE, sample.weights = NULL, sample.fraction = 0.5, mtry = min(ceiling(sqrt(ncol(X)) + 20), ncol(X)), min.node.size = 15, honesty = TRUE, honesty.fraction = 0.5, honesty.prune.leaves = TRUE, alpha = 0.05, imbalance.penalty = 0, compute.oob.predictions = TRUE, num.threads = NULL, seed = stats::runif(1, 0, .Machine$integer.max), compute.variable.importance = FALSE) { # 对X和Y进行初始检查 if (is.data.frame(X)) { if (is.null(names(X))) { stop("如果以data.frame格式提供，回归器应该命名。") } if (any(apply(X, 2, class) %in% c("factor", "character"))) { any.factor.or.character <- TRUE X.mat <- as.matrix(fastDummies::dummy_cols(X, remove_selected_columns = TRUE)) } else { any.factor.or.character <- FALSE X.mat <- as.matrix(X) } mat.col.names.df <- names(X) mat.col.names <- colnames(X.mat) } else { X.mat <- X mat.col.names <- NULL mat.col.names.df <- NULL any.factor.or.character <- FALSE } if (is.data.frame(Y)) { if (any(apply(Y, 2, class) %in% c("factor", "character"))) { stop("Y只应包含数字变量。") } Y <- as.matrix(Y) } if (is.vector(Y)) { Y <- matrix(Y,ncol=1) } #validate_X(X.mat) if (inherits(X, "Matrix") && !(inherits(X, "dgCMatrix"))) { stop("目前只支持类'dgCMatrix'的稀疏数据。") } drf:::validate_sample_weights(sample.weights, X.mat) #Y <- validate_observations(Y, X) #设置GRF的遗留参数 clusters <- vector(mode = "numeric", length = 0) samples.per.cluster <- 0 equalize.cluster.weights <- FALSE ci.group.size <- 1 num.threads <- drf:::validate_num_threads(num.threads) all.tunable.params <- c("sample.fraction", "mtry", "min.node.size", "honesty.fraction", "honesty
引用
[1] Twala, B. E. T. H., M. C. Jones, 和 David J. Hand. Good methods for coping with missing data in decision trees. Pattern Recognition Letters 29, 2008.