summaryの罠
年月を6桁の数値(YYYYMM)で表すために以下のように書いて何気なくsummaryを実行したところ、思わぬ挙動となった。
> Year <- rep(2012:2016, each=12) > Month <- rep(1:12, 5) > YM <- Year * 100 + Month > summary(YM) Min. 1st Qu. Median Mean 3rd Qu. Max. 201200 201300 201400 201400 201500 201600
上記の通り、月は1〜12の範囲なのにいずれも下2桁が00となっている。
念のため他の関数を使ってみると、ちゃんとした値が返る。
> fivenum(YM) [1] 201201.0 201303.5 201406.5 201509.5 201612.0 > min(YM) [1] 201201 > max(YM) [1] 201612
summaryはgeneric functionなので、今回の例のようなnumeric型に対するsummaryの挙動を確認するためにmethodを調べてみると:
> methods("summary") [1] summary.aov summary.aovlist* summary.aspell* [4] summary.check_packages_in_dir* summary.connection summary.data.frame [7] summary.Date summary.default summary.ecdf* [10] summary.factor summary.glm summary.infl* [13] summary.lm summary.loess* summary.manova [16] summary.matrix summary.mlm* summary.nls* [19] summary.packageStatus* summary.PDF_Dictionary* summary.PDF_Stream* [22] summary.POSIXct summary.POSIXlt summary.ppr* [25] summary.prcomp* summary.princomp* summary.proc_time [28] summary.srcfile summary.srcref summary.stepfun [31] summary.stl* summary.table summary.tukeysmooth*
色々と出てきたが、どうやらsummary.defaultがそれっぽい。
summary.defaultの中身を調べてみる。
> summary.default … 中略 … }) else if (is.numeric(object)) { nas <- is.na(object) object <- object[!nas] qq <- stats::quantile(object) qq <- signif(c(qq[1L:3L], mean(object), qq[4L:5L]), digits) names(qq) <- c("Min.", "1st Qu.", "Median", "Mean", "3rd Qu.", "Max.") if (any(nas)) c(qq, `NA's` = sum(nas)) else qq } … 中略 …
どうやらsummary.defaultの中身は{stats}のquantileで、間にmeanを挟んでいるだけの様子。そこでquantileを直接実行してみると:
> quantile(YM) 0% 25% 50% 75% 100% 201201.0 201303.8 201406.5 201509.2 201612.0
それっぽい値が返ってきている。どうやら次の処理のsignifの引数のdigitsが怪しい。
そこで今度はsummaryのdigitsを指定して実行してみる。
> summary(YM, digits = 8) Min. 1st Qu. Median Mean 3rd Qu. Max. 201201.0 201303.8 201406.5 201406.5 201509.2 201612.0
狙った通りの挙動となった。どうやら6桁の整数では大きすぎて小数点以下が丸められてしまっていたようだ。
めでたしめでたし・・・・と思ったら、よく見るとsummary(つまりquantile)の返り値とfivenumの返り値が微妙に違っている。
> summary(YM, digits = 8) Min. 1st Qu. Median Mean 3rd Qu. Max. 201201.0 201303.8 201406.5 201406.5 201509.2 201612.0 > fivenum(YM) [1] 201201.0 201303.5 201406.5 201509.5 201612.0
25パーセンタイルと75パーセンタイルの値が違う。調べてみると、fivenumは四分位数ではなくヒンジを返すのだそう。
http://cse.naro.affrc.go.jp/takezawa/r-tips/r/59.html
なのでデータが偶数個の場合に先ほどのような挙動となる様子。
こんな感じ。
> hoge <- 1:10 > summary(hoge) Min. 1st Qu. Median Mean 3rd Qu. Max. 1.00 3.25 5.50 5.50 7.75 10.00 > quantile(hoge) 0% 25% 50% 75% 100% 1.00 3.25 5.50 7.75 10.00 > fivenum(hoge) [1] 1.0 3.0 5.5 8.0 10.0 # 合わない!
> foo <- c(0, hoge) > summary(foo) Min. 1st Qu. Median Mean 3rd Qu. Max. 0.0 2.5 5.0 5.0 7.5 10.0 > quantile(foo) 0% 25% 50% 75% 100% 0.0 2.5 5.0 7.5 10.0 > fivenum(foo) [1] 0.0 2.5 5.0 7.5 10.0 # 一致した!