Zenn
🐹

Rで始める（３）データの探索

に公開
データ活用
tech

データ探索をする場合，特に「データの分布」などについて調べようとする場合，この資料[1]が歴史的に模範となる。この資料の2ページに掲載されたデータはBase Rのデフォルトに設定済みで利用可能となっている。データの内容はお花のiris(アヤメ)についてであり，その種類ごとのSepal(がく片)よりPetal(花びら)のLength(長さ)とWidth(幅)を調査をしたものである。しかしながら，その主な内容は数字羅列のデータに過ぎない。しかし，この資料には嬉しいことにデータの分布を示す手法が紹介されている。データを探索し，データの分布を示すことでDSのための良い模範となる。下の写真は日本でもよく見かけるお花・アヤメである。

図0　irisアヤメの写真（筆者が長崎市内で撮影）

irisデータ


図1　irisアヤメのSetosa, Versicolor, Virginica品種

irisアヤメ（菖蒲、文目、綾目、学名:Iris sanguinea）というお花はアヤメ科アヤメ属の多年草で，病気を治す薬としても知られているそうだ[2]。本データの作者R.A.Fisherはこのアヤメのがく片と花びらを調べて，種別(Species)にその特性を調べるためにデータの探索・分布の考え方[1:1]を提案している。なぜ，R.A.Fisherはがく片と長さに注目したのだろうか。それはお花の特徴を綿密に探索・分布を調べて気付いたものだ。品種ごとの特徴である。Setosa品種はがく片が大きく，花びらはそれに比べると小さい形をしている。他品種のVersicolorになるとがく片はSetosaよりだんだんと縮小され，代わりに花びらがSetosaより大きくなっていく。更にVirginica品種になるとがく片と花びらの大きさがSetosa品種と逆転することがわかる。Fisherの報告[1:2]では最終的にこれらの分布をヒストグラム手法を導入しており，irisアヤメの特徴をわかりやすく分析している。

本資料ではデータの特徴をわかりやすく説明できる「データ探索と分布」を紹介するためにR言語を導入されている。R言語には最初からこのirisデータが組み込まれているのが一般的になっている。

データの中身

R.A.Fisherが調査整理をしたデータは，お花iris(アヤメ)の3種(Setosa,Versicolor,Virginica)の，Sepal(がく片)のLength(長さ)とWidth(幅)およびPetal(花びら)のLength(長さ)とWidth(幅)である。ひとつの品種に対して50個のお花を調査をしているので，3品種を調査することで全部で150個(行)のデータセットとなっている。データの調査記録には1個(行)あたりに品種名を含めて「Speralがく片のWidth幅とLength長さ」「Petal花びらのWidth幅とLength長さ」の全部で5種類(列)のデータがセットである。結果，この5列のデータセットが150行続く。

以下の図2はSetosa品種の調査項目を説明している。Setosa品種の，品種名＋4種類（Sepal.Length, Sepal.Width, Petal.Lenght, Petal.Width）の合計5列となる調査項目が示されており，これが1行の合成セットデータである。このセットがSetosaだけで50個，Versicolorだけで50個，Virginicaだけで50個の全部で150個(行)が用意されている。

図2はFisherの調査内容を推定したものである。


図2　Setosa品種についての調査内容

Rによりこれらのデータセットの内容を確認するには以下の命令文で実施される。

RSource-0 : アヤメ(iris)の調査内容

{r-0}
iris
Sepal.Length Sepal.Width Petal.Length Petal.Width Species
1 5.1 3.5 1.4 0.2 setosa
2 4.9 3.0 1.4 0.2 setosa
3 4.7 3.2 1.3 0.2 setosa
- - - - - -
50 5.0 3.3 1.4 0.2 setosa
51 7.0 3.2 4.7 1.4 versicolor
52 6.4 3.2 4.7 1.5 versicolor
- - - - - -
100 5.7 2.8 4.1 1.3 versicolor
101 6.3 3.3 6.0 2.5 virginica
102 5.8 2.7 5.1 1.9 virginica
- - - - - -
150 5.9 3.0 5.1 1.8 virginica

データの探索・分布調べ

データの探索や分布を調べるにはデータをグラフ化することから始めるのが良い。そのためにはどんなデータをどのようなグラフで表現するかを決める必要がある。本資料ではそのデータ整理とグラフ化にRを用いる。このRには簡単で，尚且つ非常に便利なツールが豊富に揃っている。特に，本資料ではggplot2のライブラリーを想定しているが，所々BaseRによるグラフ化を比較説明する。

度数の導入

度数とは統計学では良く現れるもので，度を数えた数値である。一般的にはヒストグラムと呼ばれている場合が多い。例えば，アヤメの花びらの長さ(Petal.length)は1.0cmから6.9cmまで分布しているが，長さの同じ花びらが何度現れるかを数えられば度数が得られる。しかし，自然の花びらの長さは決して同じものが現れにくい。そこで花びらの長さの現れる範囲を絞り，例えば1.0cmから1.2cmまでの花びらの長さが現れる度に数えて数値化する「度数」を記録。この数を縦軸に横軸には長さの範囲を示すグラフを作成するとこれが度数分布表でありヒストグラムとも呼ばれている。

Rではこれを実現するための手順として以下の過程を踏む。

  • 花びらの長さの範囲を決める
  • 決められた長さの範囲内に現れる度数を数える
  • 範囲を横軸に，度数を縦軸に，グラフ化をする

RSource-1 : ヒストグラム表示のRソース

{r-1}

breaksP <- seq(from=min(iris[["Petal.Length"]]),to=max(iris[["Petal.Length"]]),length=31)
# Petal.length=31種類(30カテゴリ)に分類し，最小値から最高値までの値を決める
PL <- cut(iris[["Petal.Length"]], breaks=breaksP, right=TRUE, include.lowest=TRUE)
#決めた30カテゴリ，30種類の範囲を設定，30個の束を作る

# GGPLOT2無しの場合，度数分布表のヒストグラム表示
hist(iris[["Petal.Length"]], breaks=breaksP, 
     main="Petal Length Histogram", 
     xlab="Petal Length", 
     ylab="Frequency", 
     col="lightblue", 
     border="white")

# GGPLOT2有りの場合，以下libraryを読み込む，pacmanでLoadした場合は不要
library(ggplot2)
# GGPLOT2有りの場合，テーブルをフレーム型に変換
freq_data <- as.data.frame(table(PL))
colnames(freq_data) <- c("Range", "Frequency") # 列名を範囲と度数に

# GGPLOT2有りの場合，度数分布表のヒストグラム表示
ggplot(freq_data, aes(x = Range, y = Frequency)) +
  geom_bar(stat = "identity", fill = "lightblue", color = "black") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
                                # x軸ラベルを45度傾ける
  labs(title = "Petal Length Distribution",
       x = "Petal Length Range",
       y = "Frequency")


結果１：x = Petal.Lengthとして指定したGGPLOT2無しの表示

結果２：x = Petal.Lengthとして指定したGGPLOT2有りの表示

密度の導入

上記のヒストグラムは「花びらの長さ」を調べているが，例えば「がく片」と「花びら」をヒストグラムで表示するグラフ化には見栄えなどを考えて，ヒストグラムを滑らかな曲線で繋いだグラフを作ることができる。これが密度によるグラフ化である。例えば，irisアヤメの3品種の中でSetosaの品種に注目すれば，そのがく片(Sepal.Length, Sepal.Width), 花びら(Petal.Length, Petal.Width)の違いをグラフ化できる。この場合，Setosa品種に対して長さに注目したSepal.LengthとPetal.Lengthの2種類と，幅に注目したSepal.WidthとPetal.Widthの2種類の，合計4種類になる。長さに注目したSepal.LengthとPetal.Lengthの滑らかなグラフを書いてみる。

RSource-2 : アヤメ(iris)の種別グラフ作成

{r-2}

library(ggplot2)
ggplot(iris, aes(x = Sepal.Length, fill = Species)) +
  geom_density(alpha = 0.7) +
  labs(title = "種別ごとのがく片の長さの分布",
       x = "がく片の長さ (cm)",
       y = "密度") +
  theme_minimal()

ggplot(iris, aes(x = Petal.Length, fill = Species)) +
  geom_density(alpha = 0.7) +
  labs(title = "種別ごとの花びらの長さの分布",
       x = "花びらの長さ (cm)",
       y = "密度") +
  theme_minimal()


結果１：x = Sepal.Lengthとして指定した場合の結果

結果２：x = Petal.Lengthとして指定した場合の結果

ヒストグラム作成開始

上記のR言語ソースを見ながら，geom_density()の代わりにgeom_histogram()を使ってたらヒストグラムになると思いがちである。しかし，geom_density()をgeom_histogram()に入れ替えるだけではエラーが生じる[3][4]。その特性は以下のとおりである。

geom_density(), geom_histogram()

GGPLOT2で使用するgeom_density()は言葉通り密度Densityをグラフ化する関数で，geom_histogram()はデータの頻度をグラフ化する関数なのだ。「頻度」は一般的に「度数」と表現されることが多く意味は同じもの。一度，二度，三度などで数を数える「度数」は「頻度」を意味するからだ。

iris品種の特徴を調べる場合，iris品種の何かを数えて「頻度・度数」を「密度」グラフとして書いて見ることになる。上記の結果1と結果2ではiris品種ごとの長さに注目したSepal.Length及びPetal.Lengthを「密度」としてグラフ化している。Sepal.Lengthを根拠に「密度」グラフを作成したのが結果1で，Petal.Lengthを根拠に「密度」グラフを作成したのが結果2である。密度とは正確に言えば「相対度数」であり，色で塗りつぶされた範囲のグラフ面積が1になっている。密度は相対度数だが「度数」としても作成可能なのがgeom_histogram()である。以下にそのソースを示す。

RSource-3 : ヒストグラム表示のRソース

{r-3}
library(ggplot2)
ggplot(iris, aes(x = Sepal.Length, fill = Species)) +
  geom_histogram(aes(y = ..density..), alpha = 0.7, position = "identity", binwidth = 0.2) +
  labs(title = "種別ごとのがく片の長さの分布",
       x = "がく片の長さ (cm)",
       y = "密度") +
  theme_minimal()

ggplot(iris, aes(x = Petal.Length, fill = Species)) +
  geom_histogram(alpha = 0.7, position = "identity", bins = 30) +
  labs(title = "種別ごとの花びらの長さの分布",
       x = "花びらの長さ (cm)",
       y = "度数") +
  theme_minimal()

上記のソース中の，geom_histogram(aes(y = ..density..),...)は，Y軸のスケール単位を「密度」にするという意味。aes(y = ..density..)のオプションが無い場合は，Y軸の単位は命令通りのgeom_histogram()となり，自ずとY軸は「度数」の単位となる。度数は「頻度」を指すので度数分布のヒストグラム形式を取るが，中身は密度は相対度数の「確率単位」を取りたいという意味。また，ヒストグラムを書く場合，Sepal.LengthとPetal.Lengthの長さ範囲は0.2cm間隔で度数を数え指定。但し(binwidth = 0.2)を，全体のスケールから均等に30等分(bins = 30)にしたい場合も対応可能。この命令でヒストグラムの棒グラフの太さを決める。


結果３：aes(y = ..density..)有りで，binwidth = 0.2として指定した場合の結果

結果４：aes(y = ..density..)を削除して，bins = 30として指定した場合の結果

度数と相対度数の密度

「度数」とは度重ねた数を意味するものだ。下記に示すテーブルはirisアヤメのSetosa品種を示す50個のデータだ。データの数値を確認してみると，Sepalがく片は大きく，Petal花びらは小さい。Length長さやWidth幅も同様な傾向。このような傾向をグラフ化したい場合「度数や相対度数の密度」として表現することができる。まず，数える範囲を決めてその範囲内にある数を数える作業が必要である。例えば，Sepal.Length(がく片の長さ)の数え範囲を0.2cm毎に決め，該当するデータの数を数えて密度（geom_density()）として表現。別の方法としてPetal.Length(花びらの長さ)の範囲を0.2cm毎（binwidth = 0.2）に区分して，区分ごとに該当するデータの数を数えて度数（geom_histogram()）として表示できる。

度数は数を数えた数値そのものを使用するが相対度数の密度は数えた数を全体の数（50個）で割る処理をした結果の相対度数の密度である。度数は数えた数だけに意味があるが，密度はSetosa品種全体の中で区切られた範囲（0.2cm毎か，30等分の）の数が他の範囲とどのような違いがあるかを説明する相対的な度数として利用できる。

Sepal.Length Sepal.Width Petal.Length Petal.Width Species
1 5.1 3.5 1.4 0.2 setosa
2 4.9 3.0 1.4 0.2 setosa
3 4.7 3.2 1.3 0.2 setosa
4 4.6 3.1 1.5 0.2 setosa
5 5.0 3.6 1.4 0.2 setosa
- - - - - -
50 5.0 3.3 1.4 0.2 setosa

RSource-4 : 度数と密度を同時に表示するRソース

{r-4}
library(ggplot2)
ggplot(iris, aes(x = Petal.Length, fill = Species)) +
  geom_density(alpha = 0.4) +
  geom_histogram(alpha = 0.4, position = "identity", bins = 30) +
    # identityの代わりにdodgeを入れる場合は，color="black"などを追加して棒で区別
  labs(title = "種別ごとのがく片の長さの分布",
       x = "がく片の長さ (cm)",
       y = "度数と密度") +
  theme_minimal()


結果５：geom_dinsity(),geom_histogram()を併用した結果

上記の結果を見ると，区切られた範囲内の度数と，その度数を全体度数に割った密度では単位が一桁違う。結果として同じ意味合いのグラフなのに，縮尺が違うため説得力に欠ける。そこで，以下のように改変したR言語のソースを利用する。前述した（geom_histogram(aes(y = ..density..))を取り入れ，比較しやすいグラフを作成する。正確に言えば，度数を密度の単位に合わせている。

RSource-5 : 度数を密度として密度と同時に表示するRソース

{r-5}
library(ggplot2)
ggplot(iris, aes(x = Petal.Length, fill = Species)) +
  geom_density(alpha = 0.4) +
  geom_histogram(aes(y = ..density..), alpha = 0.4, position = "identity", bins = 30) +
  labs(title = "種別ごとのがく片の長さの分布",
       x = "がく片の長さ (cm)",
       y = "密度として") +
  theme_minimal()


結果６：geom_dinsity(),geom_histogram(aes(y = ..density..))の結果

練習問題

  1. 種別ごとのがく片の幅の分布
  2. 種別ごとの花びらの幅の分布
  3. がく片の長さと花びらの長さを一つのグラフに表示する

練習問題３の予想結果

脚注

  1. https://onlinelibrary.wiley.com/doi/epdf/10.1111/j.1469-1809.1936.tb02137.x ↩︎ ↩︎ ↩︎

  2. https://ja.wikipedia.org/wiki/アヤメ ↩︎

  3. https://www.jaysong.net/tutorial/R/ggplot_intro2.html ↩︎

  4. https://ggplot2.tidyverse.org/index.html ↩︎

Discussion