備忘ログ

チラシの裏的備忘録&メモ

R 3.6.3からR 4.0.2周りで変更になった因子の扱いについての気づき（ハマった）

R 生存曲線 因子型

mikuhatsune.hatenadiary.com

この記事をみて、面白いなぁと思ってやってみたらうまくいったり行かなかったりした（要するにはまった）のでメモ。本旨はR 3.6.3からR 4.0.2周りで変更になった因子の扱いについての気づき。ブックマークしてコメントにメモしようと思ったけど文量が多くなったのでこっちに書く。

3.6.3から4.0.2で色々変更になっていたのは知っていたけど、ちゃんとリファレンス読んでないからこういうのに今更ハマる。

ブログ内で指定されたパッケージ以外にも色の設定のためのjet.colors関数のために{matlab}パッケージが必要とかいうのもあったけど、そんなのは

library(matlab)

で読み込んでおけばいいのでそんなのは問題ではない。

あとは、ブログの下に書かれているCSVファイル用のデータが貼ってあるが、すでに前処理済みなのも大した問題じゃない。 コピペしたCSVファイルのsurvive、days、yearを消して、titileは「ダンキラ!!! - Boys, be DANCING! -」がコロン含まれているので処理に注意し、d1の2020-09-16を消して、CSVファイルを作りCSVファイルを読み込むときに

dat <- read.csv("sns_survival.csv", stringsAsFactors = FALSE, na.strings = "")

とna.strings = ""を指定するとおそらくコードで実施した前処理的なのがなされる前のデータが読み込める。

そいうのをクリアしてもそのままのコードでうまくいくパターンとうまく行かないパターンがあり、うまくいくパターンとしては

  • Rの3.6.3（3.6系最終バージョン）でやる

そのままのコードではうまく行かないパターンとしては

  • Rの4.0.2（現行最新バージョン）でやる

でなんでうまく行かないのかということについて見てみると、4.0.xで変更になった因子型の取り扱い周りな気がする。

dat$year <- cut(as.numeric(format(dat$d0, "%Y")), ycut, labels=tail(ycut, -1))

dat0 <- rbind(cbind.data.frame(subset(dat, select=-year), year="All"), dat)

のところで、dat$yearがラベルがlabels = tail(ycut, -1)で指定された因子型になるのは3.6.3、4.0.2ともに共通なのだけれど、dat0でdatの全データのyearにAllを入れたデータを入れた上で、datのもとのデータをくっつけようとしたときに、dat0$yearが3.6.3の場合には因子型のまま処理できているのだけれど、4.0.2では文字列型になってしまってその後の生存曲線を描く。

plot(fit, newdata=data.frame(year=levels(dat0$year)), col=cols,
     marktime=TRUE,
     axis1.at=at.t, axis1.labels=seq(at.t)-1, atrisk.at=at.t, atrisk.font=2,
     xlab="Time to service stop [year]",
     lwd=4,
     legend.x="topright", legend.cex=1, legend=FALSE
)

でnewdata = data.frame(year= levels(dat0$year))で3.6.3だとdat0$yearが因子型なのでレベルもそのままもってこれるのだけれど、4.0.2だと文字列になってしまっているためうまくいかなくなってしまう。

R 4.0.2でも生存曲線を描けるようにするために安直な解決策としては

dat0$year <- factor(dat0$year, levels = c("All", 2010:2020))

と途中でdat0$yearを因子型変えてやり、生存曲線を描画するところで

par(mar=c(4, 3, 2, 2), cex.lab=1.5, cex.axis=1.2)
plot(fit, newdata=data.frame(year=levels(dat0$year)), col=cols)

とすると

f:id:indenkun:20200921195141p:plain

と似たような生存曲線がR4.0.2でも描ける。あとはプロットの中身をもっと書けばさらに近づく。

dat0$yearのlevelsをAll、2010～2020の順に明示的にしているのに、なぜかyearが2010～2020で最後にAllになる。なぜかよく分からない。

legend("topright", legend=levels(dat0$year), col=cols, pch=15, bg="white", cex=1.5, y.intersp=1.2)

で凡例書こうとするとトラブるToLOVEる。こっちはちゃんとAllから始まって2010～2020の順に表記され、黒から始まるので凡例と生存曲線とで色が一致しなくなる。なぜかよくわからない。

この判例周りの扱いを見るともしかしたらバグもあるのかもしれないと思ったけど、よくわからない。

今回は気付きについてのメモまで。個人的な印象としては、3.6.3までのAllを追加したときにもyearで因子型が継承される挙動のほうがよくわからない。