０が多いカウントデータの回帰分析 ZIP Rで分析

０が多いようなカウントデータを分析したり、予測したりするには１手間必要だったりします。例えば　六本木で働くデータサイエンティストのブログとかではポワソン回帰・負の二項分布の使いわけを説明されています。
ただ実業務で割と多いのが、０回がほとんどで１回がちょっと多い、それからさらに少なくなるといった以下のような対数正規分布みたいなデータを扱う事が多く、なんかいいやり方ないかなぁと調べるとZIP(Zero-inflated model)というのがあるようなので備忘録がてら比較しておきます。

library(pscl)
library(boot)
data(fir)
set.seed(123)

table(fir$count)
# 学習用と検証用でデータ分割
learn_idx <- sample(1:nrow(fir),nrow(fir)*0.8)
learn_data <-  fir[learn_idx,]
eval_data <- fir[-learn_idx,]


## ポワソン回帰/ZIP/通常回帰を比較
glm1 <- glm(count ~ ., data = learn_data, family =  poisson)
lm_result <- lm(count ~ ., data = learn_data)
zip <- zeroinfl(count ~ . | ., data = learn_data)
lm_pred_result <- predict(lm_result,eval_data)
zip_pred_result <- predict(zip,eval_data)
glm_pred_result <- predict(glm1,eval_data)

mean(abs(lm_pred_result-eval_data$count))
mean(abs(zip_pred_result-eval_data$count))
mean(abs(glm_pred_result-eval_data$count))

AIC(lm_result)
AIC(glm1)
AIC(zip)

# もうちょっと０が多いデータで再チャレンジ
data(bioChemists)
table(bioChemists$art)

# 学習用と検証用でデータ分割
learn_idx <- sample(1:nrow(bioChemists),nrow(bioChemists)*0.8)
learn_data <-  bioChemists[learn_idx,]
eval_data <- bioChemists[-learn_idx,]

glm1 <- glm(art ~ ., data = learn_data, family = poisson)
lm_result <- lm(art ~ ., data = learn_data)
zip <- zeroinfl(art ~ . | ., data = learn_data,EM=TRUE)
lm_pred_result <- predict(lm_result,eval_data)
zip_pred_result <- predict(zip,eval_data)
glm_pred_result <- predict(glm1,eval_data)

mean(abs(lm_pred_result-eval_data$art))
mean(abs(zip_pred_result-eval_data$art))
mean(abs(glm_pred_result-eval_data$art))

AIC(lm_result)
AIC(glm1)
AIC(zip)

psclライブラリーのzeroinfl関数で実行できます。ポワソン回帰・線形回帰・ZIPで比較します。bootライブラリーはデータ用のライブラリーです。
最初にfirデータのカウントを予測しますが、firデータは0以外が多いときはポワソン回帰の性能がよさそうですが、AICだと微妙にZIPのほうがよさげです。bioChemistsデータは０が多いものになるとやはりZIPのほうが精度的にも良くなります。