promict.lm () с уровнем неизвестного фактора в тестовых данных
https://stackoverflow.com/questions/4285214
-
28-09-2019 - |
italiano
english
français
española
中国
日本の
العربية
Deutsch
한국어
Português
RussianВопрос
Я устанавливаю модель для факторов данных и прогнозируя. Если то newdata в predict.lm() содержит один фактор, который неизвестен моделью, все из predict.lm() Не удается и возвращает ошибку.
Есть хороший способ иметь predict.lm() Вернуть прогноз на эти уровни фактора Модель знает и Na для неизвестных уровней фактора, а не только ошибка?
Пример кода:
foo <- data.frame(response=rnorm(3),predictor=as.factor(c("A","B","C")))
model <- lm(response~predictor,foo)
foo.new <- data.frame(predictor=as.factor(c("A","B","C","D")))
predict(model,newdata=foo.new)
Я хотел бы, чтобы самая последняя команда вернуть три «реальных» прогноза, соответствующих уровню фактора «A», «B» и «C» и NA Соответствует неизвестному уровню «D».
Решение
Вставил и расширил функцию Morgenball. Отказ Это также реализовано в СПРАМЕРЕСТ Теперь.
Дополнительные возможности
- Откатывает неиспользованные уровни фактора, а не просто устанавливать недостающие значения для
NA. - выдает сообщение пользователям, что уровни фактора были сброшены
- проверки на существование переменных фактора в
test_dataи возвращает исходные данные. Сделайте, если не присутствуют - работает не только для
lm,glmи но и дляglmmPQL
Примечание. Показанная здесь функция может измениться (улучшить) со временем.
#' @title remove_missing_levels
#' @description Accounts for missing factor levels present only in test data
#' but not in train data by setting values to NA
#'
#' @import magrittr
#' @importFrom gdata unmatrix
#' @importFrom stringr str_split
#'
#' @param fit fitted model on training data
#'
#' @param test_data data to make predictions for
#'
#' @return data.frame with matching factor levels to fitted model
#'
#' @keywords internal
#'
#' @export
remove_missing_levels <- function(fit, test_data) {
# https://stackoverflow.com/a/39495480/4185785
# drop empty factor levels in test data
test_data %>%
droplevels() %>%
as.data.frame() -> test_data
# 'fit' object structure of 'lm' and 'glmmPQL' is different so we need to
# account for it
if (any(class(fit) == "glmmPQL")) {
# Obtain factor predictors in the model and their levels
factors <- (gsub("[-^0-9]|as.factor|\\(|\\)", "",
names(unlist(fit$contrasts))))
# do nothing if no factors are present
if (length(factors) == 0) {
return(test_data)
}
map(fit$contrasts, function(x) names(unmatrix(x))) %>%
unlist() -> factor_levels
factor_levels %>% str_split(":", simplify = TRUE) %>%
extract(, 1) -> factor_levels
model_factors <- as.data.frame(cbind(factors, factor_levels))
} else {
# Obtain factor predictors in the model and their levels
factors <- (gsub("[-^0-9]|as.factor|\\(|\\)", "",
names(unlist(fit$xlevels))))
# do nothing if no factors are present
if (length(factors) == 0) {
return(test_data)
}
factor_levels <- unname(unlist(fit$xlevels))
model_factors <- as.data.frame(cbind(factors, factor_levels))
}
# Select column names in test data that are factor predictors in
# trained model
predictors <- names(test_data[names(test_data) %in% factors])
# For each factor predictor in your data, if the level is not in the model,
# set the value to NA
for (i in 1:length(predictors)) {
found <- test_data[, predictors[i]] %in% model_factors[
model_factors$factors == predictors[i], ]$factor_levels
if (any(!found)) {
# track which variable
var <- predictors[i]
# set to NA
test_data[!found, predictors[i]] <- NA
# drop empty factor levels in test data
test_data %>%
droplevels() -> test_data
# issue warning to console
message(sprintf(paste0("Setting missing levels in '%s', only present",
" in test data but missing in train data,",
" to 'NA'."),
var))
}
}
return(test_data)
}
Мы можем применить эту функцию к примеру в вопросе следующим образом:
predict(model,newdata=remove_missing_levels (fit=model, test_data=foo.new))
Пытаясь улучшить эту функцию, я наткнулся на тот факт, что методы обучения SL, как lm, glm И т. Д. Нужен те же уровни в поезде и тестирование, а методы обучения ML (svm, randomForest) не удалось, если уровни удаляются. Эти методы нуждаются в всех уровнях в поезде и тест.
Общее решение довольно сложно достичь, поскольку каждая установленная модель имеет другой способ хранения их компонента коэффициента уровня (fit$xlevels за lm и fit$contrasts за glmmPQL). По крайней мере, это, кажется, соответствует lm Связанные модели.
Другие советы
Вы должны снять дополнительные уровни перед любым расчетом, например:
> id <- which(!(foo.new$predictor %in% levels(foo$predictor)))
> foo.new$predictor[id] <- NA
> predict(model,newdata=foo.new)
1 2 3 4
-0.1676941 -0.6454521 0.4524391 NA
Это более общий способ сделать это, он установит все уровни, которые не возникают в исходных данных к NA. Как упоминалось в комментариях Хадли, они могли бы выбрать это в predict() функция, но они не
Почему вы должны сделать, это становится очевидным, если вы посмотрите на сам расчет. Внутренне прогнозы рассчитаны как:
model.matrix(~predictor,data=foo) %*% coef(model)
[,1]
1 -0.1676941
2 -0.6454521
3 0.4524391
Внизу у вас есть оба модельные матрицы. Вы видите, что один для foo.new Имеет дополнительную колонну, поэтому вы не можете использовать расчет матрицы. Если вы будете использовать новый набор данных для модели, вы также получите другую модель, являющуюся одной с дополнительной фиктивной переменной для дополнительного уровня.
> model.matrix(~predictor,data=foo)
(Intercept) predictorB predictorC
1 1 0 0
2 1 1 0
3 1 0 1
attr(,"assign")
[1] 0 1 1
attr(,"contrasts")
attr(,"contrasts")$predictor
[1] "contr.treatment"
> model.matrix(~predictor,data=foo.new)
(Intercept) predictorB predictorC predictorD
1 1 0 0 0
2 1 1 0 0
3 1 0 1 0
4 1 0 0 1
attr(,"assign")
[1] 0 1 1 1
attr(,"contrasts")
attr(,"contrasts")$predictor
[1] "contr.treatment"
Вы не можете просто удалить последний столбец из матрицы модели, потому что даже если вы это сделаете, оба других уровня все еще влияют. Код для уровня A будет (0,0). За B Это (1,0), для C Это (0,1) ... и для D Это снова (0,0)! Так что ваша модель предположила бы, что A и D находятся на одном уровне, если он наивно бросит последнюю фиктивную переменную.
На более теоретической части: можно построить модель без всех уровней. Теперь, как я пытался объяснить раньше, эта модель Только Действительно для уровней, которые вы использовали при построении модели. Если вы нашли новые уровни, вы должны создать новую модель, чтобы включить дополнительную информацию. Если вы этого не сделаете, единственное, что вы можете сделать, это удалить дополнительные уровни из набора данных. Но тогда вы в основном теряете всю информацию, которая содержится в ней, так что в целом не считается хорошей практикой.
Если вы хотите иметь дело с недостающимися уровнями в ваших данных после создания модели LM, но до прогнозирования вызова (учитывая мы не знаем, какие уровни могут отсутствовать заранее) Вот функция, которую я создал, чтобы установить все уровни не в Модель на Na - прогноз также дает Na, и вы можете использовать альтернативный метод для прогнозирования этих значений.
объект Будет ваш выход LM от LM (..., Data = rundata)
данные будет кадр данных, который вы хотите создать прогнозы для
missingLevelsToNA<-function(object,data){
#Obtain factor predictors in the model and their levels ------------------
factors<-(gsub("[-^0-9]|as.factor|\\(|\\)", "",names(unlist(object$xlevels))))
factorLevels<-unname(unlist(object$xlevels))
modelFactors<-as.data.frame(cbind(factors,factorLevels))
#Select column names in your data that are factor predictors in your model -----
predictors<-names(data[names(data) %in% factors])
#For each factor predictor in your data if the level is not in the model set the value to NA --------------
for (i in 1:length(predictors)){
found<-data[,predictors[i]] %in% modelFactors[modelFactors$factors==predictors[i],]$factorLevels
if (any(!found)) data[!found,predictors[i]]<-NA
}
data
}
Похоже, вы можете понравиться случайные эффекты. Посмотрите на что-то вроде GLMER (пакет LME4). С помощью модели Bayesian вы получите эффекты, которые подходят 0, когда есть небольшая информация для использования при их оценке. Предупреждение, хотя, чтобы вам придется сделать прогноз самого себя, а не использовать прогноз ().
В качестве альтернативы, вы можете просто сделать фиктивные переменные для уровней, которые вы хотите включить в модель, например, переменную 0/1 для понедельника, один на вторник, по одной в среду и т. Д. Воскресенье будет автоматически удалено из модели, если она содержит все 0 Но имея 1 в воскресной колонке в других данных, не пройдет неудачу на шаг прогнозирования. Он просто предположил, что в воскресенье есть эффект, который в среднем в другие дни (что может или не может быть верным).
Одним из предположений линейных / логистических регрессий является практически без мультиолицереации; Поэтому, если переменные предикторов в идеале независимы друг от друга, то модель не должна видеть все возможное множество уровней фактора. Новый уровень фактора (D) является новым предиктором и может быть установлен на NA, не влияя на прогнозную способность оставшихся факторов A, B, C. Вот почему модель все еще должна иметь возможность делать прогнозы. Но добавление нового уровня D бросает ожидаемую схему. Это весь вопрос. Установка NA исправляет это.
То lme4 Пакет будет обрабатывать новые уровни, если вы установите флаг allow.new.levels=TRUE при звонке predict.
Пример: если ваш день недельного фактора в переменной dow и категорический результат b_fail, вы могли бы запустить
M0 <- lmer(b_fail ~ x + (1 | dow), data=df.your.data, family=binomial(link='logit'))
M0.preds <- predict(M0, df.new.data, allow.new.levels=TRUE)
Это пример с логистической регрессией случайных эффектов. Конечно, вы можете выполнить регулярную регрессию ... или большинство моделей GLM. Если вы хотите отправиться дальше вниз по пути Байеса, посмотрите на отличную книгу Gelman & Hill и Соревнование инфраструктура.
Быстрое и грязное решение для разделения тестирования, состоит в том, чтобы восстановить редкие ценности как «другие». Вот реализация:
rare_to_other <- function(x, fault_factor = 1e6) {
# dirty dealing with rare levels:
# recode small cells as "other" before splitting to train/test,
# assuring that lopsided split occurs with prob < 1/fault_factor
# (N.b. not fully kosher, but useful for quick and dirty exploratory).
if (is.factor(x) | is.character(x)) {
min.cell.size = log(fault_factor, 2) + 1
xfreq <- sort(table(x), dec = T)
rare_levels <- names(which(xfreq < min.cell.size))
if (length(rare_levels) == length(unique(x))) {
warning("all levels are rare and recorded as other. make sure this is desirable")
}
if (length(rare_levels) > 0) {
message("recoding rare levels")
if (is.factor(x)) {
altx <- as.character(x)
altx[altx %in% rare_levels] <- "other"
x <- as.factor(altx)
return(x)
} else {
# is.character(x)
x[x %in% rare_levels] <- "other"
return(x)
}
} else {
message("no rare levels encountered")
return(x)
}
} else {
message("x is neither a factor nor a character, doing nothing")
return(x)
}
}
Например, с data.table, звонок будет чем-то вроде:
dt[, (xcols) := mclapply(.SD, rare_to_other), .SDcol = xcols] # recode rare levels as other
куда xcols это любое подмножество colnames(dt).
