Let’s fix NAs!

Довольно часто встречаются неполные наборы данных, в которых некоторые переменные не определены. В языке R содержимое таких переменных задается как «Not Available» — или сокращенно NA. Соответственно, возникает вопрос, как поступать с неопределенными значениям: стоит ли их игнорировать или откорректировать каким-либо образом?

Исследуем некоторые аспекты этой проблемы на примере практически уже ставшего классическим набора данных, взятого из соревнования Titanic: Machine Learning from Disaster. Необходимые данные можно загрузить вручную с сайта Kaggle либо же средствами R (работает под Linux):

download.file("https://bitbucket.org/kailexx/fixnas/raw/ae65f7939974e709f10aa50c96c368120487a7f2/train.csv", destfile="train.csv", method= "wget")
train <- read.csv("train.csv", na.strings = c(NA, ""))

Посмотрим, что из себя представляет содержимое файла:

str(train)
'data.frame':	891 obs. of  12 variables:
 $ PassengerId: int  1 2 3 4 5 6 7 8 9 10 ...
 $ Survived   : int  0 1 1 1 0 0 0 0 1 1 ...
 $ Pclass     : int  3 1 3 1 3 3 1 3 3 2 ...
 $ Name       : Factor w/ 891 levels "Abbing, Mr. Anthony",..: 109 191 358 277 16 559 520 629 416 581 ...
 $ Sex        : Factor w/ 2 levels "female","male": 2 1 1 1 2 2 2 2 1 1 ...
 $ Age        : num  22 38 26 35 35 NA 54 2 27 14 ...
 $ SibSp      : int  1 1 0 1 0 0 0 3 0 1 ...
 $ Parch      : int  0 0 0 0 0 0 0 1 2 0 ...
 $ Ticket     : Factor w/ 681 levels "110152","110413",..: 525 596 662 50 473 276 86 396 345 133 ...
 $ Fare       : num  7.25 71.28 7.92 53.1 8.05 ...
 $ Cabin      : Factor w/ 147 levels "A10","A14","A16",..: NA 82 NA 56 NA NA 130 NA NA NA ...
 $ Embarked   : Factor w/ 3 levels "C","Q","S": 3 1 3 3 3 2 3 3 3 1 ...

sum(is.na(train))
[1] 866

В наборе присутствуют 891 строка и аж 866 неопределенных переменных, что наглядно демонстрирует следующий график:

Два простых решения

Первое решение — для ленивых: делать практически ничего не надо, т.к. многие функции в R имеют аргумент na.rm, выставив значение которого в TRUE, мы заставим R просто удалить NA перед тем, как что-то делать с данным. Другие же функции имеют параметр na.action, который может принимать такие значения:
na.fail — возвращать ошибку, если в данных присутствуют NA.
na.omit, na.exclude — удалить все переменные со значением NA.
na.pass — оставить данные как есть.
Собственно, если нам надо подсчитать медиану возраста пассажиров Титаника, то можно поступить так:

median(train$Age, na.rm=T)
[1] 28

Второе решение, в общем-то, вытекает из первого — заранее очистить данные от всех неопределенных значений и больше не думать о них:

train.nopain <- na.omit(train)
nrow(train.nopain)
[1] 183
sum(is.na(train.nopain))
[1] 0
median(train.nopain$Age)
[1] 36

Как видим, вместе с NA испарились и 708 строк. Не страшно, если нас не особенно интересует результат, но стоит рассмотреть и другие варианты.

Используем подручные средства

Приемлемый вариант — заменить NA на какое-то заранее выбранное значение; обычно это среднее или медиана. Для этого напишем простую функцию:

simpleFix <- function(x, imputeFn=mean){
  return(ifelse(is.na(x), imputeFn(x, na.rm=TRUE), x))
}

train.median <- train
nas.idx <- which(is.na(train.median$Age))
train.median$Age <- simpleFix(train.median$Age, median)

head(train.median$Age[nas.idx])
[1] 28 28 28 28 28 28

Все неопределенные значения Age теперь заменены на медиану. Иногда есть смысл вычислять среднее значение или медиану с каким-то условием. Например, если мы корректируем NA поля Age пассажира мужского пола, который занимает каюту первого класса, то и среднее надо вычислять с учетом этого обстоятельства:

fixAge <- function(tdf, imputeFn=mean) {
  tdf$Age[is.na(tdf$Age)] <- sapply(which(is.na(tdf$Age)), 
                                    function(i) 
                                      imputeFn(tdf$Age[tdf$Pclass == tdf$Pclass[i] & 
                                                       tdf$Sex == tdf$Sex[i]], 
                                               na.rm=T))
  return(tdf)
}

nas.idx <- which(is.na(train$Age))
train.cond <- fixAge(train, median)

head(train.cond$Age[nas.idx])
[1] 25.0 30.0 21.5 25.0 21.5 25.0

Стоит отметить, что замена неопределенных значений средним или медианой приводит к уменьшению дисперсии данных. Данный метод можно несколько усовершенствовать (и усложнить), используя сингулярное разложение (SVD) и приближение матрицей меньшего ранга. При этом неопределенные значения сначала заменяются средними, потом исходная матрица аппроксимируется матрицей меньшего ранга, и неопределенные значения, которые были заменены на средние, заменяются на значения, взятые из разложения. Процедура аппроксимации повторяется несколько раз; для этого нам надо заранее задать ранг аппроксимации и количество шагов. Рассмотрим несколько сферический пример с матрицей случайных чисел размером 10х10 с 15 неопределенными значениями.

k <- 6 # Ранг аппроксимирующей матрицы
n.iters <- 10
nrows <- 10
set.seed(100500)
train.mat <- runif(nrows * nrows) 
train.mat[sample(1:length(train.mat), 15)] <- NA
train.mat <- matrix(train.mat, nrows) # Матрица случайных чисел с 15 NA значениями
nas.idx <- which(is.na(train.mat))
train.svd <- train.mat
train.svd <- apply(train.svd, 2, simpleFix) # Заменяем NA на среднее по столбцу
for (i in 1:n.iters){
  s <- svd(train.svd, k, k)
  train.svd[nas.idx] <- (s$u %*% diag(s$d[1:k], nrow=k, ncol=k) %*%  t(s$v))[nas.idx]
}

head(train.svd[nas.idx])
[1] 0.3020229 0.4475467 0.3114711 0.7161445 0.4379184 0.6734933

Что делать с нечисловыми значениями

Если внимательно рассмотреть поле Embarked, то обнаружатся 2 неопределенных значения. Один из вариантов обработки базируется на сэмплировании:

fixSample <- function(x) {
  x[is.na(x)] <- sample(x, sum(is.na(x)), replace = T)
  return(x)
}

set.seed(111)
nas.idx <- which(is.na(train.cond$Embarked))
train.cond$Embarked <- fixSample(train.cond$Embarked)

train.cond$Embarked[nas.idx]
[1] S C
Levels: C Q S

sum(is.na(train.cond$Embarked))
[1] 0

Собственно, сэмлирование является довольно распространенным методом, и весьма популярно у современных статистиков.

Более универсальный подход

R обладает развитыми средствами для построения различных моделей — от простой линейной регрессии до техники «3B» (bagging, boosting, blending). Напишем функцию, которая использует RandomForest для вычисления неопределенных переменных; также уберем некоторые переменные PassengerId, Name, Ticket, Cabin (я, если честно, не нашел им достойного применения, а в поле Cabin столько неопределенных значений, что «довычисление» теряет смысл). Если пакет randomForest у вас не установлен, то командой install.packages("randomForest") R установит его из CRAN.

fixNA <- function(y, x) {
  require(randomForest)
  fixer <- randomForest(x[!is.na(y), ], y[!is.na(y)])
  y[is.na(y)] <- predict(fixer, x[is.na(y), ])
  return(y)
}

set.seed(111)
train.rf <- subset(train, select=-c(PassengerId, Name, Ticket, Cabin))
ageNA.idx <- which(is.na(train.rf$Age))
embNA.idx <- which(is.na(train.rf$Embarked))
sum(is.na(train.rf))
train.rf$Age <- fixNA(train.rf$Age, cbind(train.rf$Pclass, train.rf$Sex, train.rf$Parch))
train.rf$Embarked <- fixNA(train.rf$Embarked, cbind(train.rf$Pclass, train.rf$Sex, train.rf$Parch))

head(train.rf$Age[ageNA.idx])
[1] 29.65873 31.67546 26.64918 29.65873 26.64918 29.65873

head(train.rf$Embarked[embNA.idx])
[1] S S

sum(is.na(train.rf))
[1] 0

В таком виде данные практически готовы для более детального анализа.

Заключение

Существует множество подходов и методов при работе с неполными данными — проблема далеко не тривиальная. В CRAN, в частности, есть и специализированные пакеты для обработки значений NA (например, Amelia, imputation), и пакеты, которые помимо прочего позволяют проводить манипуляции с NA (функции impute в пакете Hmisc, rfImpute в randomForest). Отдельного рассмотрения требует метод, основанный на алгоритме Еxpectation Maximization.

Ссылки и литература

1. Working with missing data
2. Data Imputation
3. Missing Data and Small-Area Estimation: Modern Analytical Equipment for the Survey Statistician. Nicholas T. Longford.

Автор: kxx

Источник