🍭{ggplot2} 1.0.1.9003以降のバージョンで回帰線を引く時などに使うgeom_smooth()関数の挙動が変わるカモ

以前書いた記事に関して、気になるコメントをもらって、ちょっともやもやしていたので検証してみました。 {ggplot2}のgeom_smooth()についての挙動です。

最後のグラフで、glmの結果をggplotするときのfamilyの指定の所ですが、自分の環境ではmethod.argではなく geom_smooth(method = "glm", family = "Gamma"(link = "log"), se = FALSE) のように直指定しないとfamilyが有効になりませんでした。

というコメントに対して、まあバージョンの違いでしょう、みたいに雑に返答してしまいましたが、具体的にどう変化したのかがきちんと説明できませんでした。のでちょっと正確さに不安がある（し、{ggplot2}そのものの仕様もまた変わるかもしれない）のですが、メモとして残します。間違いがございましたらご指摘ください。助かります。

結論から言うと、やはりバージョンの問題っぽい、です。

問題になっているのは、geom_smooth()です。散布図を描いた際の直線回帰や加法モデル、一般化線形モデルなど、いくつかの統計モデリングによる推定値の当てはめに利用できる便利な関数です。

{ggplot2} 1.0.1でのgeom_smooth()

便宜上、旧geom_smoothと呼びます。

# 検証のためにCRANから最新(2015-12-06現在)のバージョンをインストール
install.packages("ggplot2")

packageVersion("ggplot2")

## [1] '1.0.1'

1.0.1 です。問題になっているgeom_smooth()について見てみましょう。

geom_smooth %>% str()

## function (mapping = NULL, data = NULL, stat = "smooth", position = "identity", ...)

1.0.1のgeom_smooth()では、一般化線形モデルによるフィッティングを行う時には、id:matsuken92 さんがおっしゃるように、

mpg %>% ggplot(aes(displ, hwy)) +
  geom_point() +
  geom_smooth(method = "glm", family = "Gamma"(link = "log"), se = FALSE)

とすれば回帰線が引かれます。method引数で回帰モデルのタイプを指定し（lmやらglmやら、gamやらloessやら...）、必要であれば（GLMなので必要）family引数を指定する、という形になりました。

f:id:u_ribo:20151206070249p:plain

{ggplot2} 1.0.1.9003でのgeom_smooth()

上記の記事を書いた自分がインストールしていたのはGitHubのバージョンだったので、以下のような仕様になっていました。

# GitHubの開発版をインストール（記事を書いた際のバージョン）
devtools::install_github("hadley/ggplot2")
packageVersion("ggplot2")

1.0.1.9003... CRANよりも新しい（開発版）です。このバージョンでのgeom_smooth()は

geom_smooth %>% str()

## function (mapping = NULL, data = NULL, stat = "smooth", method = "auto", 
##     formula = y ~ x, se = TRUE, position = "identity", na.rm = FALSE, 
##     show.legend = NA, inherit.aes = TRUE, ...)

のようになっています。いくつか引数が追加されていますね。ここでは件のmethod.argsは見られませんが、

stat_smooth %>% str()

## function (mapping = NULL, data = NULL, geom = "smooth", position = "identity", 
##     method = "auto", formula = y ~ x, se = TRUE, n = 80, span = 0.75, 
##     fullrange = FALSE, level = 0.95, method.args = list(), na.rm = FALSE, 
##     show.legend = NA, inherit.aes = TRUE, ...)

とすると現れます。また、このgeom_smooth()の違いはヘルプドキュメントのExampleを見ると明らかです。1.0.1.9003以降でのgeom_smooth()では、説明変数に対して平滑化関数などを嚙ますことを許すようですが、GLMのように目的変数に対して特定の分布を直接適用することができない、という挙動が旧geom_smoothとの変化点みたいです（回帰線を引く以外のgeom_smooth()の挙動に合わせたのかな、と推測）。なので、id:matsuken92 さんが試されたものと自分の書いたコードでの動作が異なる、という話、ではないでしょうか...

mpg %>% ggplot(aes(displ, hwy)) +
  geom_point() +
  geom_smooth(method = "glm",
              method.args = list(family = "Gamma"(link = "log")),
              se = FALSE)

元の記事同様、開発版の{ggplot2}では、method.argsでfamilyを指定して、旧geom_smooth()と同じ図を描きます（省略）。

ちなみに、旧geom_smooth()のような処理をすると

mpg %>% ggplot(aes(displ, hwy)) +
  geom_point() +
  geom_smooth(method = "glm", family = "Gamma"(link = "log"), se = FALSE)

## Error: Unknown parameters: family

となって怒られます。

もう少し詳しく見ていきます。

mpg %>% ggplot(aes(displ, hwy)) +
  geom_point() +
  geom_smooth(method = "lm")

geom_smooth()の何の引数も与えない場合の回帰線は一般化加法モデルの回帰線です（初期値）。これは説明変数に対して平滑化関数s()を与えています(formual = y ~ s(x))。これを変更して、method = "lm"とすれば直線回帰モデルの回帰線が引かれます(formual = y ~ x)。では、このmethodを"glm"や"gam"にすれば任意の回帰モデルを適用してくれるかというとそうではなくて、モデル式の情報を別途method.argsで与える必要があります。

method.argsに与えるパラメータ

method.args引数には、リスト形式でモデル式のパラメータを与えます。そのパラメータは、適用したい統計モデルの関数と一致します。一般化線形モデルであれば、glm()なので

args(glm) %>% as.list() %>% names()

##  [1] "formula"   "family"    "data"      "weights"   "subset"   
##  [6] "na.action" "start"     "etastart"  "mustart"   "offset"   
## [11] "control"   "model"     "method"    "x"         "y"        
## [16] "contrasts" "..."       ""

これらと同じパラメータを指定します（使うのはfamily）。Rでの回帰モデルを実行する関数は共通の引数名を持っていることが多いので、こうした統一は良いなと思いました。

geom_smooth()をラップする

method.argsに与える値はモデルによりけりなので、ヘルプを見ると、どうも、適用したいモデルごとにgeom_smooth()をラップした関数を作っておくのが良いようです。例えば、ロジスティック分布に従う一般化線形モデルによる回帰を行う際には

# http://www.ats.ucla.edu/stat/r/dae/logit.htm から引用
mydata <- read.csv("http://www.ats.ucla.edu/stat/data/binary.csv")
mydata %>% glm(formula = admit ~ gre,
               family = "binomial", 
               data = .)

のような引数を与えるので、これをプロットする際にはfamilyをmethod.argsに渡します（formulaはgeom_smooth()内で指定するか、省略）。その際、関数としてラップすることで異なるプロットにも適用できるのでコードを省略することが可能になります。

# ロジスティック回帰用のgeom_smooth()ラッパー関数を作っておく
binomial_smooth <- function(...) {
  geom_smooth(method = "glm", 
              method.args = list(family = "binomial"), ...)
}

# 先ほどのgeom_smooth()ラッパー関数を使用する
mydata %>% ggplot(aes(gre, admit)) + 
  geom_point() + 
  binomial_smooth()