テストデータに不明な係数レベルを持つpredict.lm()
-
28-09-2019 - |
質問
データを因数分解し、予測するためにモデルを適合させています。場合 newdata
の predict.lm()
モデルに知られていない単一因子レベルが含まれています。 すべて の predict.lm()
エラーが失敗し、返されます。
良い方法はありますか? predict.lm()
エラーのみではなく、モデルが知っている因子レベルとNAの予測を返し、モデルが知っている、NAは未知の因子レベルに対して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」に対応する3つの「実際の」予測を返したいと思います。 NA
未知のレベル「D」に対応します。
解決
整頓され、関数を拡張しました モルゲンボール. 。また実装されています Sperrorest 今。
追加機能
- 欠損値をに設定するのではなく、未使用の因子レベルをドロップします
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学習方法(列車とテストで同じレベルが必要です(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を与え、NAを与え、代替方法を使用してこれらの値を予測することができます。
物体 LMからのLM出力になります(...、data = TrainData)
データ 予測を作成するデータフレームになります
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パッケージ)のようなものを調べてください。ベイジアンモデルを使用すると、それらを推定するときに使用する情報がほとんどない場合、0に近づく効果が得られます。ただし、予測()を使用するのではなく、自分で予測する必要があるという警告。
または、モデルに含めたいレベルのダミー変数を単純に作成できます。たとえば、月曜日の可変0/1、火曜日の1つ、水曜日に1つなどです。 0です。ただし、他のデータに日曜日の列に1を持つことは、予測ステップに失敗しません。日曜日は、他の日の平均的な影響(真実である場合とそうでない場合がある場合がある場合)にあると仮定します。
線形/ロジスティック回帰の仮定の1つは、多重または複数の共有性をほとんどまたはまったくないことです。したがって、予測変数が理想的に互いに独立している場合、モデルは因子レベルのすべての種類を見る必要はありません。新しい因子レベル(d)は新しい予測因子であり、残りの因子A、B、cの予測能力に影響を与えることなくNAに設定できます。これが、モデルがまだ予測を行うことができるはずな理由です。しかし、新しいレベル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モデルを実行できます。ベイジアンパスをさらに進んで行きたい場合は、ジェルマン&ヒルの優れた本と スタン インフラストラクチャー。
スプリットテストのための迅速な解決策は、まれな値を「その他」として再調整することです。これが実装です:
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)
.