Rで学ぶ計量経済学と機械学習 3

データの可視化

This version: 2026-10-08

Author

安藤道人(立教大学) 三田匡能 (株式会社 GA technologies)

1 本講義の目的

  • 要約統計量の算出方法を学ぶ。
  • ggplotやplotlyなどを用いたデータの可視化の基礎を学ぶ。

2 パッケージの読み込み

pacmanを使うと、install.packages()やlibrary()を使わずにパッケージをまとめて読み込むことができる。

まずpacmanをインストールする。

install.packages("pacman")

今回は、readxlとtidyverseとplotlyとskimrを使う。以下のように書けば、まとめて読み込むことができる。

pacman::p_load(readxl, tidyverse, plotly, skimr)

3 ディレクトリの設定

(Posit CloudではなくローカルPCで作業する場合)作業ディレクトリを設定する。 プロジェクト機能を活用できるのであれば、そちらを活用する。

setwd(自分のディレクトリ) # ディレクトリの設定
getwd() # ディレクトリの確認

4 例示用データの読み込み

前回同様、今回も test_scores.xlsx のデータを使う

df <- read_excel("test_scores.xlsx") # エクセル読み込み
head(df)
# A tibble: 6 × 5
  class name   math english japanese
  <chr> <chr> <dbl>   <dbl>    <dbl>
1 C     安藤     69      78       73
2 A     藤村     61      60       60
3 B     大泉     77      83       61
4 B     木村     45      49       68
5 A     鈴井     52      25       52
6 C     望月     79      65       74

5 データの分布を要約する

5.1 平均・分散など

データの分布を要約するための指標(要約統計量)の代表的なものが平均(mean)と分散(variance)である。

  • 平均(mean):\(n\)個のデータ\(x_i\)の総和を、サンプルサイズ\(n\)で割ったもの。

\[ \bar{x} = \frac{1}{n}\sum_{i=1}^n x_i \]

  • 分散(variance):平均からの差の2乗和の平均

\[ \sigma^2=\frac{1}{n}\sum_{i=1}^{n}(x_i-\bar{x})^2 \]

ただし、Rのvar()関数の算出する値は、母集団の分散(母分散)の不偏推定量としての分散すなわち不偏分散

\[ s^2=\frac{1}{n-1}\sum_{i=1}^{n}(x_i-\bar{x})^2 \]

である。

# 不偏分散
var(df$english)
[1] 291.0231

分散の単位は元の単位の2乗であるため、それを元の単位に戻すために平方根をとったものが標準偏差(standard deviation)である。

# 不偏標準偏差
sd(df$english)
[1] 17.0594

5.2 要約統計量

5.2.1 summary()

summary()関数は、データフレームの各行に関する代表値などの情報を取得することができる関数である。

summary(df)
    class               name                math          english     
 Length:40          Length:40          Min.   :35.00   Min.   :25.00  
 Class :character   Class :character   1st Qu.:57.75   1st Qu.:55.75  
 Mode  :character   Mode  :character   Median :70.50   Median :65.00  
                                       Mean   :67.78   Mean   :66.05  
                                       3rd Qu.:79.00   3rd Qu.:80.00  
                                       Max.   :96.00   Max.   :95.00  
    japanese     
 Min.   : 34.00  
 1st Qu.: 55.75  
 Median : 66.50  
 Mean   : 66.10  
 3rd Qu.: 75.00  
 Max.   :100.00  

classとnameは文字列(character)型なので、Length(行数、標本数)が40であることくらいしか表示されない。

math、english、japaneseは数値(numeric)型なので、代表値などの要約統計量(summary statistics)が表示されている。

それぞれの統計量の意味は

  • Min.が最小値
  • 1st Qu.が第一四分位数(1st quantile):「データを小さいものから順に並べた時に標本数の4分の1の位置(下位25%)にくる値」
  • Medianが中央値(median):「データを小さい順に並べた時に標本数の2分の1の位置にくる値」
  • Meanが平均値(mean):「データの総計をデータの個数で割ったもの」\(\sum_{i=1}^n x_i/n\)
  • 3rd Qu.が第三四分位数(3rd quantile):「データを小さい順に並べた時に標本数の4分の3の位置(75%)にくる値」
  • Max.が最大値

である。これらのうちMean以外の5つの統計量をまとめて五数要約という。

平均値は外れ値(極端に高かったり低かったりする値)の影響を受けやすいため、これらの統計量を見ることで、中央値と平均値が大きく異なるデータは外れ値があったりデータの分布に歪みがあることが推測できるわけである。

5.2.2 skimr::skim()

skimrのskim()を使うと、小さなヒストグラムつきの要約統計量の表を出力できる。欠損値の数なども表示してくれる。

skimr::skim(df)
Data summary
Name df
Number of rows 40
Number of columns 5
_______________________
Column type frequency:
character 2
numeric 3
________________________
Group variables None

Variable type: character

skim_variable n_missing complete_rate min max empty n_unique whitespace
class 0 1 1 1 0 3 0
name 0 1 1 3 0 38 0

Variable type: numeric

skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
math 0 1 67.78 13.82 35 57.75 70.5 79 96 ▂▅▆▇▂
english 0 1 66.05 17.06 25 55.75 65.0 80 95 ▂▃▇▅▅
japanese 0 1 66.10 14.43 34 55.75 66.5 75 100 ▂▇▇▆▂

6 2変数の関係を要約する

6.1 相関係数

相関係数(correlation coefficient)は2つの量的変数の直線的な関係の強さを測るもので、次の式で定義される

\[ \begin{align} x\text{と}y\text{の相関係数}r &= \frac{x\text{と}y\text{の共分散}}{x\text{の標準偏差}\times y\text{の標準偏差}} \\ &= \frac{\frac{1}{n} \sum_{i=1}^n (x_i - \bar{x})(y_i - \bar{y})} {\sqrt{\frac{1}{n} \sum_{i=1}^n (x_i - \bar{x})^2} \sqrt{\frac{1}{n} \sum_{i=1}^n (y_i - \bar{y})^2}} \end{align} \]

Rではcor()で計算できる。相関係数は-1から+1の範囲の値をとり、絶対値で1に近いほど関係が強いことを示す。

cor(x = df$english, y = df$math)
[1] 0.5063382

6.2 相関行列

相関係数が入った行列を相関行列という

cor(df[,3:5])
              math   english  japanese
math     1.0000000 0.5063382 0.3633715
english  0.5063382 1.0000000 0.6218165
japanese 0.3633715 0.6218165 1.0000000

7 描画用パッケージ

グラフの描画で役立つパッケージを紹介する。

7.1 ggplot2

{ggplot2}というグラフの描画用パッケージは綺麗なグラフを描くことができる。

ggplotは次のように2つ以上の手順をとり、2つ以上の関数を+でつなげていき、レイヤーを重ねるように描画していく。

  1. キャンバスを用意する(ggplot(データフレーム, aes(x, y)))
  2. グラフを用意する(geom_histogram()などgeom_なんたらという関数)
  3. グラフの装飾のための関数を使う場合は、さらにレイヤーを重ねていく
# 棒グラフをつくり、gに代入する
g <- ggplot(df, aes(x = `class`)) + # dfのデータでキャンバスを作る(`で囲っているのは関数のclass()と名前が被っているため)
  geom_bar() # キャンバスに棒グラフgeom_bar()を重ねる

g # 表示

7.2 plotly

{plotly}パッケージを使うと、簡単なインタラクティブなグラフ(カーソルを合わせると情報が表示されたり、表示情報を変更したりできるグラフ)を描くことができる。

ggplotly()関数を使うと、{ggplot2}で描いたグラフをインタラクティブなものに変換できる。

plotly::ggplotly(g)   # ggplot2のグラフをplotlyのグラフへ変換

8 1変数のグラフ

8.1 ヒストグラム

ヒストグラムは量的変数の分布を確認するには極めて有用なグラフである。

# ヒストグラム
ggplot(df, aes(x = english)) +
  geom_histogram(bins = 10) # binsは棒の数。指定しない場合はbins=30になる

なおR本体の機能でヒストグラムを描く場合はhist()を使う

# ヒストグラム
hist(x = df[["math"]])

8.2 棒グラフ

# ggplotの棒グラフ
g <- ggplot(df, aes(x = `class`)) +
  geom_bar()

#グラフ表示
g

#plotlyで表示
plotly::ggplotly(g) 

9 2変数のグラフ

9.1 散布図

2変数でグラフを描く場合、キャンバスにx(グラフ横軸)とy(グラフ縦軸)の両方を指定する。

そして散布図を描く場合はgeom_point()を+でつなげればよい

# ggplot
g <- ggplot(df, aes(x = english, y = japanese)) +
  geom_point() +
  xlab("English") + # X軸のラベルを変更
  ylab("Japanese") # y軸のラベルを変更

# 散布図表示
g

# plotlyで表示
plotly::ggplotly(g) 

回帰直線を加えるには、stat_smooth()を使う。

ggplot(df, aes(x = english, y = japanese)) +
  geom_point() +
  stat_smooth(method = "lm", se =FALSE, colour = "blue") # lm: 最小二乗法、信頼区間は表示させない(FALSE)、色は青(blue)

9.2 折れ線グラフ

折れ線グラフはgeom_line()で描くことができる。

Excelで言うところの「マーカー付き折れ線グラフ」にしたい場合はgeom_point()を重ねればよい。

ggplot(df, aes(x = english, y = japanese)) +
  geom_line() + # 折れ線
  geom_point()  # マーカー

(参考)アメリカのGNP

上では説明の簡単のため散布図と同じデータを使ったが、折れ線グラフは本来は時系列データに使うのが望ましい。

そこで、以下ではRにあらかじめ収録されている練習用のデータセットのひとつであるlongleyを使った例をのせる。

このデータは1947~62年のアメリカのGNPや雇用者数などが収録されている(詳細はhelpを参照してほしい)。

head(longley)
     GNP.deflator     GNP Unemployed Armed.Forces Population Year Employed
1947         83.0 234.289      235.6        159.0    107.608 1947   60.323
1948         88.5 259.426      232.5        145.6    108.632 1948   61.122
1949         88.2 258.054      368.2        161.6    109.773 1949   60.171
1950         89.5 284.599      335.1        165.0    110.929 1950   61.187
1951         96.2 328.975      209.9        309.9    112.075 1951   63.221
1952         98.1 346.999      193.2        359.4    113.270 1952   63.639
# ggplot
g <- ggplot(longley, aes(x = Year, y = GNP)) +
  geom_line()+ # 折れ線
  geom_point()+ # マーカー
  labs(title = "1947~62年のアメリカのGNP")

# グラフ表示
g

# plotlyで表示
ggplotly(g)

9.3 箱ひげ図

箱ひげ図は五数要約の統計量(最小値、第一四分位数、中央値、第三四分位数、最大値)をプロットすることで分布の概形を表示する図である。

g <- ggplot(df, aes(x = class, y = japanese)) +
  geom_boxplot()+   # 箱ひげ図
  labs(title = "クラスごとの国語の点数の分布")


ggplotly(g)

9.4 バイオリンプロット

バイオリンプロットは、「カーネル密度推定」という手法で推定されたデータ分布の曲線を左右対称に描いたもので、異なるデータの分布を視覚的に比較できる。デフォルトでは最小値から最大値の範囲でプロットされ、geom_violin(trim = FALSE)にすると曲線の端はカットされない。

g <- ggplot(df, aes(x = class, y = japanese)) +
  geom_violin() + # バイオリンプロット
  labs(title = "クラスごとの国語の点数の分布")

ggplotly(g)

9.5 箱ひげ図とバイオリンプロットの併用

g <- ggplot(df, aes(x = class, y = japanese)) +
  geom_violin() +             # バイオリンプロット。これを先に描く。
  geom_boxplot(width = 0.3) + # 箱ひげ図を上書き。幅を0.3に指定
  labs(title = "クラスごとの国語の点数の分布")

g  # ggplotly(g)だと箱ひげ図が表示されない

10 層別化

ggplot()関数の中に入れるaes()関数のfillやcolorといった引数に質的変数を指定することで、色の塗り分けによる層別化ができる。

これにより3つ目の変数の情報を表現できる。

ggplot(df, aes(x = english, y = japanese, fill = class, color = class)) +
  geom_point() +
  labs(title = "英語と国語の得点の散布図をクラスごとに塗り分けたもの")

11 バブルチャート

散布図のサークルの大きさを、第三の変数に合わせて変更することもできる。

# ggplot
g <- ggplot(df, aes(x = english, y = japanese, size = math)) +
  geom_point(alpha = 0.5) + # alpha:バブルの透過度
  labs(title = "英語と国語の得点の散布図(バブルサイズは数学の点数を反映)") +
  scale_size(range = c(1,10)) # バブルサイズの範囲を設定

# 散布図表示
g

#plotly
ggplotly(g) 

12 グラフの保存

図の保存にはggsaveを使用する

ggsave(plot = g, file = "bubble_chart.png")

(参考)事前に集計をしてからplotを描く場合

geom_bar()はstatという引数があり、デフォルトではgeom_bar(stat="count")と、観測値の数をカウントして表示するようになっている。そのためaes()にはxを指定するだけで自動的にy軸も描画された。

事前に集計を行って自分でx軸とy軸の変数を用意し、geom_bar(stat="identity")と指定することもできる。

# クラスごとに国語の平均点を集計
mean_ja <- df |> group_by(class) |> summarise(mean_ja=mean(japanese))
mean_ja
# A tibble: 3 × 2
  class mean_ja
  <chr>   <dbl>
1 A        62.2
2 B        74.1
3 C        60.7
# プロット
mean_ja |> ggplot(aes(x=`class`, y=mean_ja)) + geom_bar(stat="identity")

# pivot_longerで縦持ちデータに変え、科目ごとに平均得点を計算して描画する
g <- df |>
  tidyr::pivot_longer(cols=c(japanese, english, math),
                      names_to="subject",
                      values_to="score") |>
  dplyr::group_by(subject) |>
  dplyr::summarise(mean_score=mean(score)) |>
  ggplot(aes(x=subject, y=mean_score)) +
  geom_bar(stat="identity")

ggplotly(g)

13 参考文献

データ可視化については、様々な情報やテキストがある。また、作りたいグラフのイメージに近い例をウェブで探し、そのRコードを参考にするのも一つの方法である。

14 (おまけ)Pythonで同様の作業を行う場合

本講義の演習はRで行うが、参考として、ここまで学んだ内容に対応するPythonのコードを簡単に紹介する。

Pythonでの可視化には、主にmatplotlib(基本的な描画)、seaborn(matplotlibをベースにした統計グラフ用のライブラリで、ggplot2に近い感覚で使える)、plotly(インタラクティブなグラフ)を使う。

# ターミナル上で実行(Rのinstall.packages()に相当)
# pip install pandas openpyxl matplotlib seaborn plotly statsmodels
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import plotly.express as px

# matplotlibはそのままだと日本語が□□□になるので、日本語フォントを指定する
plt.rcParams["font.family"] = "Hiragino Sans"  # Macの場合(Windowsなら"Yu Gothic"など)

df = pd.read_excel("test_scores.xlsx")

14.1 データの分布を要約する

pandasの.var()・.std()は、Rのvar()・sd()と同じく不偏分散・不偏標準偏差を返す。.describe()はRのsummary()に相当する。

df["english"].mean()   # 平均
df["english"].var()    # 不偏分散(Rの var() に相当)
df["english"].std()    # 不偏標準偏差(Rの sd() に相当)

df.describe()               # 要約統計量(Rの summary() に相当)
df.describe(include="all")  # 文字列の列も含めて表示(Rの skimr::skim() に近い)

14.2 2変数の関係を要約する

df["english"].corr(df["math"])            # 相関係数(Rの cor(x, y) に相当)
df[["math", "english", "japanese"]].corr()  # 相関行列(Rの cor(df[,3:5]) に相当)

14.3 描画用ライブラリ

seabornは「データフレーム・x軸・y軸を指定してグラフの関数を呼ぶ」という、ggplot2に近い書き方ができる。plotlyのplotly.expressを使うと、同じようにしてインタラクティブなグラフを描ける(Rのggplotly()に近い役割)。

# 棒グラフ(Rの ggplot(df, aes(x = class)) + geom_bar() に相当)
sns.countplot(data=df, x="class", order=["A", "B", "C"])
plt.show()

# plotlyでインタラクティブな棒グラフ
px.histogram(df, x="class")

14.4 1変数のグラフ

# ヒストグラム(Rの geom_histogram(bins = 10) に相当)
sns.histplot(data=df, x="english", bins=10)
plt.show()

14.5 2変数のグラフ

# 散布図(Rの geom_point() に相当)
sns.scatterplot(data=df, x="english", y="japanese")
plt.show()

# 回帰直線つきの散布図(Rの stat_smooth(method = "lm", se = FALSE) に相当)
sns.regplot(data=df, x="english", y="japanese", ci=None)
plt.show()

折れ線グラフの例として、Rのlongleyと同じデータはstatsmodelsというライブラリに収録されている。

import statsmodels.api as sm
longley = sm.datasets.longley.load_pandas().data

# マーカー付き折れ線グラフ(Rの geom_line() + geom_point() に相当)
sns.lineplot(data=longley, x="YEAR", y="GNP", marker="o")
plt.title("1947~62年のアメリカのGNP")
plt.show()

箱ひげ図・バイオリンプロットはsns.boxplot()・sns.violinplot()で描ける。seabornのバイオリンプロットはデフォルトで最小値・最大値の外側まで曲線が伸びるので、Rのgeom_violin()に合わせるにはcut=0を指定する。

# 箱ひげ図(Rの geom_boxplot() に相当)
sns.boxplot(data=df, x="class", y="japanese", order=["A", "B", "C"])
plt.title("クラスごとの国語の点数の分布")
plt.show()

# バイオリンプロットと箱ひげ図の併用
ax = sns.violinplot(data=df, x="class", y="japanese", order=["A", "B", "C"],
                    cut=0, inner=None, color="white")
sns.boxplot(data=df, x="class", y="japanese", order=["A", "B", "C"], width=0.3, ax=ax)
plt.show()

14.6 層別化・バブルチャート

hue(色)やsize(大きさ)に変数を指定する(Rのaes()のcolor・sizeに相当)。

# クラスごとに色分けした散布図
sns.scatterplot(data=df, x="english", y="japanese", hue="class")
plt.show()

# バブルチャート(バブルの大きさが数学の点数)
sns.scatterplot(data=df, x="english", y="japanese", size="math",
                sizes=(20, 400), alpha=0.5)
plt.show()

# plotlyでインタラクティブなバブルチャート
px.scatter(df, x="english", y="japanese", size="math")

14.7 グラフの保存

plt.savefig("bubble_chart.png")   # matplotlib/seabornのグラフ(Rの ggsave() に相当)。plt.show()の前に実行する
fig = px.scatter(df, x="english", y="japanese", size="math")
fig.write_html("bubble_chart.html")  # plotlyのグラフはhtmlとして保存できる

14.8 (参考)事前に集計をしてからplotを描く場合

# クラスごとに国語の平均点を集計(Rの group_by() |> summarise() に相当)
mean_ja = df.groupby("class", as_index=False)["japanese"].mean()
sns.barplot(data=mean_ja, x="class", y="japanese")
plt.show()

# 縦持ちデータに変え(Rの pivot_longer() に相当)、科目ごとの平均点を描画
mean_score = (
    df.melt(id_vars=["class", "name"], value_vars=["japanese", "english", "math"],
            var_name="subject", value_name="score")
    .groupby("subject", as_index=False)["score"].mean()
)
px.bar(mean_score, x="subject", y="score")