predict.lm () con un livello incognita nei dati di test
https://stackoverflow.com/questions/4285214
-
28-09-2019 - |
italiano
english
français
española
中国
日本の
العربية
Deutsch
한국어
Português
RussianDomanda
Io sono il montaggio di un modello ai dati dei fattori e predire. Se il newdata in predict.lm() contiene un unico livello di fattore che è sconosciuto al modello, tutti di predict.lm() non riesce e restituisce un errore.
C'è un buon modo per avere predict.lm() restituire una previsione per quei livelli di fattore del modello conosce e NA per i livelli incognita, invece che solo un errore?
Esempio di codice:
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)
Vorrei l'ultimo comando per tornare tre predizioni "reali" corrispondente a livelli di fattore "A", "B" e "C" ed un NA corrispondente al livello sconosciuto "D".
Soluzione
riordinato e ampliato la funzione da MorgenBall . E 'implementato anche in sperrorest ora .
Caratteristiche aggiuntive
- scende livelli di fattore non utilizzati e non solo l'impostazione dei valori mancanti per
NA. - emette un messaggio all'utente che i livelli del fattore sono state abbandonate
- controlli per l'esistenza di variabili fattore in
test_datae restituisce data.frame originale se non sono presenti - funziona non solo per
lm,glme ma anche perglmmPQL
Nota: la funzione mostrata qui può cambiare (migliorare) nel corso del tempo.
#' @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)
}
Possiamo applicare questa funzione per l'esempio nella questione come segue:
predict(model,newdata=remove_missing_levels (fit=model, test_data=foo.new))
Durante il tentativo di migliorare questa funzione, mi sono imbattuto il fatto che SL metodi come lm apprendimento, glm ecc necessità gli stessi livelli in treno e di prova, mentre ML metodi di apprendimento (svm, randomForest) non riuscire se i livelli sono stati rimossi. Questi metodi hanno bisogno di tutti i livelli in treno e di prova.
Una soluzione generale è molto difficile da realizzare, poiché ogni modello adattato ha un modo diverso di memorizzare loro componente livello fattore (fit$xlevels per lm e fit$contrasts per glmmPQL). Almeno sembra essere coerenti tra modelli correlati lm.
Altri suggerimenti
È necessario rimuovere i livelli extra prima di qualsiasi calcolo, come:
> 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
Questo è un modo più generale di farlo, sarà impostato tutti i livelli che non si verificano nei dati originali a NA. Come Hadley accennato nei commenti, che avrebbero potuto scegliere di includere questo nella funzione predict(), ma non hanno
Perché quello che dovete fare che diventa evidente se si guarda al calcolo stesso. Internamente, le previsioni sono calcolate come:
model.matrix(~predictor,data=foo) %*% coef(model)
[,1]
1 -0.1676941
2 -0.6454521
3 0.4524391
In fondo si dispone di entrambe le matrici modello. Si vede che quella per foo.new ha una colonna in più, quindi non è possibile utilizzare il calcolo della matrice più. Se si desidera utilizzare il nuovo set di dati del modello, si potrebbe anche ottenere un modello diverso, essendo uno con una variabile dummy in più per il livello supplementare.
> 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"
È possibile non solo eliminare l'ultima colonna dalla matrice del modello sia, perché, anche se si fa questo, entrambi gli altri livelli sono ancora influenzati. Il codice per A livello sarebbe (0,0). Per questo è B (1,0), per C questo (0,1) ... e per D è ancora una volta (0,0)! Così il vostro modello sarebbe supporre che A e D sono allo stesso livello se sarebbe ingenuamente cadere la variabile di ultima manichino.
In una parte più teorica: è possibile costruire un modello senza avere tutti i livelli. Ora, come ho cercato di spiegare prima, quel modello è solo valida per i livelli che hai utilizzato durante la creazione del modello. Se vi imbattete in nuovi livelli, è necessario costruire un nuovo modello per includere le informazioni supplementari. Se non lo fai, l'unica cosa che puoi fare è eliminare i livelli extra dal set di dati. Ma poi che, fondamentalmente, perde tutte le informazioni che era contenuta in esso, quindi è generalmente considerata buona pratica.
Se si vuole affrontare con i livelli mancanti nei dati dopo aver creato il vostro modello lm, ma prima di chiamare prevedere (dato che non sappiamo esattamente che cosa i livelli potrebbero mancare in anticipo) qui è la funzione che ho costruito per impostare tutti i livelli non nel modello a NA - la previsione sarà anche poi dare NA e quindi è possibile utilizzare un metodo alternativo per prevedere questi valori
.oggetto sarà la vostra uscita da lm lm (..., data = trainData)
dati sarà la cornice di dati che si desidera creare previsioni per
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
}
Suoni come te potrebbe piacere effetti casuali. Guardate in qualcosa di simile glmer (pacchetto lme4). Con un modello bayesiano, si otterrà gli effetti che si avvicinano a 0 quando c'è poca informazione da utilizzare nella stima di loro. Attenzione, però, che si dovrà fare la previsione da soli, piuttosto che usare prevedere ().
In alternativa, si può semplicemente fare fittizi variabili per i livelli che si desidera includere nel modello, per esempio una variabile 0/1 per il Lunedi, uno per il Martedì, uno per Mercoledì, ecc Domenica verrà automaticamente rimosso dal modello se contiene tutti i 0 di. Ma avere un 1 nella colonna Domenica negli altri dati non mancherà la fase di previsione. E 'solo supporre che Domenica ha un effetto che nella media degli altri giorni (che può o non può essere vero).
Una delle ipotesi di Linear / logistica regressioni è di poco o nessun multi-collinearità; quindi se le variabili predittive sono idealmente indipendenti l'uno dall'altro, quindi il modello non ha bisogno di vedere tutte le possibili varietà di livelli di fattore. Un nuovo livello fattore (D) è un nuovo predittore, e può essere impostato NA senza compromettere la capacità predire dei rimanenti fattori A, B, C. Questo è il motivo per cui il modello dovrebbe essere ancora in grado di fare previsioni. Ma aggiunta del nuovo livello D getta via lo schema previsto. Questo è l'intera questione. Impostazione correzioni NA che.
Il pacchetto lme4 gestirà nuovi livelli se si imposta la allow.new.levels=TRUE bandiera al momento della chiamata predict.
Esempio: se il giorno della settimana è il fattore in un dow variabile e un b_fail esito categorica, è possibile eseguire
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)
Questo è un esempio di un'effetti casuali regressione logistica. Naturalmente, è possibile eseguire la regressione normale ... o la maggior parte dei modelli GLM. Se si vuole la testa più in basso il percorso Bayesiano, sguardo Gelman & Hill del libro eccellente e le infrastrutture del Stan .
Una soluzione veloce-e-sporco per split test, è ricodificare valori rari come "altro". Ecco un'implementazione:
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)
}
}
Ad esempio, con data.table, la chiamata sarebbe qualcosa di simile:
dt[, (xcols) := mclapply(.SD, rare_to_other), .SDcol = xcols] # recode rare levels as other
dove xcols è un qualsiasi sottoinsieme di colnames(dt).
