# using cluster information in multiple imputation

i need to impute a dataset all categorical variables before doing analysis. I can just simply impute with mode of all data or a variable.

I belief that better option will be to classify the subjects (observations) using some short of clustering algorithm and then use this information to impute the data. The following is small data (though the real data is really big) and the idea.

md <- data.frame(V1 = c("AA", "AA", "AA", NA, "AB"), V2 = c("AB", "AB", "BB", "BB", "BB"), V3 = c("BB", NA, "BB", "BB", "BB"), V4 = c("AA", "AA", "AA", "AA", "AA"),
V5=c(rep("AB", 5)), V6 = c("BB", "BB", "AB", "AB", NA), V7 = c("AB", "AB", "BB", "BB", "BB"))

md
V1 V2   V3 V4 V5   V6   V7
1   AA AB   BB AA   AB   BB AB
2   AA AB <NA> AA   AB   BB AB
3   AA BB   BB AA   AB   AB BB
4 <NA> BB   BB AA   AB   AB BB
5   AB BB   BB AA   AB <NA> BB


In visual observation, the sample 1 and 2 are more similar in one cluster while 3,4,5 are in second cluster. The side is hand drawn dendogram (I could not due HC because of missing values).

Now I would like to impute all missing values based on the similarity. For example, missing value in column 1 is most likely to be AA as sample 4 is more similar to 3 than 1, 2, or 5. Similarly the 5 column missing value is BB as the neighbor in the cluster is also BB. Similarly column 6 value should be AB as the closest similar has AB. and so on. Thus completed data would look like: How can we perform this ?

• You're looking to do single rather than multiple imputation? Are you familiar with eg MICE? – generic_user Jul 10 '14 at 22:40
• multiple imputation might be good, but I do not know if works for categorical data – John Jul 10 '14 at 23:26
• MICE: multiple imputation by chained equations. When you have hierarchical data, I think the way it works is that you specify your models at each level. In your case, you'd probably want multinomial logits to predict the category. I.e.: pr(category|level, whatever else). It can be a pain to code, easy to miscode, and computationally intensive, but it's pretty flexible. And I guess that you'd need to address the uncertainty in your clustering algorithm as well. It might not be the best choice. – generic_user Jul 10 '14 at 23:35
• I agree with @ACD that a MICE package can be used to impute categorical data like yours. I don't think "hierarchical data" applies here, which simplifies your task. Still, doing it in R is no simple task. – rolando2 Jul 11 '14 at 0:26

The following function is based on the paper "Imputation of Missing Values for Unsupervised Data Using the Proximity in Random Forests" by Tsunenori Ishioka in eLmL 2013. Please follow the paper for methodological details. Note that this function is applicable for both categorical and numerical data.

rfunsuper <- function (x, iter=5, ntree=100){
# three different functions
#   Return k-neighbor weighted mean for numeric variables or
#   most weighted frequent factor element for factor variables.
KWmean <- function (value, weight, k=10){
if (missing(weight)){
w <- rep.int(1, length(value))
}else if (length(weight) != length(value)){
stop("'value' and 'weight' must have the same length")
}
k <- min(k, length(value))
if (is.numeric(value)){
order.weight <- order(weight, decreasing = T)
ww <- weight[order.weight]
vv <- value[order.weight]
ret <- sum(ww[1:k] * vv[1:k]) / sum(ww[1:k])

}else if(is.factor(value)){
wgt.sum <- tapply(weight, value, sum)
# most weighted frequent factor element
ret <- names(subset (wgt.sum, wgt.sum == max(wgt.sum, na.rm=T)))
}else{
stop("'value' is neither numeric nor factor")
}
return(ret)
}
#   Return relative distance between x.impute' to x.org'
relatImpOrg <- function (x.impute, x.org){
#   x.impute: imputed data
#   x.org: original data
ncol.x <- length(x.org)
x.abs.org <-  matrix(abs(as.numeric(unlist(x.org))), ncol=ncol.x)
max.x <- apply(x.abs.org, 2, max) # for normalization of features size
# x.impute' and x.org' may include factor elements
if (FALSE){ # available for only numeric
diff.x <- (x.impute - x.org) / max.x # normalize
diff.rel <- sum(diff.x^2) / sum((x.org / max.x)^2)
}else{
mat.x.impute <- matrix(as.numeric(unlist(x.impute)), ncol=ncol.x)
mat.x.org <- matrix(as.numeric(unlist(x.org)), ncol=ncol.x)
max.numx <- as.numeric(unlist(max.x))
diff.x <- sweep((mat.x.impute - mat.x.org), 2, max.numx, FUN="/")
size.org <- sweep(mat.x.org, 2, max.numx, FUN="/")
diff.rel <- sum(diff.x^2) / sum(size.org^2)
}
cat ("diff.rel =", sum(diff.x^2), "/", sum(size.org^2), "=", diff.rel, "\n")
return(diff.rel)
}
#   Impute or revise NA elements using the data proximity.
prox.nafix <- function (na.values, rough.values, x.prox){
#   na.values: data vector that includes NA; unchanged.
#   rough.values: rough data vector to be replaced; NAs cannot include.
#   x.prox: data proximity matrix; each element is positive and <= 1.
if (length(na.values) != length(rough.values)){
stop("'na.values' and 'rough.values' must have the same length");
}else if (length(rough.values) != ncol(x.prox)){
stop("'rough.values' and 'x.prox' size incorrect");
}
# NA imputation ONLY for NA data
na.list <- which(is.na(na.values))
if (length(na.list) == 0){
# no NAs
return(rough.values)
}
replaced.vales <- rough.values
for (i in 1:length(na.list)){
j <- na.list[i]
x.prox[j,j] <- 0 # ignore the weight of the data to be imputed.
replaced.vales[j] <- KWmean (rough.values, x.prox[,j])
}
return(replaced.vales)
}
require(randomForest)
x.roughfixed <- na.roughfix(x)
#For numeric variables, NAs are replaced with column medians.
#For factor variables, NAs are replaced with the most frequent levels (breaking ties at random).
rf.impute <- x
while (iter){
x.rf <- randomForest(x.roughfixed, ntree=ntree)
#randomForest implements Breiman's random forest algorithm
x.prox <- x.rf\$proximity
#a matrix of proximity measures among the input
#(based on the frequency that pairs of data points are in the same terminal nodes).
for (i in 1:ncol(x)){
rf.impute[,i] <- prox.nafix(x[,i], x.roughfixed[,i], x.prox)
}
diff.rel <- relatImpOrg(rf.impute, x.roughfixed)
if (diff.rel < 1e-5){
break
}else{
x.roughfixed <- rf.impute
iter <- iter -1
}
}
print(x.rf)
return(rf.impute)
}


Now let's implement in your situation.

library(randomForest)
md <- data.frame(V1 = c("AA", "AA", "AA", NA, "AB"),
V2 = c("AB", "AB", "BB", "BB", "BB"),
V3 = c("BB", NA, "BB", "BB", "BB"),
V4 = c("AA", "AA", "AA", "AA", "AA"),
V5=c(rep("AB", 5)),
V6 = c("BB", "BB", "AB", "AB", NA),
V7 = c("AB", "AB", "BB", "BB", "BB"))

rfunsuper (md, iter=5, ntree=100)


The imputed output exactly match to your expectations:

diff.rel = 0 / 28.11111 = 0

Call:
randomForest(x = x.roughfixed, ntree = ntree)
Type of random forest: unsupervised
Number of trees: 100
No. of variables tried at each split: 2

V1 V2 V3 V4 V5 V6 V7
1 AA AB BB AA AB BB AB
2 AA AB BB AA AB BB AB
3 AA BB BB AA AB AB BB
4 AA BB BB AA AB AB BB
5 AB BB BB AA AB AB BB