カテゴリ:統計 > R

第4章は統計を扱います。


今回「シロート統計学」のハルさんとコラボレーションすることとなりました。


ハルさん、ありがとうございます!


シロート統計学はEZRを使った統計分析をわかりやすく解説されています。




第4章はシロート統計学で使われていたEZRをRで行うとどうなるのか?といった視点で進めていきます。


今回使うデータもハルさんのサイトと同じものを使わせでいただく事になりました。それぞれ見比べることで参考にしてみてください!


今回はt検定を紹介します



まずt検定についてはハルさんのサイトをご参照ください。

 



1.準備

第4章は毎回ExcelデータをダウンロードしてRを使うのでプロジェクトで管理して行うことを勧めています。



ここではR練習というプロジェクトを作り、Excelファイルを入れるためのdataフォルダを作っています。
これを前提に次から進めていきます。


2.スクリプトファイルの作成

次にRのコードを書くためのスクリプトファイルを作ります。

スクリーンショット 2019-10-25 12.15.42

完成です。
スクリーンショット 2019-10-25 12.18.51


3.データのダウンロード

今回はハルさんのサイトのデータを使わせていただきます。

デモデータ(t検定)

これをダウンロードしてdataフォルダに入れればいいのですが実はRでできてしまいます
download.file関数を使います。" "を忘れないようにしてください。

url <- "https://haru-reha.com/wp-content/uploads/2018/03/demo-t-test.xlsx"
destfile = "data/demo-t-test.xlsx"

download.file(url, destfile)

以下説明します。

download.file(url = “ファイルのURL”,
        destfile = “保存したい場所/ファイル名”)


urlはデモデータで右クリック → リンクのアドレスをコピー

destfileは保存場所と保存のファイル名を指定します。
保存場所は今回プロジェクトを使っているのでR練習フォルダになります。加えてdata/を付け足すことでR練習フォルダ内にあるdataフォルダという意味になります。
ファイル名は自由に決められますが今回は元のファイルと同じにしました。拡張子も忘れないようにしましょう。

もしプロジェクトを使っていなければ保存場所はgetwd関数で出てきたフォルダになります。
getwd関数の()には何も入れません。

getwd()


この方法を使う最大のメリットは、次回使うExcelデータはurlの部分を変えるだけでできてしまうことです。毎回右クリックでアドレスを保存 → 保存したファイルを指定したところに移動させて・・・といった手作業必要ありません。こういった作業もRで行っていくことでRにも早く慣れてくると思います。



4.データの読み込み

データを読み込みます。
今回は【4-0】第4章を進めていく上での準備で行った方法で進めます。

View Fileでデータを確認します。
スクリーンショット 2019-10-25 14.07.12

今回は握力のデータです。A群とB群を比較します。
ただA1にデータが入っていません。実際にはB2からC62までデータが入っていることを確認します。
スクリーンショット 2019-10-25 14.11.22

次はImport Datasetを選びます。
スクリーンショット 2019-10-25 14.07.18

ポイントは2つです。
①データの名前(変数名)を付ける
何でもいいのですが今回はハルさんのサイトと同じgripにしました。

②読み込む範囲を指定する
今回A1からのデータではないので先程確認したB2からC62を指定します。
B2:C62のように左上と右下を:でつなげます
スクリーンショット 2019-10-25 14.17.12


そして右下にコードが自動的に作られます。Importを押せば完了なのですが、このコードをコピーしスクリプトに貼り付けておけば1年たった後でも同じことができます。EZRでもスクリプトを保存することができないわけではないのですが、再現性(後でしても、他の人がしても同じ事ができる)を保つためにもこういったコードを残しておく習慣をつけるようにしましょう。
スクリーンショット 2019-10-25 14.19.26


コードの一番下にあるView関数を使うことでRStudio内でもデータの確認ができます。このタブを消してもデータに影響はありません。もう1回View関数を使えばまた表示できます。
スクリーンショット 2019-10-25 14.28.00

View関数はEZRで言う表示と同じです(下図はハルさんのサイトより。比べてみてください)



これでデータの取り込みは完了です!



5.データの要約

ハルさんは次にデータの要約をしています。
EZRでのデータの要約と全く同じ機能はないですが、第2章で紹介したtidyverseパッケージのgroup_by関数とsummarize関数が使えます。group_by関数とsummarize関数に関してはこちらで紹介しています。

%>%やgroup_by関数、summarize関数はtidyverseパッケージに含まれていますのでtidyverseパッケージを呼び出します。もしtidyverseパッケージを全く使ったことが無い方はパッケージをインストールします。1度でも使ったことがあれば次の1行は必要ありません。

install.packages("readxl")

実際のコードは以下になります。イメージ図も添付します。

library(tidyverse)
grip %>% group_by(category) %>% summarize(平均 = mean(grip), 標準偏差 = sd(grip), '0%' = quantile(grip, 0), '25%' = quantile(grip, 0.25), '50%' = quantile(grip, 0.5), '75%' = quantile(grip, 0.75), '100%' = quantile(grip, 1), n = n())

スクリーンショット 2019-10-25 20.52.13

スクリーンショット 2019-10-26 7.55.30


これでA群とB群のデータのばらつきを確認することができます。
ちなみにこのコードをコピーして色がついた箇所を変更すれば他の場面でも使えます!

このままでもいいのですが、データ要約は後でグラフ作成に使うのでgrip_summaryという名前をつけます。

grip_summary <- 
grip %>% group_by(category) %>% summarize(平均 = mean(grip), 標準偏差 = sd(grip), '0%' = quantile(grip, 0), '25%' = quantile(grip, 0.25), '50%' = quantile(grip, 0.5), '75%' = quantile(grip, 0.75), '100%' = quantile(grip, 1), n = n())


6.正規性の確認

次に正規性の確認を行います。ハルさんのサイトではヒストグラムを作成しました。
ヒストグラムの作り方はこちらで紹介しています。


ハルさんのヒストグラムは棒が5本だったので同じ形にするようbins = 5とします。
A群とB群の棒を横に並べるときはposition = "dodge"を使います。
ggplot()で行を変える時は %>% ではなく + を使うので注意してください。

ggplot(data = grip) + 
  geom_histogram(aes(x = grip, fill = category), position = "dodge", bins = 5) 

スクリーンショット 2019-10-25 20.55.06



加えてシャピロ・ウィルク検定も紹介されています。シャピロ・ウィルク検定はshapriro.test関数を使います。
ただA群とB群それぞれで行いますのでgrip$gripの列をA群とB群に分ける必要があります。

shapiro.test(grip$grip[category == "A"])
shapiro.test(grip$grip[category == "B"])



データの中で特定の条件だけを抜き出すには[ ]を使います。



ハルさんのサイトのEZRで行う場合と見比べてみてください。
category == "A" と書いてあるところの意味が見えてきます。
ちなみに変数(1つ選択)のgripはgrip$gripの色のついた部分です。

スクリーンショット 2019-10-25 21.13.57


結果は以下のとおりです。
スクリーンショット 2019-10-25 19.36.44

どちらも0.5を超えているので正規分布であるという仮説は棄却されませんでした。
このあたりの解釈はハルさんのサイトをご参照ください。


7.1度に計算するsplit + map関数

shapiro.test(grip$grip[category == "A"])
shapiro.test(grip$grip[category == "B"])

上記のように1つずつ計算する方法もいいのですが群の数だけ繰り返します。プログラミングであるRは繰り返しに強いという特徴があります。


群ごとにデータを分割し、まとめて計算する方法として今回はsplit関数map関数を使います。


grip %>% 
  split(.$category) %>% 
  map(~shapiro.test(.$grip))

スクリーンショット 2019-10-26 22.42.55

まずgripのデータを
split関数を使ってA群とB群の2つのデータに分割し、map(~シャピロウィルク検定)でシャピロウィルク検定を繰り返します。

スクリーンショット 2019-10-26 23.01.43

EZRだと1つ1つ検定を繰り返す必要があります。
Rを活用すると1度に作業が終わるため、繰り返すことでのやり間違えの予防や時間の短縮にも繋がります。


8.等分散性を確認する

先程は正規性を調べましたが、今回は等分散性を調べます。
F検定と言いますがRではvar.test関数を使います。



Rで行うとこうなります。
var.test(目的変数 
~ グループ)

var.test(grip$grip ~ grip$category)

もしくはdata = を使えば$が要らなくなります。

var.test(grip ~ category, data = grip)
スクリーンショット 2019-10-27 5.57.05

p値 > 0.05だったので等分散性は棄却されませんでした。



9.t検定を行う

いままで正規分布かどうか、等分散性かどうかを確認しました。

スクリーンショット 2019-10-27 6.17.29

これでt検定が使えます。t検定はt.test関数を使います。
t検定は2種類の書き方があります。その違いはデータの形にあります。
ハルさんのサイトでEZRでは左の形しかできませんでしたがRならどちらも可能です
ちなみに左をlongデータ(縦に長い)、右をwideデータ(横に広い)といいます。
スクリーンショット 2019-10-27 21.34.00

スクリーンショット 2019-10-27 8.09.41

今回はlongデータなのでこのようになります。

t.test(grip$grip ~ grip$category, var = TRUE)

data = gripを使うとgrip$は外せます。どちらでも好きな方で大丈夫です。

t.test(grip ~ category, var = TRUE, data = grip)

スクリーンショット 2019-10-27 21.42.47
EZRと同じ結果になりました!


10.グラフの作成
EZRでは丁寧にグラフも作ってくれますがRでは自作する必要があります。
EZRでは平均と標準偏差の棒グラフを作成していました。
先程のgrip_summaryに平均と標準偏差を求めていたのでそれを使います。

棒グラフの作り方はここで紹介しています。


①最低限のグラフ(見た目は気にしない)

グラフを作る上で最低限必要な要素としては
棒グラフ:geom_bar(aes(x軸の指定、y軸の指定), stat = "identity")
エラーバー:geom_errobar(aes(x軸の指定, エラーバーの下端, エラーバーの上端))
y軸の名前が「平均」になるので「grip」に変更する

ggplot(data = grip_summary)+
  geom_bar(aes(x = category, y = 平均), stat = "identity") +
  geom_errorbar(aes(x = category, ymin = 平均 - 標準偏差, ymax = 平均 + 標準偏差)) + 
  labs(y = "grip") 
スクリーンショット 2019-10-27 22.59.56
ただこれでは見るに耐えません・・・



②見栄えを変更
次は以下を修正します。
棒グラフ:周りの線をblack、中の色をgray
エラーバー:幅を補足する(今回は0.1)

ggplot(data = grip_summary)+
  geom_bar(aes(x = category, y = 平均), color = "black", fill = "gray", stat = "identity") +
  geom_errorbar(aes(x = category, ymin = 平均 - 標準偏差, ymax = 平均 + 標準偏差), width = 0.1) + 
  labs(y = "grip")
スクリーンショット 2019-10-27 23.00.05
だいぶ見た目が良くなりました。



③EZRのグラフにできるだけ近づける
EZRのグラフに似せるためにもうひと工夫します。
背景を白にする(グラフのテーマをclassicに変える)
文字を大きくする(今回はsizeを15に変更)
棒グラフが浮いているように見えるのを修正する

ggplot(data = grip_summary)+
  geom_bar(aes(x = category, y = 平均), color = "black", fill = "gray", stat = "identity") +
  geom_errorbar(aes(x = category, ymin = 平均 - 標準偏差, ymax = 平均 + 標準偏差), width = 0.1) + 
  labs(y = "grip") +
  theme_classic(base_size = 15) +
  scale_y_continuous(expand = c(0,0), limits = c(0,50))
スクリーンショット 2019-10-27 23.00.12
これでEZRにだいぶ似ました。



11.まとめ

かなり長くなりましたが今回はt検定を紹介しました。

最初だったので1つ1つ説明しました。ボリュームが多すぎて慣れないところもあると思いますが、これから何度も出てくるものもありますので少しずつ慣れて貰えればと思います。


12.コードをまとめました
今回のコードを1つにまとめました。
パッケージはlibrary()を1度だけ行えばいいので最初に出しています。

#パッケージの読み込み
library(readxl)
library(tidyverse)

#データのダウンロード
url <- "https://haru-reha.com/wp-content/uploads/2018/03/demo-t-test.xlsx"
destfile = "data/demo-t-test.xlsx"

download.file(url, destfile)

#データの読み込み
grip <- read_excel("data/demo-t-test.xlsx", 
                   range = "B2:C62")
View(grip)  

#データの要約
grip_summary <- 
  grip %>% 
  group_by(category) %>% 
  summarize(平均 = mean(grip),
              標準偏差 = sd(grip),
              '0%' = quantile(grip, 0),
              '25%' = quantile(grip, 0.25),
              '50%' = quantile(grip, 0.5),
              '75%' = quantile(grip, 0.75),
              '100%' = quantile(grip, 1),
              n = n())
grip_summary

#正規性を調べる
#ヒストグラム
ggplot(data = grip) + 
  geom_histogram(aes(x = grip, fill = category), position = "dodge", bins = 5) 

#シャピロ・ウィルク検定
shapiro.test(grip$grip[category == "A"])
shapiro.test(grip$grip[category == "B"])

grip %>% 
  split(.$category) %>% 
  map(~shapiro.test(.$grip))


#等分散性を確認
var.test(grip ~ category, data = grip)

#t検定
t.test(grip$grip ~ grip$category, var = TRUE)

t.test(grip ~ category, var = TRUE, data = grip)

#グラフ化
ggplot(data = grip_summary)+
  geom_bar(aes(x = category, y = 平均), stat = "identity") +
  geom_errorbar(aes(x = category, ymin = 平均 - 標準偏差, ymax = 平均 + 標準偏差)) + 
  labs(y = "grip") 

ggplot(data = grip_summary)+
  geom_bar(aes(x = category, y = 平均), color = "black", fill = "gray", stat = "identity") +
  geom_errorbar(aes(x = category, ymin = 平均 - 標準偏差, ymax = 平均 + 標準偏差), width = 0.1) + 
  labs(y = "grip")

ggplot(data = grip_summary)+
  geom_bar(aes(x = category, y = 平均), color = "black", fill = "gray", stat = "identity") +
  geom_errorbar(aes(x = category, ymin = 平均 - 標準偏差, ymax = 平均 + 標準偏差), width = 0.1) + 
  labs(y = "grip") +
  theme_classic(base_size = 15) +
  scale_y_continuous(expand = c(0,0), limits = c(0,50))

                


第4章は統計を扱います。

今回「シロート統計学」のハルさんとコラボレーションすることとなりました。

ハルさん、ありがとうございます!

シロート統計学はEZRを使った統計分析をわかりやすく解説されています。



第4章はシロート統計学で使われていたEZRをRで行うとどうなるのか?といった視点で進めていきます。

今回使うデータもハルさんのサイトと同じものを使わせでいただく事になりました。それぞれ見比べることで参考にしてみてください!



今回は4章を進めるにあたって抑えておきたいRの基本的な操作について復習します。

1.プロジェクトを作る

第4章では色々なデータを扱うため、フォルダやRのファイルが煩雑になりやすいです。そこで今回はプロジェクトを使って管理します。

RStudioでプロジェクトを作成する方法は第1章で紹介しています。




今回は「R練習」というプロジェクトを作りました。スクリーンショット 2019-10-20 7.54.14


更にデータを入れるフォルダを作ります。分析に使うExcelファイルはデスクトップやマイドキュメントに置かずにここに置くようにします!
スクリーンショット 2019-10-20 23.07.53

そしてRのコードを書くスクリプトファイルを作ります。
今回は【4-1】で使うためのスクリプトファイルを作成しました。
スクリーンショット 2019-10-20 23.12.27

これで完成です!
スクリーンショット 2019-10-20 7.57.37




2.Excelのデータを読み込む

Excelを読み込む方法はいくつかありますが、ここではコードを書かない方法を紹介します。


まずデータが入ってるdataフォルダに進みRStudioのfilesタブからImport DatasetでExcelファイルをRで読み込みます。


スクリーンショット 2019-10-23 8.08.00

View Fileではエクセルファイルを直接開くことができます
(エクスプローラやFinderを使わなくていいので便利!)

ここではImport Datasetを選びます。

スクリーンショット 2019-10-23 8.11.09

Name:Rで使うときのデータの名前(自由に決められる)
sheet:Excelに複数タブあればどれを使うか指定できる
Range:データの読み込む範囲を指定できる

主には上記3つで読み込む範囲を指定します。
もし1枚目のタブでデータの読み込む範囲の左上がA1だったらNameの指定だけでも大丈夫です。

そしてこの方法の良いところは右にコードを自動で書いてくれることです。これを忘れずにスクリプトにコピーします。

そうすることで2回目以降は1秒もかからずに行えるのでかなりの時短になります!




まとめ

今後色々な検定やデモデータを使うことになります。EZRを使うとあまり意識されない部分ですが、プロジェクトファイルを作ることでデータの管理やスクリプトファイルの管理が進めやすくなります。

次回はt検定を紹介します。





第3章ではggplot2を使ったグラフの作り方について説明してきました。

【3-1】ExcelにはないRでグラフを作るメリットと特徴

【3-2】ggplot2でグラフを作る流れを説明します

【3-3】Rのggplot2で散布図を作るgeom_point関数

【3-4】Rのggplot2でヒストグラムを作るgeom_histogram関数

【3-5】Rのggplot2で密度曲線を作るgeom_density関数

【3-6】Rのggplot2で箱ひげ図を作るgeom_boxplot関数

【3-7】棒グラフの基本とRのggplot2で棒グラフを作るgeom_bar関数

【3-8】ggplot2で折れ線グラフを作るgeom_line関数

【3-9】ggplot2でヒートマップを作るgeom_tile関数

【3-10】Rのggplot2でグループ毎にグラフを作りまとめて表示するfacet_wrap関数

【3-11】Rのggplot2で作った複数のグラフを1つにまとめるgridExtraパッケージ


今まで基本的なグラフの作り方について紹介しましたが、ggplot2は他にも様々なな調整やができます。


今回はグラフを作るときに参考になるサイトを紹介します。


質・量共にこのサイトより充実しているものばかりですが、ある程度慣なれていないと読み解けない箇所もあります。(初心者に向けた)説明の細かさはこのサイトが1番だと思いますので、困ったらこのサイトで確認してください。



チートシート
チートシートとは直訳するとカンニングペーパーという意味で、スライド1〜2枚でその機能を確認できる便利シートになります。基本的な使い方でこまればまずはチートシートを確認しています。

Rにもいろいろなチートシートがあり、ggplot2のチートシートもあります。
ggplot2のチートシートはRStudioのヘルプ→Cheetsheetsにあります(英語版)。

スクリーンショット 2019-10-08 21.44.51

また日本語版もありますのでダウンロードしてお使いください。



書籍

RユーザのためのRStudio[実践]入門
−tidyverseによるモダンな分析フローの世界

松村 優哉 (著), 湯谷 啓明 (著), 紀ノ定 保礼 (著), 前田 和寛 (著)

通称「宇宙本」
この本は初めてRを勉強する方に最初に進めている本です。
Rの基本的な操作から第2章で行ったデータハンドリング、第3章で行ったggplot2によるグラフ作成も紹介されています。他にもウェブスクレイピングやレポートの作成などRを使ってデータ分析を行いたい方にはピッタリはまると思います。


Rグラフィックスクックブック ―ggplot2によるグラフ作成のレシピ集

この本ではggplot2の細かい使い方が紹介されている辞書的な1冊です。
ただ2013年の本なのでRやggplotのバージョンの問題で上手く動かないものがあるというコメントもあります。
実はこの本の第2版がオンラインで無料で読むことができます(英語)。



サイト

グラフ描画ggplot2の辞書的まとめ20のコード

ggplot2による可視化入門


biostatistics

これらのサイトはggplot2の基本的な使い方や代表的なグラフ、オプションの説明があります。
まだ慣れない時期はこれらのサイトを見比べながら少しずつイメージを掴むようにしてました。


From Data to Viz
スクリーンショット 2019-10-08 22.59.39



海外のサイトですが用途を選べばどのグラフを使えばいいか、そのコードサンプルを表示してくれます。そもそもどんなグラフを使えばいいか分からない場合は参考になります。




R Graphics Cookbook, 2nd edition
先程の「Rグラフィックスクックブック」の第2版です。
英語版では書籍がありますが、オンライ版は無料となっています。
ggplot2の辞書的な使い方ができますので、困ったことがあるとよくお世話になっています。
英語ではありますが、第3章を一通り試す+google翻訳でも参考になる部分は大きいと思っています。

Data Visualization
こちらも英語版ではありますが、Rを使ったグラフ作成について考え方から実際のコードまで紹介されています。中上級変ではありますが勉強になります。

Introduction to Data Science
私がRを勉強するのに利用したMOOC(Edx)のテキストです。英語ですが無料で読めます。
ggplot2に関しても詳しい説明があります。
こちらではデータ視覚化の原則(これだけは絶対にするな集)もあるので一読を進めます!



r-wakalang
Rに関するどんな質問でも気軽に相談できるslackになります。
Rは実はコミュニティーが平和で、初心者的な質問でも誰かがわかりやすく回答してくれる非常に珍しい(貴重な!)場になっています。
ggplot2に関してのスレッドもあります。


まとめ
今回はRのggplot2を使うときに参考になる書籍やサイトを紹介しました。
これからも参考になる情報があれば更新したいと思います。

(2019.09.17 更新)

前回はmutate関数を使って列の追加や修正を行いました。
そして最後に以下のコメントで終了しました。
mutate関数を使うと「FIMの合計点」とかだけでなく、年齢を「年代ごとに分類する」だったり、ある評価の「カットオフ値以上とそれ以下」に分けた列をExcelを使わずに作成することができます。
【2-2】Rのmutate関数を使って列の追加や修正を行うより)

今回は連続変数を条件に応じて複数のカテゴリーに変換していきます。

Rではいくつか方法があるのですが、今回はifelse関数、case_when関数、cut関数を紹介します。


今回使う架空のデータは年齢、性別、病名が入っています。
下のコードをまとめてスクリプトにコピペして実行してください。
library(tidyverse)
set.seed(1) 年齢 <- floor(rnorm(100,60,20)) 性別 <- sample(c("男性","女性"),100,replace = TRUE) 病名 <- sample(c("脳梗塞", "脳出血", "骨折", "靭帯損傷","心筋梗塞"), 100, replace = TRUE) data <- data.frame(年齢,性別,病名) head(data)
スクリーンショット 2019-02-10 22.48.49



今回は%>%やmutate関数を使って作業します。
これらの使い方に不安がある場合は前回の記事をご参照ください



1.ifelse関数をつかって分類する

条件から2つに分類する場合

ifelse関数は条件からTRUEとFALSEの2つに分けることができます。

ifelse(条件TRUEの場合FALSEの場合)



まず例としてdata_20という変数名に以下のように分類していきます。
・20歳未満を「未成年」、20歳以上を「成人」
・性別の男性を「male」、女性を「female」

data_20 <- data %>% 
  mutate(成人 = ifelse(.$年齢 < 20 , "未成年", "成人")) %>% 
  mutate(性別 = ifelse(.$性別 == "男性", "male", "female"))
head(data_成人)

スクリーンショット 2019-02-11 0.44.03

スクリーンショット 2019-02-11 0.46.48



条件から3つ以上分類する場合

3つ以上分類するはifelse関数のFALSEにifelseを入れ込みます。

ifelse(条件, TRUE, ifelse(条件2, TRUE, ifelse(条件3, TRUE, FALSE)))

このような形を「入れ子」と呼ぶそうです。繰り返せばいくらでも条件を増やすことができます。


今度は「年代」という変数名で20歳ごとに区切ってみます。
20歳未満なら「0〜」、
そうでなければ40歳未満なら「20〜」
そうでなければ60歳未満なら「40〜」・・・

上のイメージをコードにするとこうなります。

data_age <- data %>% 
  mutate(年代 = ifelse(.$年齢 < 20,"0〜", 
                      		ifelse(.$年齢 <40 ,"20〜", 
                        		ifelse(.$年齢 < 60,"40〜",
                        			ifelse(.$年齢 < 80, "60〜",
                        				ifelse(.$年齢 < 100, "80〜", "100〜")))))) %>% 
  mutate(年代 = factor(.$年代, levels = c("0〜", "20〜", "40〜", "60〜", "80〜", "100〜"), ordered = TRUE))
head(data_age)
スクリーンショット 2019-02-12 21.22.24

 mutate(年代 = factor(.$年代, levels = c("0〜", "20〜", "40〜", "60〜", "80〜", "100〜"), ordered = TRUE)の部分はfactorの順番を並べ替えています。

factorは50音順に並ぶので、もしこれをしなかったらグラフなど作る時に「0〜, 100〜, 20〜・・・」と100歳以上の場所がずれてしまいます。それを修正するためにlevels = でfactorの並べ替えをしています。
参照:【1-10】Rでよく使われる型について説明します。


個人的に感じる「ifelse関数+入れ子」のメリットは「他の関数を覚えなくていい」で、デメリットは終わりの))))))))の数が多すぎて間違えやすく、エラーに悩まされることです(笑)



2.case_when関数を使って分類する

3つ以上の分類で入れ子構造対策で慣れるとわかりやすいのがcase_when関数かもしれません。

case_when(条件A ~ 結果A, 条件B ~ 結果B, 条件C ~ 結果C)

今回は0歳以上 かつ 20歳未満といったように「かつ」が入る時は & を使います

年齢が0歳以上   かつ 年齢が20歳未満 を「0〜」
年齢が20歳以上 かつ 年齢が40歳未満 を「20〜」
年齢が40歳以上 かつ 年齢が60歳未満 を「60〜」・・・

上のイメージをコードにするとこうなります。

data_age <- data %>% 
  mutate(年代 = case_when(
    .$年齢 >= 0  & .$年齢 <20 ~ "0〜",
    .$年齢 >= 20 & .$年齢 <40 ~ "20〜",
    .$年齢 >= 40 & .$年齢 <60 ~ "40〜",
    .$年齢 >= 60 & .$年齢 <80 ~ "60〜",
    .$年齢 >= 80 & .$年齢 <100 ~ "80〜",
    .$年齢 >= 100 ~ "100〜")) %>% 
  mutate(年代 = factor(.$年代, levels = c("0〜", "20〜", "40〜", "60〜", "80〜", "100〜"), ordered = TRUE))

head(data_age)
スクリーンショット 2019-02-12 23.58.42


入れ子よりも( )の数が少ないのがありがたいところです。


数値でなくカテゴリーを条件で更に分類する場合

上記の場合は数値でしたが次のような場合はどうでしょう?

最初に行ったdata_20という変数名に以下のように分類していきます。
・成人という変数名で20歳未満を「未成年」、20歳以上を「成人」とする
・性別の男性を「male」、女性を「female」に変更
・疾患分類という変数名で、骨折・靭帯損傷を「運動器疾患」、脳梗塞・脳出血を「脳血管疾患」、心筋梗塞を「心疾患」とする

data_age <- data %>% 
  mutate(年代 = case_when(
    .$年齢 < 20 ~ "未成年",
    .$年齢 <= 20 ~ "成人")) %>% 
  mutate(性別 = case_when(
    .$性別 == "男性" ~ "male",
    .$性別 == "女性" ~ "female")) %>%
  mutate(疾患分類 = case_when(
    .$病名 %in% c("骨折", "靭帯損傷") ~ "運動器疾患",
    .$病名 %in% c("脳梗塞", "脳出血") ~ "脳血管疾患",
    .$病名 %in% "心筋梗塞" ~ "心疾患"))
head(data_age)
スクリーンショット 2019-02-13 0.58.43


ここで新たな ==%in% が出てきました。

スクリーンショット 2019-02-11 2.04.40
A == Bと A %in% Bは似ていますが、Bにあたる部分で==であれば1つしか入りません。
%in%は1つでも複数でも大丈夫です。

上記のコードで言うと .$性別 == "男性" ~ "male" の==は %in%に変えても大丈夫ですし、.$病名 %in% "心筋梗塞" ~ "心疾患" の %in% は == に変えても大丈夫です。

しかし.$病名 %in% c("骨折", "靭帯損傷") ~ "運動器疾患" に関しては ==に変えるとエラーが出ます。

上記図は今回の記事だけでなく条件式を作る時に必要になる考え方です。


case_when関数のメリットは入れ子構造にならないことで、デメリットとしては2つに分類するならifelseの方がコードが短くなります。



その他を作る(2019.09.17追記)

Exploratoryというソフトを開発しているKan Nishidaさんのツイートでcase_whenの解説を見つけました。ExploratoryはRをベースにしたデータ分析のソフトです。Rはプログラムを勉強してからでしか使えないところがありますが、ExploratoryはRコマンダーのようにプログラムを使わずにデータ分析ができます。

すべての人にデータサイエンスを

最新のデータサイエンスでは、さまざまな種類のデータを扱うことができ、機械学習や統計アルゴリズムを使用して、深い統計の知識がなくてもデータに隠れているパターンや傾向を見つけることができます。しかし、それはデータサイエンティストとプログラマーに限られています。

Exploratoryは誰もが、プログラミングを必要とせずに、モダンで最先端のデータサイエンスのアルゴリズムに簡単にアクセスできるようにします。私たちは誰もが最新のテクノロジーにアクセスし、データを通じて世界をより深く理解していくべきだと強く信じています。
(https://exploratory.io/より)





Exploratoryのサイトでもcase_when紹介がありました。

 


もし運動器疾患脳血管疾患以外をその他にしたければ以下のようになります。
TRUEはここでは残り全部といった意味になります。
data_age <- data %>% 
  mutate(年代 = case_when(
    .$年齢 < 20 ~ "未成年",
    .$年齢 <= 20 ~ "成人")) %>% 
  mutate(性別 = case_when(
    .$性別 == "男性" ~ "male",
    .$性別 == "女性" ~ "female")) %>%
  mutate(疾患分類 = case_when(
    .$病名 %in% c("骨折", "靭帯損傷") ~ "運動器疾患",
    .$病名 %in% c("脳梗塞", "脳出血") ~ "脳血管疾患",
    TRUE ~ "その他"

勉強になりました!


3.cut関数を使って分類する


数値を分類するに限って言うとcut関数も使えます。
cut関数はコードはスッキリしますが、少しクセがあり、おまじないが必要です。

cut(目的の列,
      breaks = (下限, カットする数値, 上限),
      right = FALSE,
      include.lowest = TRUE,
      labels =c("カテゴリー毎の名前"))

スクリーンショット 2019-02-13 1.49.04


「10代・20代」などでは0〜19のような分け方をするので、right = FALSEが必要になります。
加えてinclude.lowest = TRUEを加えないと端の値を読み込まず<NA>とデータなしとみなされてしまいます。


これらをふまえると以下のコードになります。

 data_age <- data %>% 
  mutate(年代 = cut(.$年齢, 
                  breaks = c(0, 20, 40, 60, 80, 100, 120),
                  right = FALSE, 
                  include.lowest = TRUE,
                  labels = c("0〜", "20〜", "40〜", "60〜", "80〜", "100〜")))
head(data_age)
スクリーンショット 2019-02-13 1.51.09


cut関数を使えば細かい条件が必要にならない分、right = や include.lowest = を忘れないように注意が必要です。


補足.mutate関数の組み合わせについて

mutate関数を複数回行うにはいくつかの方法があります。

スクリーンショット 2019-02-16 19.31.02
(追記)2019/02/16更新
すみません。mutate関数では②の方法は使えませんでした。
②の方法は次で紹介するfilter関数で使うことができます。


本来は,でつなぐと1つの行で済むのですが、初心者だと、()や,の数がうまくあわなかったりします。プログラミングの上級者からは意見があるかもしれませんが、プログラミング未経験・初心者はまず自分の覚えやすいものから始めていいと思っています。



まとめ

今回は条件によって分類する方法についてifelse関数、case_when関数、cut関数を紹介しました。

それぞれに特徴がありますので、必要に応じて使い分けてください。

また今回は演算子(<, ==, %in%など)の紹介も行いました。

演算子はデータを集計するときに必要になりますのでまた紹介していきたいと思います。

第3章ではggplot2を使ったグラフの作り方について説明してきました。

【3-1】ExcelにはないRでグラフを作るメリットと特徴

【3-2】ggplot2でグラフを作る流れを説明します

【3-3】Rのggplot2で散布図を作るgeom_point関数

【3-4】Rのggplot2でヒストグラムを作るgeom_histogram関数

【3-5】Rのggplot2で密度曲線を作るgeom_density関数

【3-6】Rのggplot2で箱ひげ図を作るgeom_boxplot関数

【3-7】棒グラフの基本とRのggplot2で棒グラフを作るgeom_bar関数

【3-8】ggplot2で折れ線グラフを作るgeom_line関数

【3-9】ggplot2でヒートマップを作るgeom_tile関数



【3-10】Rのggplot2でグループ毎にグラフを作りまとめて表示するfacet_wrap関数では3つのグラフを作成しました。今回はその3つのグラフを並べて1つのグラフにしてしまうgridExtraパッケージを紹介します。


1.使用するデータ

今回は【3-10】Rのggplot2でグループ毎にグラフを作りまとめて表示するfacet_wrap関数で取り込んだデータを使用します。データやコードについては以下をご参照ください。


今回はtidyverseパッケージに加えgridExtraパッケージを使用します。

上記パッケージを一度も使ったことがなければインストールします。
install.packages("tidyverse")
既にインストールしていればlibrary関数で呼び出します。
library(tidyverse) 


macの日本語文字化け対策にヒラギノ角ゴ Pro W3 ("HiraKakuPro-W3"を使っています。
windowの場合はエラーが出るかもしれませんのでtheme_gray(base_family = "HiraKakuPro-W3") +をの色をつけた箇所を消してください。

url <- "https://github.com/mitti1210/myblog/blob/master/fim.csv?raw=true" 
dat <- read.csv(url)

dat_stat <-
  dat %>% 
    gather(key = 項目, value = 点数, 食事:合計, factor_key = TRUE) %>% 
    group_by(項目, sheet) %>% 
    summarize(平均 = mean(点数)) %>% 
    mutate(平均 = round(平均, 1)) 

p_1 <-
dat_stat %>% 
  filter(項目 == "合計") %>% 
  ggplot(aes(x = sheet, y = 平均)) +
  theme_gray(base_family = "HiraKakuPro-W3") +
  geom_bar(aes(fill = sheet), stat = "identity") +
  geom_text(aes(label = 平均), position = position_stack(vjust = 0.5)) +
  facet_wrap( ~ 項目) +
  labs(x = "", y = "") +
  guides(fill = "none")

p_2 <-
dat_stat %>% 
  filter(項目 %in% c("運動合計", "認知合計")) %>% 
  ggplot(aes(x = sheet, y = 平均)) +
  theme_gray(base_family = "HiraKakuPro-W3") +
  geom_bar(aes(fill = sheet), stat = "identity") +
  geom_text(aes(label = 平均), position = position_stack(vjust = 0.5)) +
  facet_wrap( ~ 項目) +
  labs(x = "", y = "") +
  guides(fill = "none")

p_3 <-
dat_stat %>% 
  filter(!項目 %in% c("運動合計", "認知合計", "合計")) %>% 
  ggplot(aes(x = sheet, y = 平均)) +
  theme_gray(base_family = "HiraKakuPro-W3") +
  geom_bar(aes(fill = sheet), stat = "identity") +
  geom_text(aes(label = 平均), position = position_stack(vjust = 0.5)) +
  facet_wrap( ~ 項目) +
  labs(x = "", y = "") +
  guides(fill = "none")
スクリーンショット 2019-08-24 1.07.31
スクリーンショット 2019-08-24 1.07.41
スクリーンショット 2019-08-24 1.07.55


2.gridExtraパッケージについて


gridExtraパッケージのgrid.arrange関数は複数のグラフを並べて1つのグラフにしてくれます。



gridExtraパッケージを一度も使ったことがなければインストールします。
install.packages("gridExtra")
既にインストールしていればlibrary関数で呼び出します。
library(gridExtra) 



3.基本的な使い方

grid.arrange(グラフ1, グラフ2, ...)とつなげるだけです。そうすると縦に並びます。

grid.arrange(p_1, p_2)
スクリーンショット 2019-08-24 1.13.21



4.行数、列数の指定

行数と列数を指定することもできます。

grid.arrange(p_1, p_3, ncol = 1)
ncol = 1は横には1つ並べるという意味になります(何もつけないとこれになる)。
スクリーンショット 2019-08-27 7.17.26



grid.arrange(p_1, p_3, ncol = 2)
ncol = 2は横に2つ並べるとなります。
縦には制限なく並びます。

スクリーンショット 2019-08-27 7.17.34




grid.arrange(p_1, p_2, p_3, ncol = 2, nrow = 2)
ncol = 2は横に2つ並べる、nrow = 2は縦に2つ並べるという意味になります。
スクリーンショット 2019-08-27 7.17.45






5.大きさの比率を変える

縦と横の比率を変えることもできます。

横の比率を指定するwidths =と縦の比率を指定するheights =があります。
また1:2というように2つの数を指定する必要がありますので、width = c(1, 2)のように指定します。
もし横に3つ並べるならwidth = c(1, 3, 2)のように3つ数を指定します。

grid.arrange(p_1, p_2, p_3, ncol = 2, nrow = 2, widths = c(2, 1), heights = c(2,3))

スクリーンショット 2019-08-28 7.51.45





6.自分でレイアウトを作る

さらに細かいレイアウトも可能です。この方法は以下のサイトを参考にさせていただきました。




①自分でレイアウトを作る

3パターンレイアウトを作ってみました。
スクリーンショット 2019-08-28 21.55.43


ポイントはmatrix型にすることです。rbind関数c関数で縦につなぎます。
layout1 <- 
rbind(c(1, 2), c(1, 3)) layout2 <-
rbind(c(1, 3), c(2, 3)) layout3 <-
rbind(c(1, 2), c(3, 3))



②関数で作ったレイアウトを指定する

作ったレイアウトを指定するにはlayout_matrix = を指定します。
加えてwidths =heights =で微調整しています。


レイアウト1
grid.arrange(p_1, p_2, p_3, layout_matrix = layout1, widths = c(1, 3), heights = c(1, 2.5))
スクリーンショット 2019-08-28 22.00.35


レイアウト2
grid.arrange(p_1, p_2, p_3, layout_matrix = layout2, widths = c(1, 2))
スクリーンショット 2019-08-28 22.00.49



レイアウト3
grid.arrange(p_1, p_2, p_3, layout_matrix = layout3, heights = c(1, 2))
スクリーンショット 2019-08-28 22.01.01


7.まとめ

今回はパッケージを使って複数のグラフを1つにまとめる方法を紹介しました。

工夫次第でいろいろな見せ方ができます。gridExtraは便利なパッケージですのでぜひ試してみてください。










↑このページのトップヘ