Help us understand the problem. What is going on with this article?

ggplot2のグラフ表現を改良する

奈良文化財研究所で9月10日〜11日に開催された「考古学・文化財のためのデータサイエンス研究集会」で発表した資料考古学のためのデータビジュアライゼーション(ggplot2を使用したデータ可視化)について、Ben Marwick先生から指導をいただきましたので、備忘録として公開します。

alphaはaes()の外に出す

下記のコードではalpha=0.7に設定していますが、alphaの指定をaes()の中で行っているので、無意味な凡例がついてしまっています。
あと、pirnt(p)は明示しなくても良いのでは?という提案もいただきました。

p<-pot%>%
    ggplot(aes(x=口径,fill=分類,alpha=0.7))+
        geom_density()+
        scale_fill_ptol()+
        theme_minimal()
print(p)

density_befor.png

修正したのがこれ。alphaaes()の中に入れるのではなく、geom_density()の引数として指定しました。コードもスッキリしています。

pot%>%
    ggplot(aes(x=口径,fill=分類))+
        geom_density(alpha=0.7)+
        scale_fill_ptol()+
        theme_minimal()

density_after.png

箱ひげ図を改良する

次の指摘はこちらの箱ひげ図。Marwick先生のご指摘は2点。

  • ggforceパッケージのgeom_sina()を使って箱ひげ図とドットを重ねて表示する。
  • カラーパレットは viridis colour palettesを使ったほうが良い。

ggforceパッケージは初めて知りました。viridis colourは連続量のイメージがありますが、離散的なデータにも対応しているようです。

pot%>%
    ggplot(aes(x=分類,y=口径,fill=分類))+
        geom_boxplot()+
        scale_fill_ptol()+
        coord_flip()+
        theme_minimal()

boxplot_bifore.png

ということで改良してみました。
viridisカラーパレットを使用する理由として、単に見やすいということだけではなく、ユニバーサルなカラーパレットであるということが挙げられます。brewerやptolも優れたカラーパレットですが、viridisカラーが現時点では一日の長があるようです。

pot%>%
    ggplot(aes(x=分類,y=口径,fill=分類))+ 
        geom_boxplot(alpha = 0.2)+      #不透明度を0.2
        geom_sina(aes(colour = 分類),     #geom_sina()関数でaes()の引数にcolour=分類を指定
              alpha = 0.4, 
              size = 3) +
        scale_fill_viridis_d() +        #viridis_d()は離散量、連続量ならviridis_c()を指定する
        scale_colour_viridis_d() +
        coord_flip()+
        theme_minimal()

boxplot_after.png

統計データは小数点2〜3桁で表示すれば十分

分散分析の結果を表示した表に対してのコメントです。こうした統計表の表示には小数点2〜3桁で十分とのこと。

# aov関数の結果をTukeyHSD関数に渡す
tkh<-aov(口径~分類,data=pot)%>%TukeyHSD()
tkh$分類%>%kable(format="markdown")
diff lwr upr p adj
B型-A型 -6.58 -8.1885528 -4.971447 0.0000000
C型-A型 -4.54 -6.1485528 -2.931447 0.0000000
C型-B型 2.04 0.4314472 3.648553 0.0087802

改良したのがこちら。小数点3桁でまるめられて良い感じ。ただ、もともとあった変数の組み合わせ("B型-A型"の部分)が消えてしまっているのが少し残念です。

tkh <- 
    aov(口径 ~ 分類, data = pot) %>% 
    TukeyHSD() %>% 
    .$分類 %>%    #TukeyHSD関数の結果から$分類を選択 
    as_tibble() %>%     #tibble_df形式に変換
    mutate_if(is.numeric, round,3)      #mutate_if()でnumericクラスのカラムにround関数を適用する。
tkh%>%kable(format="markdown")
diff lwr upr p adj
-6.58 -8.189 -4.971 0.000
-4.54 -6.149 -2.931 0.000
2.04 0.431 3.649 0.009

集計表タイプのデータはgeom_barのかわりにgeom_colを使う

geom_bar(stat="identity")を指定していますが、geom_col(position="fill")とするべき、との提案です。結果はどちらも同じですが、geom_colstat="identity"がデフォルトとなっていますので、集計表タイプのデータ描画に利用しやすくなっています。

befor

toj%>%
    ggplot(aes(x=遺跡名,y=点数,fill=器種))+
        geom_bar(stat="identity",position="fill")+
        coord_flip()+
        scale_fill_ptol()+
        theme_minimal()

after

toj%>%
    ggplot(aes(x=遺跡名,y=点数,fill=器種))+
        geom_bar(stat="identity",position="fill")+
        coord_flip()+
        scale_fill_ptol()+
        theme_minimal()

bar_before.png

tidyverseだけでデータ整形

細別器種を大別器種に再分類していますが、コードが冗長です。Marwick先生は次のように指摘します。

  • grepl()を多用しているが、もっとシンプルに記述できる。
  • tidyverseで全てを処理したいならglepl()ではなくstringr::str_detect()を使っても良い。
toj2<-toj%>%
    mutate(大別器種 = case_when(
        grepl("碗",器種)|grepl("皿",器種)|grepl("その他の食膳具",器種)    == TRUE ~ "食膳具",
        grepl("壺・甕・瓶",器種) == TRUE ~ "貯蔵具",
        grepl("灯明皿・油注",器種)|grepl("その他",器種)|grepl("すり鉢",器種)|
            grepl("鍋",器種)|grepl("土瓶",器種) == TRUE ~ "その他",
        ))

次のように修正します。シンプル。

toj2<-toj%>%
    mutate(
        大別器種 = case_when(
            str_detect(器種,"碗 | 皿 | その他食膳具") ~ "食膳具",
            str_detect(器種,"壺・甕・瓶") ~ "貯蔵具",
            str_detect(器種,"灯明皿・油注 | その他 | すり鉢 | 鍋 | 土瓶") ~ "その他",
        )
    )

大別器種というカラムが作られます。

遺跡名 器種 点数 大別器種
弁天貝塚 134 食膳具
弁天貝塚 84 食膳具
弁天貝塚 その他食膳具 34 食膳具
弁天貝塚 土瓶 6 その他
弁天貝塚 0 その他

散布図行列に有意性の評価を入れる

GGally::ggpairs()は見栄えのよい表現ですが、有意性の評価を図中に入れるためPerformanceAnalytics::chart.Correlation()を使用してはどうか、という提案をいただきました。

GGally::ggpairs()による描画
pairs_before.png

PerformanceAnalytics::chart.Correlation()による描画です。文字サイズの調整のやり方などがうまくできていませんが、有意性の評価が星印で示されています。

library(PerformanceAnalytics)
iron %>%
    select(全長, 刀身長, 茎長, 刀身先幅, 刀身元幅, 刀身元厚, 茎先幅) %>%
        chart.Correlation(histogram = TRUE, pch = 19)

pairs_after.png

図中に数式を描画する

研修では以下のようにlm()を適用して必要なパラメーターを算出して手動で回帰式を表示していました。

icoe<-lm(刀身長 ~ 刀身元幅,data=iron)%>%summary()
icoe$coefficients%>%kable()
Estimate Std. Error t value Pr(>|t|)
(Intercept) -6.280888 1.9780720 -3.175257 0.0022892
刀身元幅 10.723991 0.7889999 13.591878 0.0000000

y=10.72399X-6.280888

ggpmisc::stat_poly_eq()ggpmisc::stat_fit_glance()を利用して回帰式、相関係数、p値をプロットに記述する手法を教えていただきました。
複雑な記述になるので、引数やパラメーターの意味が理解できていませんが、再現性を保つためにはこの手法が欠かせないところです。

library(ggpmisc)
iron %>%
    ggplot(aes(x=刀身元幅,y=刀身長))+
        geom_point()+
        geom_smooth(method="lm")+
        theme_minimal() +
        stat_poly_eq(formula = y ~ x,
            eq.with.lhs = "italic(hat(y))~`=`~",
            aes(label = paste(stat(eq.label), 
                stat(rr.label), 
                sep = "~~~")
            ), parse = TRUE
        ) +
        stat_fit_glance(label.y = 0.9,
            method = "lm",
            method.args = list(formula = y ~ x),
            aes(label = sprintf(
                '~~italic(P)~"="~~%.25f',
                stat(p.value)
                )
            ),parse = TRUE
        )

lm_after.png

Why do not you register as a user and use Qiita more conveniently?
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away
Comments
Sign up for free and join this conversation.
If you already have a Qiita account
Why do not you register as a user and use Qiita more conveniently?
You need to log in to use this function. Qiita can be used more conveniently after logging in.
You seem to be reading articles frequently this month. Qiita can be used more conveniently after logging in.
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away