PCA – Breast cancerデータセット

2020-11-20 / tau / コメントする

概要

scikit-learnの主成分分析モデル(PCA)をBreast cancerデータセットに適用して、その挙動を確認する。

30個の特徴量（全て連続量）を持つ569個の腫瘍データについて、悪性(marignant)／良性(benign)がターゲットとして与えられている。PCAによって特徴量のみの分析で、少ない主成分によってある程度明確な分離が可能なことが示される。

手順

以下の手順・コードで計算した。

パッケージをインポート
Breast cancerデータセットを準備
データセットをスケーリング
- StandardScalerで特徴量データを標準化している
PCAモデルのインスタンスを生成
- 引数n_components=3で3つの主成分まで計算させている
fit()メソッドによって、モデルにデータを学習させる
成分やその寄与率を確認
- 主成分はPCA.comonents_を、寄与率はPCA.explained_variance_ratio_を確認
transform()メソッドによって、主成分に沿ってデータを変換
3つの主成分について3次元可視化
2つの主成分について2次元可視化

主成分と寄与率

以下に主成分と寄与率を計算するまでのコードを示す。

import pandas as pd
import matplotlib.pyplot as plt
from sklearn.datasets import load_breast_cancer
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
from mpl_toolkits.mplot3d import Axes3D

ds = load_breast_cancer()
X_scaled = StandardScaler().fit_transform(ds.data)
df = pd.DataFrame(X_scaled, columns=ds.feature_names)

pca = PCA(n_components=3)
pca.fit(df)

print("Components:\n{}".format(pca.components_))
# Components:
# [[ 0.21890244  0.10372458  0.22753729  0.22099499  0.14258969  0.23928535
#    0.25840048  0.26085376  0.13816696  0.06436335  0.20597878  0.01742803
#    0.21132592  0.20286964  0.01453145  0.17039345  0.15358979  0.1834174
#    0.04249842  0.10256832  0.22799663  0.10446933  0.23663968  0.22487053
#    0.12795256  0.21009588  0.22876753  0.25088597  0.12290456  0.13178394]
#  [-0.23385714 -0.0597061  -0.21518137 -0.23107672  0.18611299  0.15189161
#    0.06016536 -0.03476749  0.19034877  0.36657548 -0.10555215  0.0899797
#   -0.08945723 -0.15229263  0.20443047  0.23271588  0.19720726  0.13032158
#    0.18384799  0.28009201 -0.21986638 -0.0454673  -0.19987843 -0.21935186
#    0.17230435  0.14359318  0.09796412 -0.00825721  0.14188335  0.27533948]
#  [-0.00853135  0.06454967 -0.00931432  0.02869944 -0.10429256 -0.07409161
#    0.00273383 -0.0255634  -0.0402399  -0.02257393  0.26848142  0.37463394
#    0.26664551  0.21600645  0.30883931  0.15477941  0.17646341  0.22465797
#    0.28858416  0.21150348 -0.04750699 -0.0422979  -0.04854648 -0.0119023
#   -0.25979771 -0.23607553 -0.17305727 -0.17034366 -0.27131257 -0.23279102]]

print("Explained :{}".format(pca.explained_variance_ratio_))
# Explained :[0.44272026 0.18971182 0.09393163]

import pandas as pd

import matplotlib.pyplot as plt

from sklearn.datasets import load_breast_cancer

from sklearn.preprocessing import StandardScaler

from sklearn.decomposition import PCA

from mpl_toolkits.mplot3d import Axes3D

ds = load_breast_cancer()

X_scaled = StandardScaler().fit_transform(ds.data)

df = pd.DataFrame(X_scaled, columns=ds.feature_names)

pca = PCA(n_components=3)

pca.fit(df)

print("Components:\n{}".format(pca.components_))

# Components:

# [[ 0.21890244 0.10372458 0.22753729 0.22099499 0.14258969 0.23928535

# 0.25840048 0.26085376 0.13816696 0.06436335 0.20597878 0.01742803

# 0.21132592 0.20286964 0.01453145 0.17039345 0.15358979 0.1834174

# 0.04249842 0.10256832 0.22799663 0.10446933 0.23663968 0.22487053

# 0.12795256 0.21009588 0.22876753 0.25088597 0.12290456 0.13178394]

# [-0.23385714 -0.0597061 -0.21518137 -0.23107672 0.18611299 0.15189161

# 0.06016536 -0.03476749 0.19034877 0.36657548 -0.10555215 0.0899797

# -0.08945723 -0.15229263 0.20443047 0.23271588 0.19720726 0.13032158

# 0.18384799 0.28009201 -0.21986638 -0.0454673 -0.19987843 -0.21935186

# 0.17230435 0.14359318 0.09796412 -0.00825721 0.14188335 0.27533948]

# [-0.00853135 0.06454967 -0.00931432 0.02869944 -0.10429256 -0.07409161

# 0.00273383 -0.0255634 -0.0402399 -0.02257393 0.26848142 0.37463394

# 0.26664551 0.21600645 0.30883931 0.15477941 0.17646341 0.22465797

# 0.28858416 0.21150348 -0.04750699 -0.0422979 -0.04854648 -0.0119023

# -0.25979771 -0.23607553 -0.17305727 -0.17034366 -0.27131257 -0.23279102]]

print("Explained :{}".format(pca.explained_variance_ratio_))

# Explained :[0.44272026 0.18971182 0.09393163]

寄与率は第1主成分が44%、第2主成分が19%、第3主成分が9%。第3成分まで3/4の情報を説明していることになる。

また、第1主成分は全ての特徴量がプラス方向で寄与している。

主成分をヒートマップで視覚化してみると、各主成分の符号や大きさが直感的に把握しやすくなるが、第2～第3主成分がmeanとworst系の特徴量が小さい方が影響が大きい点、3つの主成分についてerrorが大きいほど影響が大きい点など、意味づけは難しい。

plt.matshow(pca.components_, cmap='seismic')
plt.yticks([0, 1, 2], ["Component-1", "Component-2", "Comoponent-3"])
plt.colorbar()
plt.xticks(range(len(ds.feature_names)),
    ds.feature_names, rotation=60, ha='left')
plt.xlabel("Feature")
plt.ylabel("Principal components")

plt.show()

plt.matshow(pca.components_, cmap='seismic')

plt.yticks([0, 1, 2], ["Component-1", "Component-2", "Comoponent-3"])

plt.colorbar()

plt.xticks(range(len(ds.feature_names)),

ds.feature_names, rotation=60, ha='left')

plt.xlabel("Feature")

plt.ylabel("Principal components")

plt.show()

可視化

3次元

3つの主成分について3次元で可視化してみると、悪性／良性がかなりはっきりと分離されている。

2次元

2つの主成分のみでも、悪性／良性がよく区分されている。

まとめ

Irisデータの場合と同じく、特徴量分析のみでクラスの別がよくあぶりだされている。

PCA – Irisデータセット

2020-11-20 / tau / コメントする

概要

scikit-learnの主成分分析モデル(PCA)をIrisデータに適用して、その挙動を確認する。

クラス分類のターゲットを用いていないにもかかわらず、少ない主成分でクラスがかなり明確に分類されることがわかる。

計算の手順

以下の手順・コードで計算した。

必要なパッケージをインポート
Irisデータセットを準備
データセットをスケーリング
- StandardScalerで特徴量データを標準化している
PCAモデルのインスタンスを生成
- 引数n_componentsを指定せず、4つの特徴量全てを計算
モデルにデータを学習させる
- fit()メソッドのみでよいが、後のグラフ化のためにfit_transform()メソッドを実行
- X_transに主成分によって変換したデータを格納
主成分やその寄与率を確認
- 主成分はPCA.comonents_を、寄与率はPCA.explained_variance_ratio_を確認
3つの主成分について3次元可視化
2つの主成分について2次元可視化

主成分と寄与率

以下に、主成分と寄与率を計算するまでのコードを示す。

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.datasets import load_iris
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
from mpl_toolkits.mplot3d import Axes3D

iris_data = load_iris()

X = iris_data['data']
y = iris_data['target']
feature_names = iris_data['feature_names']
target_names = iris_data['target_names']

X_scaled = StandardScaler().fit_transform(X)

pca = PCA()
X_trans = pca.fit_transform(X_scaled)

print(pd.DataFrame(pca.components_, columns=feature_names))
#    sepal length (cm)  sepal width (cm)  petal length (cm)  petal width (cm)
# 0           0.521066         -0.269347           0.580413          0.564857
# 1           0.377418          0.923296           0.024492          0.066942
# 2          -0.719566          0.244382           0.142126          0.634273
# 3          -0.261286          0.123510           0.801449         -0.523597

print(pca.explained_variance_ratio_)
# [0.72962445 0.22850762 0.03668922 0.00517871]

import numpy as np

import pandas as pd

import matplotlib.pyplot as plt

from sklearn.datasets import load_iris

from sklearn.preprocessing import StandardScaler

from sklearn.decomposition import PCA

from mpl_toolkits.mplot3d import Axes3D

iris_data = load_iris()

X = iris_data['data']

y = iris_data['target']

feature_names = iris_data['feature_names']

target_names = iris_data['target_names']

X_scaled = StandardScaler().fit_transform(X)

pca = PCA()

X_trans = pca.fit_transform(X_scaled)

print(pd.DataFrame(pca.components_, columns=feature_names))

# sepal length (cm) sepal width (cm) petal length (cm) petal width (cm)

# 0 0.521066 -0.269347 0.580413 0.564857

# 1 0.377418 0.923296 0.024492 0.066942

# 2 -0.719566 0.244382 0.142126 0.634273

# 3 -0.261286 0.123510 0.801449 -0.523597

print(pca.explained_variance_ratio_)

# [0.72962445 0.22850762 0.03668922 0.00517871]

寄与率の方を見てみると、第1主成分で約73%、第2主成分で23 %と、2つの主成分で特徴をほぼ説明しきっている（第3、第4主成分の寄与はほとんど無視できる）。

第1主成分の各要素の符号を見てみる。萼の長さ、花弁の長さと幅は同程度でプラス方向に効いていて、萼の幅はマイナス方向の効果を持っている。このことから、萼の細長さと花弁の全体的な大きさによって、アヤメの花が特徴づけられていると考えられる。また第2主成分は、萼の幅で殆ど特徴が決まっている。

可視化

3次元

4つの主成分のうち3つについて3次元で可視化してみると、3つのアヤメの種類がかなりきれいに分離されているのがわかる。

fig, ax = plt.subplots(subplot_kw=dict(projection='3d'))

X0 = X_trans[y==0]
X1 = X_trans[y==1]
X2 = X_trans[y==2]
ax.scatter(X0[:, 0], X0[:, 1], X0[:, 2])
ax.scatter(X1[:, 0], X1[:, 1], X1[:, 2])
ax.scatter(X2[:, 0], X2[:, 1], X2[:, 2])

ax.set_xlabel("Component-1")
ax.set_ylabel("Component-2")
ax.set_zlabel("Component-3")

plt.show()

fig, ax = plt.subplots(subplot_kw=dict(projection='3d'))

X0 = X_trans[y==0]

X1 = X_trans[y==1]

X2 = X_trans[y==2]

ax.scatter(X0[:, 0], X0[:, 1], X0[:, 2])

ax.scatter(X1[:, 0], X1[:, 1], X1[:, 2])

ax.scatter(X2[:, 0], X2[:, 1], X2[:, 2])

ax.set_xlabel("Component-1")

ax.set_ylabel("Component-2")

ax.set_zlabel("Component-3")

plt.show()

2次元

第2主成分まででほとんどの特徴を説明できそうなので、2次元の散布図で表示してみる。

実際、2つの主成分だけでかなりきれいに3つのクラスが分かれている。少し重なっている部分があるが、先の主成分を3つの3次元グラフで傾きを調整すると、より明確にクラスが分けられる。

pca = PCA(n_components=2)
X_trans = pca.fit_transform(X_scaled)

print(pd.DataFrame(pca.components_, columns=feature_names))
print(pca.explained_variance_ratio_)

fig, ax = plt.subplots()

X0 = X_trans[y==0]
X1 = X_trans[y==1]
X2 = X_trans[y==2]
ax.scatter(X0[:, 0], X0[:, 1])
ax.scatter(X1[:, 0], X1[:, 1])
ax.scatter(X2[:, 0], X2[:, 1])

ax.set_xlabel("Component-1")
ax.set_ylabel("Component-2")

plt.show()

pca = PCA(n_components=2)

X_trans = pca.fit_transform(X_scaled)

print(pd.DataFrame(pca.components_, columns=feature_names))

print(pca.explained_variance_ratio_)

fig, ax = plt.subplots()

X0 = X_trans[y==0]

X1 = X_trans[y==1]

X2 = X_trans[y==2]

ax.scatter(X0[:, 0], X0[:, 1])

ax.scatter(X1[:, 0], X1[:, 1])

ax.scatter(X2[:, 0], X2[:, 1])

ax.set_xlabel("Component-1")

ax.set_ylabel("Component-2")

plt.show()

なお今回の計算では、PCAのモデルインスタンス生成時にn_components=2としている。その結果は以下の通りで、1つ前の結果と同じ値になっている。

   sepal length (cm)  sepal width (cm)  petal length (cm)  petal width (cm)
0           0.521066         -0.269347           0.580413          0.564857
1           0.377418          0.923296           0.024492          0.066942
[0.72962445 0.22850762]

sepal length (cm) sepal width (cm) petal length (cm) petal width (cm)

0 0.521066 -0.269347 0.580413 0.564857

1 0.377418 0.923296 0.024492 0.066942

[0.72962445 0.22850762]

主成分分析の特徴

IrisデータセットへのPCAの適用結果から、以下のようにまとめられる。

主成分分析の計算において、ターゲットのクラス分類は全く用いていない（特徴量データのみを用いている）
ターゲットのクラス分類は、散布図を描くときの色分けにのみ利用している
それにも関わらず、散布図において3つのクラスがかなりきれいに分離されている
特徴量の線形和に沿った分散の最大化、という問題設定で、その背後にあるアヤメの種類がうまく分類されている

主成分分析の定式化

2020-11-20 / tau / コメントする

概要

主成分分析では、複数の特徴量を持つデータセットから、そのデータセットの特徴を最もよく表す特徴量軸を発見していく。

ここで「特徴を最もよく表す」ことを数学的に「最も分散が大きくなる」と定義する。そして、分散が最も大きくなるような方向を探すことを目的とする。

ある軸に沿った分散が大きくなるということは、その軸に沿った性質のバリエーションが多いことになる。逆に分散が小さい場合は、その性質を表す数量によっては各データの特徴の違いが判別しにくい。

主成分分析では、分散が最大となるような軸の方向を発見することが目的となる。この軸は元の特徴量の線形和で表現されるもので、各特徴量の係数は、それぞれの特徴量の寄与を表す。

(1) $\begin{align*} \boldsymbol{v} &= a_1 \boldsymbol{x}_1 + \ldosts + a_m \boldsymbol{x_m} \\ &= a_1 \left( \begin{array}{c} \x_1 \\ 0 \\ \vdots \\ 0 \end{array} \right) + \cdots + a_m \left( \begin{array}{c} 0 \\ \vdots \\ 0 \\ x_m \end{array} \right) \\ v &= | \boldsymbol{v} | = a_1 x_1 + \cdots + a_m x_m \end{align*}$

以後、複数の特徴量を持つデータを、特徴量を成分とするベクトルでx表し、多数のベクトルデータx_iがデータセットを構成しているとする。

最大化すべき分散の導出

多数のデータの中のデータiが空間内の点に対応し、その位置ベクトルをx_iであるとする。このx_iの成分が特徴量に対応する。長さが1のあるベクトルdが与えられたとき、x_iのdへの射影の長さは以下のように計算される。

(2) $\begin{align*} x_{i | \boldsymbol{d}} = {\boldsymbol{x}_i}^T \boldsymbol{d} = \boldsymbol{d}^T \boldsymbol{x}_i \quad (| \boldsymbol{d} | = 1) \end{align*}$

たとえば特徴量が2つなら、2次元で以下のような計算になる。

(3) $\begin{align*} \boldsymbol{x}_i = \left( \begin{array}{C} x_{i1} \\ x_{i2} \end{array} \right) , \quad \boldsymbol{d} = \left( \begin{array}{C} d_1 \\ d_2 \end{array} \right) \end{align*}$

(4) $\begin{align*} x_{i | \boldsymbol{d}} = ( d_1 , d_2 ) \left( \begin{array}{c} x_{i1} \\ x_{i2} \end{array} \right) = ( d_1 x_{i1} + d_2 x_{i2} ) \end{align*}$

n個のデータ(i = 1～n)について、射影の平均は以下のように計算される。これは全データのベクトルdの方向に沿った値の平均となる。

(5) $\begin{align*} E( x_{i | \boldsymbol{d}} ) = E \left( \boldsymbol{d}^T \boldsymbol{x}_i \right) = \boldsymbol{d}^T E \left( \boldsymbol{x}_i \right) = \boldsymbol{d}^T \boldsymbol{\mu}_i \end{align*}$

これも2次元の場合で確認すると以下の通り。

(6) $\begin{align*} E(x_{i | \boldsymbol{d}} ) &= E\left[ (d1, d2) \left( \begin{array}{c} x_{i1} \\ x_{i2} \end{array} \right) \right] = (d_1, d_2) \left( \begin{array}{c} E(x_{i1}) \\ E(x_{i2}) \end{array} \right) \\ &= (d_1, d_2) \left( \begin{array}{c} \mu_{i1} \\ \mu_{i2} \end{aray} \right) \end{align*}$

式(5)を使ってベクトルdの方向に沿ったデータの分散を計算する。

(7) $\begin{align*} V( x_{i | \boldsymbol{d}} ) &= V \left( \boldsymbol{d}^T \boldsymbol{x}_i \right) \\ &= E \left[ \left( \boldsymbol{d}^T \boldsymbol{x}_i - E \left( \boldsymbol{d}^T \boldsymbol{x}_i \right) \right)^2 \right] \\ &= E \left[ \left( {\boldsymbol{d}}^T \left( \boldsymbol{x}_i - E(\boldsymbol{x}_i) \right) \right)^2 \right] \\ &= E \left[ {\boldsymbol{d}}^T (\boldsymbol{x}_i - \boldsymbol{\mu}_i ) (\boldsymbol{x}_i - \boldsymbol{\mu}_i )^T \boldsymbol{d} \right] \\ &= \boldsymbol{d}^T E\left[ (\boldsymbol{x}_i - \boldsymbol{\mu}_i ) (\boldsymbol{x}_i - \boldsymbol{\mu}_i )^T \right] \boldsymbol{d} \\ &= \boldsymbol{d}^T \boldsymbol{\Sigma} \boldsymbol{d} \end{align*}$

中央の平均の項が共分散行列Σとなっていることに留意。これより、あるベクトルが与えられたとき、その方向に沿った全データの成分の分散が、そのベクトルと元のデータの共分散行列を使って求めることができる。

こちらを2次元で確認すると以下の通り。

(8) $\begin{align*} &E\left[ (\boldsymbol{x}_i - \boldsymbol{\mu}_i ) (\boldsymbol{x}_i - \boldsymbol{\mu}_i )^T \right] \\ &= E \left[ \left( \begin{array}{c} x_{i1} - \mu_1 \\ x_{i2} - \mu_2 \end{array} \right) (x_{i1} - \mu_1, x_{i2} - \mu_2) \right] \\ &= \left[ \begin{array}{cc} (x_{i1} - \mu_1)^2 & (x_{i1} - \mu_1)(x_{i2} - \mu_2) \\ (x_{i2} - \mu_2)(x_{i1} - \mu_1) & (x_{i2} - \mu_2)^2 \end{array} \right] \end{align*}$

分散の最大化

式(8)で計算された分散が最大となるようにベクトルdの方向を決定する。このとき、dの大きさが1であるという制約条件があるため、問題は制約条件付きの最大化問題となる。

(9) $\begin{gather*} {\rm max} \quad \boldsymbol{d}^T \boldsymbol{\Sigma} \boldsymbol{d} \quad \rm{s.t.} \; | \boldsymbol{d} | = 1 \end{gather*}$

これをLagrangeの未定乗数法で解いていく。。

(10) $\begin{gather*} L( \boldsymbol{d}, \lambda ) = \boldsymbol{d}^T \boldsymbol{\Sigma} \boldsymbol{d} - \lambda (|\boldsymbol{d}|^2 - 1) = 0 \\ \frac{\partial L}{\partial d_i} = 0 \quad ( {\rm for \; all} \; i ) \end{gather*}$

Lagrange関数の第1項については、

(11) $\begin{align*} \boldsymbol{d}^T \boldsymbol{\Sigma} &= \begin{array}{ccc} ( & d_1 V_1 + \cdots + d_n C_{n1} & , \\ & \vdots & ,\\ & d_1 C_{1j} + \cdots + d_n C_{n, j} & , \\ & \vdots & ,\\ & d_1 C_{1n} + \cdots + d_n V_n & ) \end{array} \end{align*}$

より、以下のような長い式になる。

(12) $\begin{align*} \boldsymbol{d}^T \boldsymbol{\Sigma d} &= \begin{array}{c} {d_1}^2 V_1 + \cdots + d_j d_1 C_{j1} + \cdots + d_n d_1 C_{n1} + \\ \vdots \\ d_1 d_j C_{1j} + \cdots + {d_j}^2 V_j + \cdots + d_n d_j C_{nj} + \\ \vdots \\ d_1 d_n C_{1n} + \cdots + d_j d_n C_{jn} + \cdots + {d_n}^2 V_n \end{array} \end{align*}$

また第2項の括弧の中については以下のようになる。

(13) $\begin{align*} | \boldsymbol{d} |^2 - 1 = ( {d_1}^2 + \cdots + {d_j}^2 + \cdots + {d_n}^2 ) -1 \end{align*}$

これらを前提に、Lをd_jで微分すると以下のようになる。

(14) $\begin{align*} 2 d_1 C_{1j} + \cdots + 2 d_j V_{j} + \cdots 2 d_n C_{2j} - 2 \lambda d_j = 0 \end{align*}$

全てのd_jについて考慮した連立方程式を行列形式で表すと以下のようになる。

(15) $\begin{gather*} \boldsymbol{\Sigma d} = \lambda \boldsymbol{d} \\ | \boldsymbol{d} | = 1 \end{gather*}$

1つ目の式は共分散行列に関する固有値問題の式で、d_i (i=1～n)とλのn+1個の変数に対してn個の式となる。これに先ほど脇に置いていたdの大きさに関する制約式を加えて式の数もn+1個となり、dとλが求められる。

特徴量が2つの場合

特徴量が2つの場合を考え、以下のように記号を定義する。

(16) $\begin{align*} \boldsymbol{\Sigma} = \left( \begin{array}{cc} \sigma_{11} & \sigma_{12} \\ \sigma_{21} & \sigma_{22} \end{array} \right) , \quad \boldsymbol{d} = \left( \begin{array}{c} d_1 & d_2 \end{array} \right) \end{align*}$

このとき、分散を最大化する方向の単位ベクトルdを求める方程式は以下のようになる。

(17) $\begin{equation*} \left\{ \begin{array}{l} \left( \begin{array}{cc} \sigma_{11} & \sigma_{12} \\ \sigma_{21} & \sigma_{22} \end{array} \right) \left( \begin{array}{c} d_1 & d_2 \end{array} \right) = \lambda \left( \begin{array}{c} d_1 & d_2 \end{array} \right) \\ {d_1}^2 + {d_2}^2 = 1 \end{array} \right. \end{equation*}$

1つ目の式を解くと、

(18) $\begin{equation*} \left\{ \begin{array}{l} \sigma_{11} d_1 + \sigma_{12} d_2 = \lambda d_1 \\ \sigma_{21} d_1 + \sigma_{22} d_2 = \lambda d_2 \end{array} \rhight. \end{equation*}$

この方程式は不定なのでd₁、d₂それぞれは求められないが、μ = d₂/d₁は計算できる。これは固有ベクトルの方向が定まる。具体的には下記の通り。

(19) $\begin{gather*} \left\{ \begin{array}{l} \sigma_{11} + \sigma_{12} \mu = \lambda \\ \sigma_{21} + \sigma_{22} \mu = \lambda \mu \end{array} \right. \\ \lambda = \sigma_{11} + \sigma_{12} \mu = \frac{\sigma_{21}}{\mu} + \sigma_{22} \\ \sigma_{12} \mu^2 + ( \sigma_{11} - \sigma_{22} ) \mu - \sigma_{21} = 0 \end{gather*}$

これを解いてベクトルdの方向が定まる。これに制約条件|d|² ＝１を加味することで、大きさ1の単位ベクトルとしてdが決定される。

この解き方は最大化問題ではないので、連立方程式から2つの固有ベクトルと固有値が求まる。

第2主成分以降

一般的な固有値問題では、元の変数と同じ数の固有ベクトルと固有値のセットが求まるが、最大化問題として解いた場合には主成分が1つだけ求まる。

scikit-learnのPCAインスタンス生成時にn_componentsで主成分の数に制約をかけることができるが、このことから、PCA.fit()の実行時には連立方程式を解いているのではなく、最大化問題で1つずつ主成分を計算しているのではないかと思われる。

第2主成分以降の計算についての紹介はあまり見られないが、以下の手順と考えらえれる。

各データについて、第1主成分の方向への射影を計算
その射影の符号を逆にしたベクトルを各データに加える
これで第1主成分に沿ったばらつきが全てゼロになるので、残りの成分の中で最大となるベクトルの方向を計算し、第2主成分とする
以上を繰り返し、順次最大主成分沿いの情報を消しながら、各主成分を計算

主成分の意味

主成分の意味の一つに、元のデータを構成する成分という捉え方がある。

たとえば特徴量の数がnである元データXがあり、主成分の数をm(<= n)でモデルを構築するとする。scikit-learnでPCAのインスタンスを生成するのにn_components=mと指定し、fit(X)を実行すると、m個の主成分が生成される。この主成分は共分散行列に対する固有ベクトルであり、要素数n個（特徴量数に等しい）の1次元配列がm行（主成分の数に等しい）並んだ2次元配列として、PCAインスタンスのプロパティーcomponents_に保存される。

(20) $\begin{equation*} \tt{components\_} = \left[ \begin{array}{ccc} (p_{0, 0} & \cdots & p_{0, n-1} ) \\ & \vdots &\\ (p_{m-1, 0} & \cdots & p_{m-1, n-1}) \end{array} \right] = \left[ \begin{array}{c} \boldsymbol{p}_0 \\ \vdots \\ \boldsymbol{p}_m \end{array} \right] \end{equation*}$

元のデータは、各主成分（固有ベクトル）の重み付き和として表現される。

(21) $\begin{equation*} \boldsymbol{x} = (x_0, ..., x_n) = a_0 \boldsymbol{p}_0 + a_1 \boldsymbol{p}_1 + a_2 \boldsymbol{p}_2 + \cdots \end{equation*}$

この様子を2次元で示したのが以下の図で、直行する2つの主成分から元データの１つxが定まる。

xの主成分1、2の方向の大きさはxの各主成分に対する射影で、それらの長さはxと各主成分の内積で得られる。

(22) $\begin{equation*} \boldsymbol{x} = (x_0, ..., x_n) = ( \boldsymbol{x} \cdot \boldsymbol{p}_0 ) \boldsymbol{p}_0 + ( \boldsymbol{x} \cdot \boldsymbol{p}_1 ) \boldsymbol{p}_1 + ( \boldsymbol{x} \cdot \boldsymbol{p}_2 ) \boldsymbol{p}_2 + \cdots \end{equation*}$

DataFrame – get_dummies – One-hot

2020-11-19 / tau / コメントする

概要

DataFrameのget_dummies()メソッドは、属性データ（カテゴリーデータ）をone-hot-encodingの形に変換してくれる。

Scikit-learnにもOneHotEncoderがあるが、get_dumies()はデータの切り貼りをせずにダイレクトに属性変数だけをone-hotの形にしてくれるので便利。

基本

get_dumies()の引数にDataFrameを指定すると、文字列で属性指定されたデータが自動で認識されてon-hotの形に変換される。

import pandas as pd

data = [
    ["male", 1200],
    ["male", 1000],
    ["male", 800],
    ["female", 1500],
    ["female", 800]
]

df = pd.DataFrame(data, columns=["gender", "price"])
print(df)
#    gender  price
# 0    male   1200
# 1    male   1000
# 2    male    800
# 3  female   1500
# 4  female    800

print(pd.get_dummies(df))
#    price  gender_female  gender_male
# 0   1200              0            1
# 1   1000              0            1
# 2    800              0            1
# 3   1500              1            0
# 4    800              1            0

import pandas as pd

data = [

["male", 1200],

["male", 1000],

["male", 800],

["female", 1500],

["female", 800]

]

df = pd.DataFrame(data, columns=["gender", "price"])

print(df)

# gender price

# 0 male 1200

# 1 male 1000

# 2 male 800

# 3 female 1500

# 4 female 800

print(pd.get_dummies(df))

# price gender_female gender_male

# 0 1200 0 1

# 1 1000 0 1

# 2 800 0 1

# 3 1500 1 0

# 4 800 1 0

分解された列名は、"元の列名_属性名"となり、それぞれに対応する属性の列のみが1、その他の列は0となる。列の並びは、属性名の辞書順。数値データの列は無視される。

属性データが複数列の場合

文字列の属性データが複数列ある場合も、自動的にone-hotに分解してくれる。

data = [
    ["male", "Tokyo", 10000],
    ["male", "Tokyo", 70000],
    ["female", "Tokyo", 80000],
    ["male", "Nagoya", 7000],
    ["female", "Nagoya", 8000],
    ["female", "Nagoya", 6000],
    ["male", "Osaka", 9000],
    ["female", "Osaka", 79000],
]

df = pd.DataFrame(data, columns=["gender", "city", "price"])
print(df)
#    gender    city  price
# 0    male   Tokyo  10000
# 1    male   Tokyo  70000
# 2  female   Tokyo  80000
# 3    male  Nagoya   7000
# 4  female  Nagoya   8000
# 5  female  Nagoya   6000
# 6    male   Osaka   9000
7  female   Osaka  79000

print(pd.get_dummies(df))
#    price  gender_female  gender_male  city_Nagoya  city_Osaka  city_Tokyo
# 0  10000              0            1            0           0           1
# 1  70000              0            1            0           0           1
# 2  80000              1            0            0           0           1
# 3   7000              0            1            1           0           0
# 4   8000              1            0            1           0           0
# 5   6000              1            0            1           0           0
# 6   9000              0            1            0           1           0
# 7  79000              1            0            0           1           0

data = [

["male", "Tokyo", 10000],

["male", "Tokyo", 70000],

["female", "Tokyo", 80000],

["male", "Nagoya", 7000],

["female", "Nagoya", 8000],

["female", "Nagoya", 6000],

["male", "Osaka", 9000],

["female", "Osaka", 79000],

]

df = pd.DataFrame(data, columns=["gender", "city", "price"])

print(df)

# gender city price

# 0 male Tokyo 10000

# 1 male Tokyo 70000

# 2 female Tokyo 80000

# 3 male Nagoya 7000

# 4 female Nagoya 8000

# 5 female Nagoya 6000

# 6 male Osaka 9000

7 female Osaka 79000

print(pd.get_dummies(df))

# price gender_female gender_male city_Nagoya city_Osaka city_Tokyo

# 0 10000 0 1 0 0 1

# 1 70000 0 1 0 0 1

# 2 80000 1 0 0 0 1

# 3 7000 0 1 1 0 0

# 4 8000 1 0 1 0 0

# 5 6000 1 0 1 0 0

# 6 9000 0 1 0 1 0

# 7 79000 1 0 0 1 0

属性が数値表現の場合

属性値が文字列ではなく数値表現の場合、get_dummies()の引数に単にDataFrameを渡すだけでは変換してくれない（通常の数量データとして認識される）。

そこで、変換したい列をcolumns引数で指定する。

data = [
    [1, 10, 100],
    [1, 20, 100],
    [2, 10, 200],
    [1, 20, 200],
    [1, 30, 100],
]

df = pd.DataFrame(data, columns=["id1", "id2", "value"])
print(df)
#    id1  id2  value
# 0    1   10    100
# 1    1   20    100
# 2    2   10    200
# 3    1   20    200
# 4    1   30    100

print(pd.get_dummies(df, columns=["id1"]))
#    id2  value  id1_1  id1_2
# 0   10    100      1      0
# 1   20    100      1      0
# 2   10    200      0      1
# 3   20    200      1      0
# 4   30    100      1      0

data = [

[1, 10, 100],

[1, 20, 100],

[2, 10, 200],

[1, 20, 200],

[1, 30, 100],

]

df = pd.DataFrame(data, columns=["id1", "id2", "value"])

print(df)

# id1 id2 value

# 0 1 10 100

# 1 1 20 100

# 2 2 10 200

# 3 1 20 200

# 4 1 30 100

print(pd.get_dummies(df, columns=["id1"]))

# id2 value id1_1 id1_2

# 0 10 100 1 0

# 1 20 100 1 0

# 2 10 200 0 1

# 3 20 200 1 0

# 4 30 100 1 0

複数の属性データの列がある場合、columns引数でリスト指定する。

print(pd.get_dummies(df, columns=["id1", "id2"]))

#    value  id1_1  id1_2  id2_10  id2_20  id2_30
# 0    100      1      0       1       0       0
# 1    100      1      0       0       1       0
# 2    200      0      1       1       0       0
# 3    200      1      0       0       1       0
# 4    100      1      0       0       0       1

print(pd.get_dummies(df, columns=["id1", "id2"]))

# value id1_1 id1_2 id2_10 id2_20 id2_30

# 0 100 1 0 1 0 0

# 1 100 1 0 0 1 0

# 2 200 0 1 1 0 0

# 3 200 1 0 0 1 0

# 4 100 1 0 0 0 1

属性名の指定

prefix引数で文字列を指定すると、属性名がその文字列で置き換えられる。ただし複数の属性列が全て同じ文字列になる。

print(pd.get_dummies(df, columns=["id1", "id2"], prefix="attr"))

#    value  attr_1  attr_2  attr_10  attr_20  attr_30
# 0    100       1       0        1        0        0
# 1    100       1       0        0        1        0
# 2    200       0       1        1        0        0
# 3    200       1       0        0        1        0
# 4    100       1       0        0        0        1

print(pd.get_dummies(df, columns=["id1", "id2"], prefix="attr"))

# value attr_1 attr_2 attr_10 attr_20 attr_30

# 0 100 1 0 1 0 0

# 1 100 1 0 0 1 0

# 2 200 0 1 1 0 0

# 3 200 1 0 0 1 0

# 4 100 1 0 0 0 1

属性列ごとにprefixを変えて指定したい場合はリストで指定。

print(pd.get_dummies(df, columns=["id1", "id2"], prefix=["attr", "cat"]))

#    value  attr_1  attr_2  cat_10  cat_20  cat_30
# 0    100       1       0       1       0       0
# 1    100       1       0       0       1       0
# 2    200       0       1       1       0       0
# 3    200       1       0       0       1       0
# 4    100       1       0       0       0       1

print(pd.get_dummies(df, columns=["id1", "id2"], prefix=["attr", "cat"]))

# value attr_1 attr_2 cat_10 cat_20 cat_30

# 0 100 1 0 1 0 0

# 1 100 1 0 0 1 0

# 2 200 0 1 1 0 0

# 3 200 1 0 0 1 0

# 4 100 1 0 0 0 1

DataFrame – データの概観

2020-11-19 / tau / コメントする

概要

DataFrameの規模、格納されているデータの概要や基礎統計量を概観する各種の手順。Scikit-learnのBoston housingデータセットを例にする。

DataFrameの規模・形状

sizeプロパティーで全データ数、shapeプロパティーで行数と列数を確認。

import pandas as pd
from sklearn.datasets import load_boston

ds = load_boston()
df = pd.DataFrame(ds.data, columns=ds.feature_names)

print(df.size)
# 6578

print(df.shape)
# (506, 13)

import pandas as pd

from sklearn.datasets import load_boston

ds = load_boston()

df = pd.DataFrame(ds.data, columns=ds.feature_names)

print(df.size)

# 6578

print(df.shape)

# (506, 13)

データの先頭部分と末尾部分

head()メソッド／tail()メソッドで先頭／末尾の5行分が得られる。引数で抜き出す行数を指定。

print(df.head())
#       CRIM    ZN  INDUS  CHAS    NOX  ...  RAD    TAX  PTRATIO       B  LSTAT
# 0  0.00632  18.0   2.31   0.0  0.538  ...  1.0  296.0     15.3  396.90   4.98
# 1  0.02731   0.0   7.07   0.0  0.469  ...  2.0  242.0     17.8  396.90   9.14
# 2  0.02729   0.0   7.07   0.0  0.469  ...  2.0  242.0     17.8  392.83   4.03
# 3  0.03237   0.0   2.18   0.0  0.458  ...  3.0  222.0     18.7  394.63   2.94
# 4  0.06905   0.0   2.18   0.0  0.458  ...  3.0  222.0     18.7  396.90   5.33
# 
# [5 rows x 13 columns]

print(df.tail(3))
#         CRIM   ZN  INDUS  CHAS    NOX  ...  RAD    TAX  PTRATIO       B  LSTAT
# 503  0.06076  0.0  11.93   0.0  0.573  ...  1.0  273.0     21.0  396.90   5.64
# 504  0.10959  0.0  11.93   0.0  0.573  ...  1.0  273.0     21.0  393.45   6.48
# 505  0.04741  0.0  11.93   0.0  0.573  ...  1.0  273.0     21.0  396.90   7.88
# 
# [3 rows x 13 columns]

print(df.head())

# CRIM ZN INDUS CHAS NOX ... RAD TAX PTRATIO B LSTAT

# 0 0.00632 18.0 2.31 0.0 0.538 ... 1.0 296.0 15.3 396.90 4.98

# 1 0.02731 0.0 7.07 0.0 0.469 ... 2.0 242.0 17.8 396.90 9.14

# 2 0.02729 0.0 7.07 0.0 0.469 ... 2.0 242.0 17.8 392.83 4.03

# 3 0.03237 0.0 2.18 0.0 0.458 ... 3.0 222.0 18.7 394.63 2.94

# 4 0.06905 0.0 2.18 0.0 0.458 ... 3.0 222.0 18.7 396.90 5.33

# [5 rows x 13 columns]

print(df.tail(3))

# CRIM ZN INDUS CHAS NOX ... RAD TAX PTRATIO B LSTAT

# 503 0.06076 0.0 11.93 0.0 0.573 ... 1.0 273.0 21.0 396.90 5.64

# 504 0.10959 0.0 11.93 0.0 0.573 ... 1.0 273.0 21.0 393.45 6.48

# 505 0.04741 0.0 11.93 0.0 0.573 ... 1.0 273.0 21.0 396.90 7.88

# [3 rows x 13 columns]

`info()`～各列の基本情報の表示

info()メソッドは、DataFrameの概要に関する概略情報を出力する。直接標準出力にプリントする点に注意。

df.info()

# <class 'pandas.core.frame.DataFrame'>
# RangeIndex: 506 entries, 0 to 505
# Data columns (total 13 columns):
#  #   Column   Non-Null Count  Dtype  
# ---  ------   --------------  -----  
#  0   CRIM     506 non-null    float64
#  1   ZN       506 non-null    float64
#  2   INDUS    506 non-null    float64
#  3   CHAS     506 non-null    float64
#  4   NOX      506 non-null    float64
#  5   RM       506 non-null    float64
#  6   AGE      506 non-null    float64
#  7   DIS      506 non-null    float64
#  8   RAD      506 non-null    float64
#  9   TAX      506 non-null    float64
#  10  PTRATIO  506 non-null    float64
#  11  B        506 non-null    float64
#  12  LSTAT    506 non-null    float64
# dtypes: float64(13)
# memory usage: 51.5 KB

df.info()

# <class 'pandas.core.frame.DataFrame'>

# RangeIndex: 506 entries, 0 to 505

# Data columns (total 13 columns):

# # Column Non-Null Count Dtype

# --- ------ -------------- -----

# 0 CRIM 506 non-null float64

# 1 ZN 506 non-null float64

# 2 INDUS 506 non-null float64

# 3 CHAS 506 non-null float64

# 4 NOX 506 non-null float64

# 5 RM 506 non-null float64

# 6 AGE 506 non-null float64

# 7 DIS 506 non-null float64

# 8 RAD 506 non-null float64

# 9 TAX 506 non-null float64

# 10 PTRATIO 506 non-null float64

# 11 B 506 non-null float64

# 12 LSTAT 506 non-null float64

# dtypes: float64(13)

# memory usage: 51.5 KB

たとえば一部にNaNが含まれる場合の出力は以下のようになる。

import numpy as np
df1 = pd.DataFrame(np.arange(15).reshape(-1, 3), columns=["a", "b", "c"])
df1.loc[1, "b"] = np.nan
df1.loc[3, "c"] = np.nan

df1.info()

# <class 'pandas.core.frame.DataFrame'>
# RangeIndex: 5 entries, 0 to 4
# Data columns (total 3 columns):
#  #   Column  Non-Null Count  Dtype  
# ---  ------  --------------  -----  
#  0   a       5 non-null      int32  
#  1   b       4 non-null      float64
#  2   c       4 non-null      float64
# dtypes: float64(2), int32(1)
# memory usage: 228.0 bytes

import numpy as np

df1 = pd.DataFrame(np.arange(15).reshape(-1, 3), columns=["a", "b", "c"])

df1.loc[1, "b"] = np.nan

df1.loc[3, "c"] = np.nan

df1.info()

# <class 'pandas.core.frame.DataFrame'>

# RangeIndex: 5 entries, 0 to 4

# Data columns (total 3 columns):

# # Column Non-Null Count Dtype

# --- ------ -------------- -----

# 0 a 5 non-null int32

# 1 b 4 non-null float64

# 2 c 4 non-null float64

# dtypes: float64(2), int32(1)

# memory usage: 228.0 bytes

`dscribe()`～基本的な統計量

describe()メソッドは、各列のデータについて、個数や平均といった基本的な統計量を計算する。

print(df.describe())

#              CRIM          ZN       INDUS  ...     PTRATIO           B       LSTAT
# count  506.000000  506.000000  506.000000  ...  506.000000  506.000000  506.000000
# mean     3.613524   11.363636   11.136779  ...   18.455534  356.674032   12.653063
# std      8.601545   23.322453    6.860353  ...    2.164946   91.294864    7.141062
# min      0.006320    0.000000    0.460000  ...   12.600000    0.320000    1.730000
# 25%      0.082045    0.000000    5.190000  ...   17.400000  375.377500    6.950000
# 50%      0.256510    0.000000    9.690000  ...   19.050000  391.440000   11.360000
# 75%      3.677083   12.500000   18.100000  ...   20.200000  396.225000   16.955000
# max     88.976200  100.000000   27.740000  ...   22.000000  396.900000   37.970000

print(df.describe())

# CRIM ZN INDUS ... PTRATIO B LSTAT

# count 506.000000 506.000000 506.000000 ... 506.000000 506.000000 506.000000

# mean 3.613524 11.363636 11.136779 ... 18.455534 356.674032 12.653063

# std 8.601545 23.322453 6.860353 ... 2.164946 91.294864 7.141062

# min 0.006320 0.000000 0.460000 ... 12.600000 0.320000 1.730000

# 25% 0.082045 0.000000 5.190000 ... 17.400000 375.377500 6.950000

# 50% 0.256510 0.000000 9.690000 ... 19.050000 391.440000 11.360000

# 75% 3.677083 12.500000 18.100000 ... 20.200000 396.225000 16.955000

# max 88.976200 100.000000 27.740000 ... 22.000000 396.900000 37.970000

特定の列の統計量を見たいときは列を指定。

なおstd（標準偏差）については、ddof=1を指定した結果と同じであり、n−1で割った不偏分散。

print(df["CRIM"].describe())
# count    506.000000
# mean       3.613524
# std        8.601545
# min        0.006320
# 25%        0.082045
# 50%        0.256510
# 75%        3.677083
# max       88.976200
# Name: CRIM, dtype: float64

print(df["CRIM"].std(ddof=1))
# 8.60154510533249

print(df["CRIM"].describe())

# count 506.000000

# mean 3.613524

# std 8.601545

# min 0.006320

# 25% 0.082045

# 50% 0.256510

# 75% 3.677083

# max 88.976200

# Name: CRIM, dtype: float64

print(df["CRIM"].std(ddof=1))

# 8.60154510533249

全ての列を見たいときには、set_option()メソッドの引数でdisplay.max_columnsパラメーターを指定する。初期値に戻すときはreset_option()メソッド。

pd.set_option('display.max_columns', 100)
print(df.describe())
pd.reset_option('display.max_columns')

#              CRIM          ZN       INDUS        CHAS         NOX          RM  \
# count  506.000000  506.000000  506.000000  506.000000  506.000000  506.000000   
# mean     3.613524   11.363636   11.136779    0.069170    0.554695    6.284634   
# std      8.601545   23.322453    6.860353    0.253994    0.115878    0.702617   
# min      0.006320    0.000000    0.460000    0.000000    0.385000    3.561000   
# 25%      0.082045    0.000000    5.190000    0.000000    0.449000    5.885500   
# 50%      0.256510    0.000000    9.690000    0.000000    0.538000    6.208500   
# 75%      3.677083   12.500000   18.100000    0.000000    0.624000    6.623500   
# max     88.976200  100.000000   27.740000    1.000000    0.871000    8.780000   
# 
#               AGE         DIS         RAD         TAX     PTRATIO           B  \
# count  506.000000  506.000000  506.000000  506.000000  506.000000  506.000000   
# mean    68.574901    3.795043    9.549407  408.237154   18.455534  356.674032   
# std     28.148861    2.105710    8.707259  168.537116    2.164946   91.294864   
# min      2.900000    1.129600    1.000000  187.000000   12.600000    0.320000   
# 25%     45.025000    2.100175    4.000000  279.000000   17.400000  375.377500   
# 50%     77.500000    3.207450    5.000000  330.000000   19.050000  391.440000   
# 75%     94.075000    5.188425   24.000000  666.000000   20.200000  396.225000   
# max    100.000000   12.126500   24.000000  711.000000   22.000000  396.900000   
# 
#             LSTAT  
# count  506.000000  
# mean    12.653063  
# std      7.141062  
# min      1.730000  
# 25%      6.950000  
# 50%     11.360000
# 75%     16.955000
# max     37.970000

pd.set_option('display.max_columns', 100)

print(df.describe())

pd.reset_option('display.max_columns')

# CRIM ZN INDUS CHAS NOX RM \

# count 506.000000 506.000000 506.000000 506.000000 506.000000 506.000000

# mean 3.613524 11.363636 11.136779 0.069170 0.554695 6.284634

# std 8.601545 23.322453 6.860353 0.253994 0.115878 0.702617

# min 0.006320 0.000000 0.460000 0.000000 0.385000 3.561000

# 25% 0.082045 0.000000 5.190000 0.000000 0.449000 5.885500

# 50% 0.256510 0.000000 9.690000 0.000000 0.538000 6.208500

# 75% 3.677083 12.500000 18.100000 0.000000 0.624000 6.623500

# max 88.976200 100.000000 27.740000 1.000000 0.871000 8.780000

# AGE DIS RAD TAX PTRATIO B \

# count 506.000000 506.000000 506.000000 506.000000 506.000000 506.000000

# mean 68.574901 3.795043 9.549407 408.237154 18.455534 356.674032

# std 28.148861 2.105710 8.707259 168.537116 2.164946 91.294864

# min 2.900000 1.129600 1.000000 187.000000 12.600000 0.320000

# 25% 45.025000 2.100175 4.000000 279.000000 17.400000 375.377500

# 50% 77.500000 3.207450 5.000000 330.000000 19.050000 391.440000

# 75% 94.075000 5.188425 24.000000 666.000000 20.200000 396.225000

# max 100.000000 12.126500 24.000000 711.000000 22.000000 396.900000

# LSTAT

# count 506.000000

# mean 12.653063

# std 7.141062

# min 1.730000

# 25% 6.950000

# 50% 11.360000

# 75% 16.955000

# max 37.970000

属性変数のカウント

属性変数の属性値のカウントには、value_counts()メソッドを使う。このメソッドは、ユニークな値の数をカウントして集計する。

print(df["CHAS"].value_counts())
# 0.0    471
# 1.0     35
# Name: CHAS, dtype: int64

print(df["RAD"].value_counts())
# 24.0    132
# 5.0     115
# 4.0     110
# 3.0      38
# 6.0      26
# 8.0      24
# 2.0      24
# 1.0      20
# 7.0      17
# Name: RAD, dtype: int64

print(df["CHAS"].value_counts())

# 0.0 471

# 1.0 35

# Name: CHAS, dtype: int64

print(df["RAD"].value_counts())

# 24.0 132

# 5.0 115

# 4.0 110

# 3.0 38

# 6.0 26

# 8.0 24

# 2.0 24

# 1.0 20

# 7.0 17

# Name: RAD, dtype: int64

2つの属性変数をファンシーインデックスで指定すると、「2つの属性のユニークな組み合わせ」の数が集計される。

print(df[["CHAS", "RAD"]].value_counts())

# CHAS  RAD 
# 0.0   24.0    124
#       5.0     104
#       4.0     102
#       3.0      36
#       6.0      26
#       2.0      24
#       8.0      19
#       1.0      19
#       7.0      17
# 1.0   5.0      11
#       24.0      8
#       4.0       8
#       8.0       5
#       3.0       2
#       1.0       1
# dtype: int64

print(df[["CHAS", "RAD"]].value_counts())

# CHAS RAD

# 0.0 24.0 124

# 5.0 104

# 4.0 102

# 3.0 36

# 6.0 26

# 2.0 24

# 8.0 19

# 1.0 19

# 7.0 17

# 1.0 5.0 11

# 24.0 8

# 4.0 8

# 8.0 5

# 3.0 2

# 1.0 1

# dtype: int64

Python – 平均・分散・共分散

2020-11-17 / tau / コメントする

1次元配列の場合

平均

平均はmean()メソッドで計算。

import numpy as np

x = np.array([1, 1, 2, 3, 5, 8, 13, 21, 34, 55])

n = x.size
s = x.sum()
m = x.mean()

print(m, s / n)
# 14.3 14.3

import numpy as np

x = np.array([1, 1, 2, 3, 5, 8, 13, 21, 34, 55])

n = x.size

s = x.sum()

m = x.mean()

print(m, s / n)

# 14.3 14.3

分散・標準偏差

標本分散・標本標準偏差

分散はvar()メソッド、標準偏差はstd()メソッドで計算。デフォルトでは標本分散と標本標準偏差。

v = x.var()
sv = x.std()

print(v, np.sum((x - m)**2) / n)
# 285.01000000000005 285.01000000000005

print(sv, np.sqrt(v))
# 16.88223918797504 16.88223918797504

v = x.var()

sv = x.std()

print(v, np.sum((x - m)**2) / n)

# 285.01000000000005 285.01000000000005

print(sv, np.sqrt(v))

# 16.88223918797504 16.88223918797504

不偏推定量

分散の算出では、引数ddof (delta degrees of freedom)で指定した引数をデータ数nから引いた数で偏差の二乗和を割る。デフォルトではddof=0なので標本分散、標本標準偏差になるが、ddof=1を指定すると、不偏分散、標準偏差の不偏推定量が計算される。

u = x.var(ddof=1)
su = x.std(ddof=1)

print(u, np.sum((x - m)**2) / (n - 1))
# 316.6777777777778 316.6777777777778

print(su, np.sqrt(u))
# 17.795442612584207 17.795442612584207

u = x.var(ddof=1)

su = x.std(ddof=1)

print(u, np.sum((x - m)**2) / (n - 1))

# 316.6777777777778 316.6777777777778

print(su, np.sqrt(u))

# 17.795442612584207 17.795442612584207

2次元配列の場合

平均

平均はmean()メソッドで計算。引数axis=0/1で列方向／行方向の計算方向を指定。計算結果は1次元配列で、要素数は列方向の平均なら行数、行方向の平均なら列数と同じになる。

X = np.array([1, 1, 2, 3, 5, 8, 13, 21, 34, 55]).reshape(-1, 2)

n = X.shape[0]
s = X.sum(axis=0)
m = X.mean(axis=0)

print(m, s / n)
# [11.  17.6] [11.  17.6]

X = np.array([1, 1, 2, 3, 5, 8, 13, 21, 34, 55]).reshape(-1, 2)

n = X.shape[0]

s = X.sum(axis=0)

m = X.mean(axis=0)

print(m, s / n)

# [11. 17.6] [11. 17.6]

分散・標準偏差

標本分散・標本標準偏差

分散はvar()メソッド、標準偏差はstd()メソッドで計算。デフォルトでは標本分散、標本標準偏差を計算。デフォルトでは標本分散、標本標準偏差。

v = X.var(axis=0)

print(v, np.sum((X - m)**2, axis=0) / n)
# [150.   398.24] [150.   398.24]

v = X.var(axis=0)

print(v, np.sum((X - m)**2, axis=0) / n)

# [150. 398.24] [150. 398.24]

不偏推定量

引数ddof=1を指定すると、データ数nからddof=1を引いた数で偏差の二乗和が割られ、不偏分散、標準偏差の不偏推定量が計算される。

u = X.var(axis=0, ddof=1)

print(u, np.sum((X - m)**2, axis=0) / (n - 1))
# [187.5 497.8] [187.5 497.8]

u = X.var(axis=0, ddof=1)

print(u, np.sum((X - m)**2, axis=0) / (n - 1))

# [187.5 497.8] [187.5 497.8]

共分散

cov()メソッドで二つのデータの共分散行列を計算する。結果は2次元配列で得られ、対角要素は各データの分散、対角要素以外の要素は行数・列数に対応する共分散。

デフォルトでは引数bias=Falseとなっていて、偏差の積和をn−1で割った不偏推定量が計算される。

print(np.cov(X[:, 0], X[:, 1]))
# [[187.5 305.5]
#  [305.5 497.8]]

Xd = X - m
print(np.sum(Xd[:, 0] * Xd[:, 1]) / (n - 1))
# 305.5

print(np.cov(X[:, 0], X[:, 1]))

# [[187.5 305.5]

# [305.5 497.8]]

Xd = X - m

print(np.sum(Xd[:, 0] * Xd[:, 1]) / (n - 1))

# 305.5

bias=Trueを指定すると、偏差の積和を割る数がnとなり、標本に対する分散・共分散が計算される。

print(np.cov(X[:, 0], X[:, 1], bias=True))
# [[150.   244.4 ]
#  [244.4  398.24]]

print(np.sum(Xd[:, 0] * Xd[:, 1]) / n)
# 244.4

print(np.cov(X[:, 0], X[:, 1], bias=True))

# [[150. 244.4 ]

# [244.4 398.24]]

print(np.sum(Xd[:, 0] * Xd[:, 1]) / n)

# 244.4

感度＝陽性的中率の特性

2020-11-16 / tau / コメントする

概要

機械学習のモデルの性能や感染症検査の確からしさを検証する際、陽性的中率（適合度）や陰性的中率を確認すべきだが、これらの値が、そもそものデータの特性やモデル／検査の性能によってどのように変化するかを確認する。

具体的には、注目事象の率と真陽性率（感度）・真陰性率（特異度）を変化させたときの、陽性的中率・陰性的中率の変化を見る。

これらの値の意味や計算方法については、Confusing matrixを参照。

その結果から、以下のようなことがわかった。

予測モデルや検査において、単に感度のみを向上させても適合度（陽性的中率）は大きく変化しない
特異度を向上させることで適合度は大きく向上する
ターゲット比率がとても小さい場合、感度・特異度をかなり大きくしても、適合度は小さな値になる

2020年現在、世界的に大きな影響を及ぼしているCOVID-19（新型コロナウィルス）感染症のPCR検査では、一般に感度が70％程度、特異度が90％以上、陽性的中率が数%程度という値が多い。感度が7割程度というのは少し低く、陽性的中率がそもそも小さすぎるという気がしていたが、上記のことと符合することがわかった。

指標

以下の指標を、目的として計算する指標とする。

PPV(Positive Predicted Value)：陽性的中率、適合度、Precision
NPV(Negative Predicted Value)：陰性的中率

これらの指標を計算するために用いる指標は以下の通り。

TR(Target Rate)：注目事象の全体比率（ターゲット比率）
TPR(True Positive Rate)：真陽性率、感度(Sencitivity)
TNR(True Negative Rate)：真陰性率、特異度(Specificity)

例えば感染症の例で言うと、有病率(TR)、検査の感度(TPR)、特異度(TNR)がわかっているときに、陽性的中率(PPV)、陰性的中率(NPV)を求めることに相当する。

PPV・NPVの計算式の導出

元データの構成

まず、confusing matrixを以下のように表現する。これは、データ数で表現されたテーブルの各要素を全データ数で割った率で表すことに相当する。

$\begin{align*} \begin{array}{cc|cc|c} & & \mathrm{Prediction}\\ & & \mathrm{Positive} & \mathrm{Negative} & \mathrm{Sum} \\ \hline \mathrm{Fact} & \mathrm{Positive} & tp & fn & r_1 \\ & \mathrm{Negative} & fp & tn & r_2 \\ \hline & \mathrm{Sum}& c_1 & c_2 & 1 \end{array} \end{align}$

PPV・NPVの計算式

まず、事実(Fact)がpositiveである率がr₁に相当し、これはTR (target rate)に等しい。このTRと率TPRを使って、Positiveの行のtp(true positive)とfn (false negative)の率を計算。

(1) $\begin{align*} r_1 &= TR \\ tp &= r_1 \cdot TPR = TR \cdot TPR \\ fn &= r_1 \cdot (1 - TPR) = TR (1 - TPR) \end{align*}$

2行目の合計r₂については、行和の合計が1になることから以下のように計算される。

(2) $\begin{align*} r_2 &= 1 - r_1 = 1 - TR \end{align*}$

このr₂と率TNRからNegativeの行のtn(true negative)とfp (false positive)を計算。

(3) $\begin{align*} tn &= r_2 \cdot TNR = (1 - TR) TNR \\ fp &= r_2 (1 - TNR) = (1 - TR) (1 - TNR) \end{align*}$

tpとfpからc₁を、tnとfnからc₂を計算。

(4) $\begin{align*} c_1 &= tp + fp = TR \cdot TPR + (1 - TR) (1 - TNR) \\ c_2 &= tn + fn = (1 - TR) TNR + TR (1 - TPR) \end{align*}$

PPV（陽性的中率、感度）はc₁に対するtpの率で計算される。以下の式は分数の分数で若干ややこしいが、3つの指標が1回ずつ現れ、整った形になる。

(5) $\begin{align*} PPV &= \frac{tp}{c_1} = \frac{TR \cdot TPR}{TR \cdot TPR + (1 - TR) (1 - TNR) } \\ &= \frac{1}{1 + \left(\dfrac{1}{TR} - 1 \right) \dfrac{1 - TNR}{TPR}} \end{align*}$

NPV（陰性的中率、特異度）はc₂に対するtnの率で計算される。以下の式とPPVの式を比べると、はTRの分数項ついて逆数であり、TPRとTNRが入れ替わっていて、PPVとNPVで対称性がある。

(6) $\begin{align*} NPV &= \frac{tn}{c_2} = \frac{(1 - TR) TNR}{(1 - TR) TNR + TR (1 - TPR)}\\ &= \frac{1}{1 + \dfrac{TR}{1 - TR} \dfrac{1 - TPR}{TNR}} \end{align*}$

パラメーターに応じたPPV・NPVの変化

PPV

上記の結果を用いて、ターゲット比率、真陽性率（感度）、真陰性率（特異度）の様々な値に対するPPV（陽性的中率）、NPV（陰性的中率の変化を観察する。

まず、ターゲット比率が1に近い（ほとんどがターゲットとなるような）状態から、ターゲットが0に近いような（ターゲットとなるデータがほとんどないような）状態の間で、PPVがどのように変化するか確認してみる。

TPR（感度）の値によって曲線の形に若干の変化はあるがあまり大きくは変わらず、むしろTNR（特異度）の値による曲線の形状の変化が大きい。ここでTRが0.1～0と小さい範囲のところを見てみる。

やはり感度の影響はあまり大きくないようである。TNRを大きくするにしたがって曲線の形状は大きく変化し、ターゲット比率が小さいところでの適合度が向上するが、ターゲット比率が0に近いところではPPVがかなり小さくなる。

次に、TPRを変化させたときの曲線の違いが分かるように、表示させる変数を入れ替えてみる。まずTRが1～0の全域。

やはり感度による曲線の変化は小さく、特異度の影響が大きい。以下のようにTRが0.1～0の範囲を拡大しても同様の傾向。

以上の結果から、以下のことが言える。

ターゲット比率が低くなるほどNPVは小さくなる（適合度が低くなり、予測／検査の信頼性が下がる）
予測モデルや検査のTPR（感度）を上げることによるPPVの向上効果はあまり大きくない（いたずらに感度を上げても顕著な効果はない）
TNR（特異度）の向上によって、適合度は大きく向上する
ターゲット比率がとても小さい場合、その率の現象に従って適合度は急激に低下する

さらにこれを一般的な表現でまとめると、

予測モデルや検査において、単に感度のみを向上させても適合度（陽性的中率）は大きく変化しない
特異度を向上させることで適合度は大きく向上する
ターゲット比率がとても小さい場合、感度・特異度をかなり大きくしても、適合度は小さな値になる

NPV

PPVと同様にNPVについても計算してみた。

まずいくつかのTPRに対して、TNRを変化させて曲線を描いたもの。PPVの場合と比べて形状が左右逆で、TNRを固定してTPRを変化させたときの図と同じ傾向。

次に、いくつかのTNRを固定してTPRを変化させたもの。これもPPVと形状、TPR、TNRの関係が逆になっている。

PPVとNPVの関係

PPVとNPVが同じTPR、TNRに対してどのように変化するか重ねてみる。

TPRとTNRを同程度とすることでターゲット比率0.5付近で双方が等しくなり、その値を高くすることで、より広い範囲でPPVが向上する。

シミュレーションによる挙動確認

これまでの結果は、confusion matrixの各要素にTR、TPRなどの比率を適用してPPV、NPVを計算した。この方法は、ある予測／判定が理論通りに再現された場合だが、実際にはターゲットとなる事象の割合も、予測がpositive/negativeになる割合も確率事象である。

そこで念のため、多数の二値(True/False)正解データをランダムに生成し、これに対してTPR、TNRの設定に従った答えを出す疑似的なモデルで「予測」する。その結果を整理したconfusion_matrixからPPVを計算したのが以下の図である。

その結果は計算式による場合と同じで、理論上の挙動と実世界で起こるであろう挙動が一致している。

処理内容は以下の通り。

与えられたTrue/Falseに対して、あらかじめ設定したTPR/TNRと一様乱数に従ってTrue/Falseを「予測」する疑似予測モデルを準備
TR=1～0の間で100個のデータについてPPVを計算する
- 1つのTRについて10万個の2値正解データを生成
- 正解データセットを疑似予測モデルに適用して予測データセットを得る
- 予測データセットからconfusion matrixを構成し、その要素からPPVを計算し、配列に格納
以上の結果をプロット

import numpy as np
import random as rnd
import matplotlib.pyplot as plt
from sklearn.metrics import confusion_matrix

def generate_data(target_rate=0.5, data_size=1000):
    y = np.empty(data_size, dtype=bool)
    for i in range(data_size):
        y[i] = rnd.random() < target_rate
    return y

class ImitationLearningModel:
    def __init__(self, tpr, tnr):
        self.tpr = tpr
        self.tnr = tnr

    def predict(self, y_fact):
        y_pred = np.empty(y_fact.size, dtype=bool)
        for i in range(y_fact.size):
            if y_fact[i]:
                y_pred[i] = rnd.random() < self.tpr
            else:
                y_pred[i] = rnd.random() > self.tnr
        return y_pred

fig, ax = plt.subplots()

tpr = 0.9
tnrs = [0.5, 0.7, 0.9, 0.99]

n_data = 100
trs = np.linspace(0, 1, n_data)
ppv = np.empty((len(tnrs), trs.size))

for j, tr in enumerate(trs):
    y_fact = generate_data(target_rate=tr, data_size=100000)
    for i, tnr in enumerate(tnrs):
        model = ImitationLearningModel(tpr, tnr)
        y_pred = model.predict(y_fact)
        cfmat = confusion_matrix(y_fact, y_pred, labels=[True, False])
        ppv[i, j] = cfmat[0, 0] / (cfmat[0, 0] + cfmat[1, 0])

for i, tnr in enumerate(tnrs):
    ax.plot(trs, ppv[i], label=tnr)

ax.set_xlim(0, 1)
ax.set_ylim(0, 1)
ax.invert_xaxis()
ax.set_xlabel("Target Rate")
ax.set_ylabel("PPV")
ax.legend()
ax.legend(title="TNR", loc='lower left')
ax.set_title("TPR={}".format(tpr))

plt.show()

import numpy as np

import random as rnd

import matplotlib.pyplot as plt

from sklearn.metrics import confusion_matrix

def generate_data(target_rate=0.5, data_size=1000):

y = np.empty(data_size, dtype=bool)

for i in range(data_size):

y[i] = rnd.random() < target_rate

return y

class ImitationLearningModel:

def __init__(self, tpr, tnr):

self.tpr = tpr

self.tnr = tnr

def predict(self, y_fact):

y_pred = np.empty(y_fact.size, dtype=bool)

for i in range(y_fact.size):

if y_fact[i]:

y_pred[i] = rnd.random() < self.tpr

else:

y_pred[i] = rnd.random() > self.tnr

return y_pred

fig, ax = plt.subplots()

tpr = 0.9

tnrs = [0.5, 0.7, 0.9, 0.99]

n_data = 100

trs = np.linspace(0, 1, n_data)

ppv = np.empty((len(tnrs), trs.size))

for j, tr in enumerate(trs):

y_fact = generate_data(target_rate=tr, data_size=100000)

for i, tnr in enumerate(tnrs):

model = ImitationLearningModel(tpr, tnr)

y_pred = model.predict(y_fact)

cfmat = confusion_matrix(y_fact, y_pred, labels=[True, False])

ppv[i, j] = cfmat[0, 0] / (cfmat[0, 0] + cfmat[1, 0])

for i, tnr in enumerate(tnrs):

ax.plot(trs, ppv[i], label=tnr)

ax.set_xlim(0, 1)

ax.set_ylim(0, 1)

ax.invert_xaxis()

ax.set_xlabel("Target Rate")

ax.set_ylabel("PPV")

ax.legend()

ax.legend(title="TNR", loc='lower left')

ax.set_title("TPR={}".format(tpr))

plt.show()

Confusion Matrix～混同行列

2020-11-12 / tau / コメントする

概要

精度が高いのに性能が悪い？

クラス分類の機械学習の結果、全体の精度のほか、注目しているクラスの分類性能などについて確認しておく必要がある。

たとえば製造部品の良／不良を判別するケースで不良品の確率が1/1000などとても小さい場合や、疾病の判定をするケースで罹患する率が1万人に1人と非常に低い場合を考えてみる。

求めているのは、僅かに発生する不良品を選りだすことや、稀に罹患している人を特定することだ。このとき、ターゲットでない（正常品や罹患していない）多数のクラスを正確に分類できれば全体の正解率は上がる。ところがその一方で、求めている事象（不良品や罹患者の）ターゲットの分類精度が低いと、正解率には影響しないが本来求めているターゲットの分類機能としては低くなる。

誤判定の度合い

このほか、疾病に罹患していると判定したのに実際には罹患していない場合や、罹患していないと判定したのに実は罹患している場合など、分類器の誤判定の度合いも重要だ。

間違ってターゲットを特定してもいいから漏れがないようにしたいのか、誤って特定するリスクを避けたいのか、それらをもちいるケースに応じて分類器の性能がどうあるべきかを検討する必要がある。

confusion matrixの活用

このような場合にもちいられるのがconfusion matrixである。それは機械学習において用いられるテーブルで、クラス分類のターゲットクラスと予測されたクラスを行と列にとり、各々がどのように一致しているか／異なっているかを示したものである。

その要素と行／列の合計から、予測モデルの性能を示す様々な指標を計算することができる。

Confusion matrixの構成

3クラス分類の例

Confusion matrixは以下のように構成される。

行（または列）に正解のクラス列を、（または行）に予測されたクラスの列を、同じ順番で並べる
各正解クラスに対して、予測されたクラスの数を入れていく

例えば画像認識で果物を判別する予測モデルを考え、りんご、梨、洋梨の3つのクラスを分類するものとする。このとき、ある予測を行った結果として得られたconfusion matrixの一例を示す。

行の側に正解（事実）、列の側に予測（判定）を置いているが、この定義は場合によって入れ替わる（Scikit-learnのライブラリーではこれと同じ配置だが、WikipediaのConfusion matrixの解説では逆になっている）。

		予測
		りんご	梨	洋梨	計
正解	りんご	80	15	5	100
	梨	27	70	3	100
	洋梨	2	3	95	100
	計	109	88	103	300

2クラス分類での一般化

左記の果物のconfusion matrixは3クラス分類の例だが、以下は2クラス分類で考えていく。様々な2クラス分類におけるconfusion matrixの共通した構成を一般化したのが以下の図である。2つのクラスのうち”Positive”と表現しているのが「特に捕捉したい／注目している事象」に分類されるもので、たとえば不良品や疾病の発見、成長企業の特定など。他方はそれ以外で、製造品が正常、被検者が罹患していないといった捕捉の対象としない事象に対応する。

Positiveな事象（注目している事象）について、英語では”relevant instances”などの表現が使われているが、relevantの意味には、直訳の「関連～」だけではなく「重要な」というニュアンスもあるようなので、ここでは「注目」という用語を使う。

まず表の左上と右下について。予測した事象と実際の事象が一致している場合で色が同じになっている。この場合は予測が正しいという意味でTrueとする。左上はPositive（注目事象）と予測してそれが正しいのでTrue Positive (TP)と呼ぶ。また右下はNegative（注目していない事象）と予測してそれが正しいのでTrue Negative (TN)と呼ぶ。

次に左下と右上について。今度は予測結果の色と事実の色が異なっている。この場合は予測が誤っているという意味でFalseとする。左下はPositiveと予測したが誤りなのでFalse Positive (FP)、一方右上は予測がNegativeだがそれが誤りなのでFalse Negative (FN)と呼ぶ。

True/FalseとPositive/Negativeの順番とテーブル上の位置がややこしいが、常に予測・判定結果から見てそれが事実に対して正しいか誤りかと考えて「正しい／誤った、Positive判定／Negative判定」と定義されている。

Positive/Negativeを疾病検査の陽性(potitive)／陰性(negative)にあてはめると、TP：真陽性、TN：真陰性、FP：偽陽性、FN：偽陰性とも呼ばれる。

2クラス分類の例

疾患検査

次に2クラス分類の実世界での例を見てみる。まず、よくある例として、ある疾患にかかっているかどうかを検査する例。この場合はまさしく予測が陽性(Positive)か陰性(Negative)かに相当する。FPならば罹患していない人が不要な治療・対応を受けることになり、FNならば罹患している人を見逃すことになる。

犯人特定

次に、カメラ画像や様々な証拠などから犯人を見つけ出すような問題。対象者が犯人であるという事象に注目して、これをpositiveな判定としている。FPの場合は無実の人の誤認逮捕に結びつき、FNならば犯人を取り逃がすことになる。

ヒット商品予測

これまでの2つは、どちらかと言えば注目事象がよくない影響を及ぼすものだったが、これがよい効果をもたらす例を考えてみる。以下の例は、開発しようとしている商品がヒットするかどうか、いろいろな情報に基づいて予測しようとするものである。FPならヒットしない商品に無駄な投資をすることになり、FNならばヒット商品の開発の機会を逃すことになる。

4つの象限の意味・結果

これらの例も見ながら、confusion matrixの4つの象限がどのような意味を持つか、以下のように整理してみる。

TP: 注目対象を正しく分類する。対処すべき事象が特定できる。
FP: 注目すべきでないものを誤って注目対象に分類してしまう。注目対象が好ましくない事象の場合はその対策に余計なコストがかかったり、場合によっては謂れのない差別などの対象となったりする。好ましい事象の場合は、無駄なコストをかけることになる。統計学で言う第2種の過誤にあたる。
TN: 注目対象以外のものを正しく分類する。注目対象を誤って見逃すことがなく、被害の拡大や機会損失を避けられる。
FN: 注目対象を誤って注目対象以外に分類してしまう。捕捉すべき望ましくないものを見逃して影響が拡大したり、望ましいものを見逃して利得を得る機会を逃したりする。統計学で言う第1種の過誤にあたる。

指標

Confusion matrixの4つの象限の値から、複数の指標が導かれる。それぞれの和名には、時々異なるものを指している例があるので、英語表現を基本にする。

全体に対する率

まず、4象限全体（すなわち全データ数）に対する率を考える。これらは注目事象か非注目事象かに関わらない、モデル全体の正確さを表す。

Accuracy（正解率・正確度）

予測結果が正しく注目対象と非注目対象を言い当てた率。4象限の対角要素の合計の、総計に対する率を計算する。様々な機械学習モデルのスコアとして計算される値に相当する。Accuracyは「（ばらつきはともかく）予測が真値をどれだけ（平均的に）言い当てているか」という意味。このAccuracyを「精度」と呼んでいる場合があるが、科学的な表現としては少しずれている。

(1) $\begin{equation*} \frac{TP + TN}{TP + FP + FN + TN} \end{equation*}$

Error Rate（不正解率）

Accuracyと逆で、予測結果が実際の注目・非注目対象から外れた率。4象限の非対角要素の合計の、総計に対する率として計算する。

(2) $\begin{equation*} \frac{FP + FN}{TP + FP + FN + TN} \end{equation*}$

正解・事実に対する率

表の横方向の、各行それぞれの合計に対する率。正解に対するモデルの正確さを表す。

Sensitivity/Recall/TPR
（感度・再現率・検出率・真陽性率）

正解が注目事象の場合に、モデルも注目事象と分類する率。疾病検査を例にすると、その検査が疾病をとらえる「感度・検出率」となる。TNR(True Positive Rate)：真陽性率については、この後も真～率、偽～率が出てくるが、これらはすべて行方向に対する（正解・事実に関する）予測・分類の正確さと定義される。”recall”のニュアンス（呼び戻す、想起するなど）は、この指標の意味に繋がりにくい。むしろ無理やり日本語にしたような「再現率」の方がまだ本来の意味に近いと感じられる。

(3) $\begin{equation*} TPR = \frac{TP}{TP + FN} \end{equation*}$

Specificity/TNR（特異度・真陰性率）

正解が注目していない事象の場合に、モデルがそれを間違いなく分類した率。疾病検査なら、罹患していない人の結果が陰性となる率。「特異度」という訳は”specific”～「特別な」というあたりから名付けたのかもしれないがセンスが悪い。むしろこれは問題ないものを問題ないと分類する率だから、「特異」ではないはずだ。せめてspecify～特定するで「特定率」くらいならまだしもか。TNRはTrue Negative Rate。

(4) $\begin{equation*} TNR = \frac{TN}{FP + TN} \end{equation*}$

FNR（偽陰性率）

FNR(False Negative Rate)はTPRの裏で、正解が注目事象なのにそうでないと分類してしまった率。罹患しているのに検査で陰性となる率に相当する。

(5) $\begin{equation*} FNR = 1 - TPR = \frac{FN}{TP + FN} \end{equation*}$

FPR（偽陽性率）

FPR(False Positive Rate)はTNRの裏で、正解が注目していない事象なのに注目事象だと判定してしまった率。罹患していないのに検査で陽性だと判定されてしまう率に相当する。

(6) $\begin{equation*} FPR = 1 - TNR = \frac{FP}{FP + TN} \end{equation*}$

予測・判定結果に対する率

表の縦方向、各列の合計に対する率。分類結果がどの程度信頼できるかを表す。日本のサイトではPPV、高々NPVまでしか紹介されていないが、英語版のWikipediaではすべて図入りで説明されている。

Precision/PPV（適合度・精度・陽性的中率）

Precision(適合度）はモデルが注目事象と予測した場合に、実際にそれが注目事象である率。疾病検査で陽性判定の場合に実際に罹患している率に相当する。なお、科学上の表現でのprecision（精度）は、本来ばらつきの小ささを意味する。PPVはPositive Predictive Value。

(7) $\begin{equation*} PPV = \frac{TP}{TP + FP} \end{equation*}$

NPV（陰性的中率）

NPV(Negative Predictive Value)はモデルが注目事象ではないと分類して、それがあたっている率。疾病検査で陰性の場合に罹患していない率に相当する。敢えて日本語で言うなら「適正排除率」くらいか。

(8) $\begin{equation*} NPV = \frac{TN}{FN + TN} \end{equation*}$

FDR（陽性誤り率？）

FDR(False Discovery Rate)はモデルが注目事象であると分類したのに、実際には非注目事象である率。英語表現の直訳なら「間違って発見する率」。検査で陽性判定だが罹患していない率に相当する。日本語の訳はないが、敢えて言うなら「過剰陽性判定率」とか「陽性失中率」くらいか。

(9) $\begin{equation*} FDR = 1 - PPV = \frac{FP}{TP + FP} \end{equation*}$

FOR（陰性誤り率？）

FOR(False Omission Rate)はモデルが注目事象でないと分類したのに注目事象である率。疾病検査で陰性判定だが、実は罹患している率に相当する。英語表現の直訳なら「間違って無視してしまう率」だが、日本語なら「逸失率」くらいか。

(10) $\begin{equation*} FOR = 1 - NPV = \frac{FN}{FN + TN} \end{equation*}$

指標間のトレードオフに対する疑問

一般に、Sensitivity（感度・検出率）とPrecision（適合度・陽性的中率）はトレードオフの関係にある、と述べられることが多い。これは単純な仕組みで感度を上げようとするときに、注目対象以外でも多めに陽性と判定すれば率は上がるが、その場合は陽性判定でも注目対象以外のものが多くなって適合率が下がる、ということから来ている。

ここが少しおかしい。感度を上げるにはTPを大きくしFNを小さくしなければならない。このとき適合度の側から見れば、TPが大きくなるなら適合度も上がるし、FNを小さくしたときにFPが大きくなるという相関関係がなければならない。

実際には、見落としを少なくしようとすれば、無関係なケースを陽性と判定する「濡れ衣(FP)」は増えるだろう。しかしこの「濡れ衣(FP)」は、いくら増えても感度には寄与しない。これは感度が上がっていないのに（不安なので）陽性が多めに出るようにしているに過ぎないと思われる。だとすると、このような方針は単に適合度を下げているだけで感度は向上せず、「トレードオフ」とは言えない。

その他の指標

F値

感度と適合度のトレードオフにには疑問があるが、そのバランスを保って双方向上させるというのは重要だ。このような指標がF値(F value)と呼ばれるもので、感度と適合度の調和平均として定義されている。

(11) $\begin{equation*} F = \left( \frac{1}{2} \left( \frac{1}{TPR} + \frac{1}{PPV} \right) \right) ^{-1} = \frac{2 TPR \cdot PPV}{TPR + PPV} \end{equation*}$

これを4象限のパラメーターを使って書き直してみる。

(12) $\begin{align*} F = \left( \frac{1}{2} \left( \frac{TP + FN}{TP} + \frac{TP + FP}{TP} \right) \right) ^{-1} &= \frac{2 TP^2}{2TP + FN + FP} \\ &= \frac{TP}{1 + \dfrac{FN + FP}{TP}} \end{align*}$

一般にF値は感度と適合率のトレードオフを想定して、双方を加味した指標とされているが、双方のバランスがとれた状態がF値を最大化するというわけでもなさそうだ。

用語について

“confusion”は、LONGMAN、Cambridgeなどの英英辞典を見ると、（不明瞭な状況、人や物事などの誤認による）混乱・混迷、（不快な状況下での）困惑というニュアンスで、confusion matrixを的確に表現できるものがない。英語サイトで「confusion matrixの語源は何か？」という問いかけがいくつか見られた。どうも心理学にその元があるようだが、その中で言及されている”classification matrix”の方が明快に思われる。実際、TP、FNなどのタームや指標の名称がかなりconfusingなことからみると、アメリカ流のジョークとも思えてしまう。

和名は「混同行列」とされているが、これも何と何を混同するのか不明瞭だ。”confusion”の的確な訳ではないので、何となくそれに近い言葉を一生懸命にあてたのかもしれない。それならいっそのこと、より的確な用語（判定行列とか）をあてればよかったのにと感じる。

Python – 行・列単位の合計・率の計算

2020-11-11 / tau / コメントする

概要

2次元のndarrayやDataFrameで、行単位や列単位での合計を計算したり、それを使って行単位／列単位の率を計算する方法。

2次元のndarrayの場合は、

合計はsum()メソッドの引数にaxisを指定
- 列和ならaxis=0、行和ならaxis=1
- 結果は1次元配列で得られる
率の計算はこれらの合計の配列を使うが、列和に対する各列要素の率なら1次元配列の行ベクトルのまま、行和に対する各行要素の率なら2次元の列ベクトルに変換して除算

DataFrameの場合、まず合計を求めるには、

合計はsum()メソッドの引数にaxisを指定
- 列和ならaxis=0、行和ならaxis=1
- 結果はSeriesオブジェクトで得られる

その上で率の計算には2通りある。1つ目はSeriesオブジェクトの内容をndarrayとして取り出して計算する方法で、

Series.valuesで列和／行和の配列を取り出し、ndarrayの場合と同じ方法で計算する

もう1つの方法はSeriesオブジェクトのままでdiv()メソッドにaxisを指定する方法で、

列和に対する各列要素の率を計算するには、div(列和Series, axis=1)
行和に対する各行要素の率を計算するには、div(列和Series, axis=0)

`ndarray`の場合

確認

まず確認のために、以下の配列を準備する。

import numpy as np

a = np.arange(1, 10).reshape(3, 3)
vrow = np.arange(1, 4)
vcol = np.arange(1, 4).reshape(-1, 1)

print(a)
# [[1 2 3]
#  [4 5 6]
#  [7 8 9]]

print(vrow)
# [1 2 3]

print(vcol)
# [[1]
#  [2]
#  [3]]
# ]

import numpy as np

a = np.arange(1, 10).reshape(3, 3)

vrow = np.arange(1, 4)

vcol = np.arange(1, 4).reshape(-1, 1)

print(a)

# [[1 2 3]

# [4 5 6]

# [7 8 9]]

print(vrow)

# [1 2 3]

print(vcol)

# [[1]

# [2]

# [3]]

# ]

行ベクトルを2次元配列に加えると、配列の各行に対して行ベクトルが加えらえる。

print(a + vrow)

# [[ 2  4  6]
#  [ 5  7  9]
#  [ 8 10 12]]

print(a + vrow)

# [[ 2 4 6]

# [ 5 7 9]

# [ 8 10 12]]

列ベクトルを2次元配列に加えると、配列の各列に対して列ベクトルが加えられる。

print(a + vcol)

# [[ 2  3  4]
#  [ 6  7  8]
#  [10 11 12]]

print(a + vcol)

# [[ 2 3 4]

# [ 6 7 8]

# [10 11 12]]

つまり、ndarrayの2次元配列に行または列のベクトルを加えると、加える方のベクトルの形状に合わせて各行／列に演算が実行される。これは他のオペレーターについても同じ。

この演算は直感的にも分かりやすく、列ごと／行ごとの小計に対する比率の計算も思い浮かべることができる。

合計

ndarrayの各列／行に沿った合計を計算する。合計計算はndarrayのsum()メソッドを使うが、引数を省略すると全要素の和となる。引数にaxis=0を指定すると列方向に処理がされ（つまり配列の各列の要素が列方向に足され）、axis=1を指定すると行方向に処理がされる（つまり配列の各行の要素が行方向に足される）。

以下の例は、最初に使った2次元配列の列方向の和（の行ベクトル）と行方向の和（の列ベクトル）を計算している。

sums_along_cols = a.sum(axis=0)
sums_along_rows = a.sum(axis=1).reshape(-1, 1)

print(sums_along_cols)
# [12 15 18]

print(sums_along_rows)
# [[ 6]
#  [15]
#  [24]]

sums_along_cols = a.sum(axis=0)

sums_along_rows = a.sum(axis=1).reshape(-1, 1)

print(sums_along_cols)

# [12 15 18]

print(sums_along_rows)

# [[ 6]

# [15]

# [24]]

率の計算

2次元配列を行合計ベクトルで割ると、各行の要素が合計ベクトルの各要素で割られる（各列の行要素の合計は1になる）。

print(a / sums_along_cols)

# [[0.08333333 0.13333333 0.16666667]
#  [0.33333333 0.33333333 0.33333333]
#  [0.58333333 0.53333333 0.5       ]]

# NOTE
# 1/12 2/15 3/18
# 4/12 5/15 6/18
# 7/12 8/15 9/18

print(a / sums_along_cols)

# [[0.08333333 0.13333333 0.16666667]

# [0.33333333 0.33333333 0.33333333]

# [0.58333333 0.53333333 0.5 ]]

# NOTE

# 1/12 2/15 3/18

# 4/12 5/15 6/18

# 7/12 8/15 9/18

また、2次元配列を列合計ベクトルで割ると、各列の要素が合計ベクトルの各要素で割られる。繰り返しになるが、この場合の合計ベクトルは2次元の列ベクトルになっている。

print(a / sums_along_rows)

# [[0.16666667 0.33333333 0.5       ]
#  [0.26666667 0.33333333 0.4       ]
#  [0.29166667 0.33333333 0.375     ]]

# NOTE
#  1/6  2/6  3/6
# 4/15 5/15 6/15
# 7/24 8/24 9/24

print(a / sums_along_rows)

# [[0.16666667 0.33333333 0.5 ]

# [0.26666667 0.33333333 0.4 ]

# [0.29166667 0.33333333 0.375 ]]

# NOTE

# 1/6 2/6 3/6

# 4/15 5/15 6/15

# 7/24 8/24 9/24

以下の図のように、この

`DataFrame`の場合

確認

まず確認のために、以下のDataFrameとSeriesを準備する。2次元配列はndarrayの例と同じものを流用。

import pandas as pd

df = pd.DataFrame(a)
sr = pd.Series(np.arange(1, 4))

print(df)
#    0  1  2
# 0  1  2  3
# 1  4  5  6
# 2  7  8  9

import pandas as pd

df = pd.DataFrame(a)

sr = pd.Series(np.arange(1, 4))

print(df)

# 0 1 2

# 0 1 2 3

# 1 4 5 6

# 2 7 8 9

合計ベクトルが`ndarray`の場合

演算をほどこすベクトルがndarrayの場合、2次元配列の時と同じように、ベクトルが行／列によって自動的に加えられる方向が決められる。

vrow = np.arange(1, 4)
vcol = np.arange(1, 4).reshape(-1, 1)

print(vrow)
# [1 2 3]

print(vcol)
# [[1]
#  [2]
#  [3]]

print(df + vrow)
#    0   1   2
# 0  2   4   6
# 1  5   7   9
# 2  8  10  12

print(df + vcol)
#     0   1   2
# 0   2   3   4
# 1   6   7   8
# 2  10  11  12

vrow = np.arange(1, 4)

vcol = np.arange(1, 4).reshape(-1, 1)

print(vrow)

# [1 2 3]

print(vcol)

# [[1]

# [2]

# [3]]

print(df + vrow)

# 0 1 2

# 0 2 4 6

# 1 5 7 9

# 2 8 10 12

print(df + vcol)

# 0 1 2

# 0 2 3 4

# 1 6 7 8

# 2 10 11 12

合計ベクトルが`Series`の場合

DataFrameのsum()メソッドで行や列の合計を計算するとSeriesオブジェクトで結果が得られるため、その挙動を確認しておく。

sr = pd.Series(np.arange(1, 4))
print(sr)

0    1
1    2
2    3
dtype: int32

sr = pd.Series(np.arange(1, 4))

print(sr)

0 1

1 2

2 3

dtype: int32

Seriesを単純にDataFrameと演算子で結ぶと、行ベクトルとして扱われる。

print(df + sr)

#    0   1   2
# 0  2   4   6
# 1  5   7   9
# 2  8  10  12

print(df + sr)

# 0 1 2

# 0 2 4 6

# 1 5 7 9

# 2 8 10 12

Seriesは行・列の概念を持たないが、演算の方向を明示するのに以下の方法をとる。

演算子の代わりに演算メソッドを使う
演算メソッドの引数axisで演算の方向を指定する

print(df.add(sr, axis=0))

#     0   1   2
# 0   2   3   4
# 1   6   7   8
# 2  10  11  12

print(df.add(sr, axis=1))

#    0   1   2
# 0  2   4   6
# 1  5   7   9
# 2  8  10  12

print(df.add(sr, axis=0))

# 0 1 2

# 0 2 3 4

# 1 6 7 8

# 2 10 11 12

print(df.add(sr, axis=1))

# 0 1 2

# 0 2 4 6

# 1 5 7 9

# 2 8 10 12

演算子に対応するメソッドは、add、sub、mul、div、mod、powが準備されている。

合計

列／行ごとの合計は、ndarrayと同じくDataFrameのsum()メソッドで引数axisを指定して計算する。

sums_along_cols = df.sum(axis=0)
sums_along_rows = df.sum(axis=1)

print(sums_along_cols)
# 0    12
# 1    15
# 2    18
# dtype: int64

print(sums_along_rows)
# 0     6
# 1    15
# 2    24
# dtype: int64

sums_along_cols = df.sum(axis=0)

sums_along_rows = df.sum(axis=1)

print(sums_along_cols)

# 0 12

# 1 15

# 2 18

# dtype: int64

print(sums_along_rows)

# 0 6

# 1 15

# 2 24

# dtype: int64

率の計算

`Series`を`ndarray`で取り出して計算する方法

Seriesのvaluesプロパティーでその内容をndarrayとして取り出せる。

列和に対しては、それを行ベクトルのまま除算すれば、各列の要素を行ベクトルの対応する要素で割った値となる。

print(df / sums_along_cols.values)

#           0         1         2
# 0  0.083333  0.133333  0.166667
# 1  0.333333  0.333333  0.333333
# 2  0.583333  0.533333  0.500000

print(df / sums_along_cols.values)

# 0 1 2

# 0 0.083333 0.133333 0.166667

# 1 0.333333 0.333333 0.333333

# 2 0.583333 0.533333 0.500000

行和に対しては、ndarrayを1列の列ベクトルに変換して除算すれば、各行の要素を列ベクトルの対応する要素で割った値となる。

print(df / sums_along_rows.values.reshape(-1, 1))

#           0         1      2
# 0  0.166667  0.333333  0.500
# 1  0.266667  0.333333  0.400
# 2  0.291667  0.333333  0.375

print(df / sums_along_rows.values.reshape(-1, 1))

# 0 1 2

# 0 0.166667 0.333333 0.500

# 1 0.266667 0.333333 0.400

# 2 0.291667 0.333333 0.375

`Series`のままで計算する方法

合計Seriesをそのまま使って除算する場合はDataFrameのdiv()メソッドを使いaxisを指定するが、axisの指定の仕方に注意が必要。

列和で除算する場合は、その各要素が各行の各要素に対応するため、列和を行ベクトルとみて各行に除算を適用する(axis=1)。

print(df.div(sums_along_cols, axis=1))

#           0         1         2
# 0  0.083333  0.133333  0.166667
# 1  0.333333  0.333333  0.333333
# 2  0.583333  0.533333  0.500000

print(df.div(sums_along_cols, axis=1))

# 0 1 2

# 0 0.083333 0.133333 0.166667

# 1 0.333333 0.333333 0.333333

# 2 0.583333 0.533333 0.500000

分かりにくいので図示すると以下のようになる。まず合計を求めるのにsum()メソッドでaxis=0として列和を求める。この合計で各要素を割るのに1列目の要素は合計の1つ目の要素、2列目の要素は合計の2つ目の要素・・・で割る必要があるので、div()メソッドでaxis=1とする。こうすると合計のSeriesは行ベクトルとみなされて、それが各行の要素に適用される。

行和で除算する場合は、その各要素が各列の各要素に対応するため、行和を列ベクトルとみて各列に除算を適用する(axis=0)。

print(df.div(sums_along_rows, axis=0))

#           0         1      2
# 0  0.166667  0.333333  0.500
# 1  0.266667  0.333333  0.400
# 2  0.291667  0.333333  0.375

print(df.div(sums_along_rows, axis=0))

# 0 1 2

# 0 0.166667 0.333333 0.500

# 1 0.266667 0.333333 0.400

# 2 0.291667 0.333333 0.375

これも分かりにくいので以下のように図示する。合計を求めるのにsum()メソッドでaxis=1として行和を求める。この合計で各要素を割るのに1行目の要素は合計の1つ目の要素、2行目の要素は合計の2つ目の要素・・・で割る必要があるので、div()メソッドでaxis=1とする。こうすると合計のSeriesは列ベクトルとみなされて、それが各列の要素に適用される。

このように、div()のような演算子メソッドでaxisを使う方法はややこしい（少なくとも私には）。

実行速度

各計算方法の実行速度には、あまり大きな差は出なかった。

import numpy as np
import pandas as pd
import time

counter = 100

a = np.arange(1000000).reshape(1000, 1000)
df = pd.DataFrame(a)

sum_of_cols = df.sum(axis=0)
sum_of_rows = df.sum(axis=1)

t0 = time.time()

for n in range(counter):
    df.values / sum_of_cols.values
t_array_c = time.time()
print(t_array_c - t0)

for n in range(counter):
    df.values / sum_of_rows.values.reshape(-1, 1)
t_array_r = time.time()
print(t_array_r - t_array_c)

for n in range(counter):
    df / sum_of_cols.values
t_df_array_c = time.time()
print(t_df_array_c - t_array_r)

for n in range(counter):
    df / sum_of_rows.values.reshape(-1, 1)
t_df_array_r = time.time()
print(t_df_array_r - t_df_array_c)

for n in range(counter):
    df.div(sum_of_cols, axis=1)
t_df_div_c = time.time()
print(t_df_div_c - t_df_array_r)

for n in range(counter):
    df.div(sum_of_rows, axis=0)
t_df_div_r = time.time()
print(t_df_div_r - t_df_div_c)

import numpy as np

import pandas as pd

import time

counter = 100

a = np.arange(1000000).reshape(1000, 1000)

df = pd.DataFrame(a)

sum_of_cols = df.sum(axis=0)

sum_of_rows = df.sum(axis=1)

t0 = time.time()

for n in range(counter):

df.values / sum_of_cols.values

t_array_c = time.time()

print(t_array_c - t0)

for n in range(counter):

df.values / sum_of_rows.values.reshape(-1, 1)

t_array_r = time.time()

print(t_array_r - t_array_c)

for n in range(counter):

df / sum_of_cols.values

t_df_array_c = time.time()

print(t_df_array_c - t_array_r)

for n in range(counter):

df / sum_of_rows.values.reshape(-1, 1)

t_df_array_r = time.time()

print(t_df_array_r - t_df_array_c)

for n in range(counter):

df.div(sum_of_cols, axis=1)

t_df_div_c = time.time()

print(t_df_div_c - t_df_array_r)

for n in range(counter):

df.div(sum_of_rows, axis=0)

t_df_div_r = time.time()

print(t_df_div_r - t_df_div_c)

実行時間は以下の通りで、各計算手法の間に差はない。敢えて言えば、DataFrameを使った場合に僅かに時間がかかっている。

0.4757249355316162
0.4577751159667969
0.524827241897583
0.5440847873687744
0.5674364566802979
0.5242717266082764

0.4757249355316162

0.4577751159667969

0.524827241897583

0.5440847873687744

0.5674364566802979

0.5242717266082764

なお、この計算はpandasのバージョン1.1.4で実行したが、upgrade前のバージョン0では、3つ目と5つ目、DataFrameで行単位の演算を行うときに20秒台と2桁長い時間がかかっていた。

axisの方向

2020-11-08 / tau / コメントする

概要

配列などのメソッドの引数で指定するaxis=0/1について確認。

axis=0は配列やDataFrameを列単位で捉えて、その列の中で処理を行いながら、すべての列に対して処理が行われる
- SerieseオブジェクトがDataFrameの処理の対象となる場合は列として扱われ、DataFrameの各列を処理しながらすべての列に適用される
axis=1は配列やDataFrameを行単位で捉えて、その行の中で処理を行いながら、すべての行に対して処理が行われる。
- SerieseオブジェクトがDataFrameの処理の対象となる場合は行として扱われ、DataFrameの各行を処理しながらすべての行に適用される

`ndarray`の場合

まずndarrayの2次元配列で確認する。

ary = np.arange(1, 7).reshape(2, 3)
print(ary)

# [[1 2 3]
#  [4 5 6]]

ary = np.arange(1, 7).reshape(2, 3)

print(ary)

# [[1 2 3]

# [4 5 6]]

`max()`メソッド

axis=0は列単位で各列の最大値を探し、それらを要素とする配列（要素数＝列数の1次元配列）
axis=1は行単位で各行の最大値を探し、それらを要素とする配列（要素数＝行数の1次元配列）

print(ary.max(axis=0))
# [4 5 6]

print(ary.max(axis=1))
# [3 6]

print(ary.max(axis=0))

# [4 5 6]

print(ary.max(axis=1))

# [3 6]

`sum()`メソッド

axis=0は列単位で各列の合計を要素とする配列（要素数＝列数の1次元配列）
axis=1は行単位で各行の合計を要素とする配列（要素数＝行数の1次元配列）

print(ary.sum(axis=0))
# [5 7 9]

print(ary.sum(axis=1))
# [ 6 15]

print(ary.sum(axis=0))

# [5 7 9]

print(ary.sum(axis=1))

# [ 6 15]

`repeat()`メソッド

axis=0は列単位で各列の要素が指定回数繰り返される
axis=1は行単位で各行の要素が指定回数繰り返される

print(ary.repeat(2, axis=0))

# [[1 2 3]
#  [1 2 3]
#  [4 5 6]
#  [4 5 6]]

print(ary.repeat(2, axis=1))

# [[1 1 2 2 3 3]
#  [4 4 5 5 6 6]]

print(ary.repeat(2, axis=0))

# [[1 2 3]

# [1 2 3]

# [4 5 6]

# [4 5 6]]

print(ary.repeat(2, axis=1))

# [[1 1 2 2 3 3]

# [4 4 5 5 6 6]]

図による理解

sum()メソッドを例に、axis=0/1に対する挙動を図にすると、以下のようになる。

`DataFrame`の場合

以下のDataFrameとSeriesオブジェクトで確認する。Seriesオブジェクトは行として扱われ、array_like、1次元の配列でも同じ結果になる。

df = pd.DataFrame(np.arnge(1, 9).reshape(3, 3))
sr = pd.Series(np.arange(1, 3))

print(df)

#    0  1  2
# 0  1  2  3
# 1  4  5  6
# 2  7  8  9

print(sr)

# 0    1
# 1    2
# 2    3
# dtype: int32

df = pd.DataFrame(np.arnge(1, 9).reshape(3, 3))

sr = pd.Series(np.arange(1, 3))

print(df)

# 0 1 2

# 0 1 2 3

# 1 4 5 6

# 2 7 8 9

print(sr)

# 0 1

# 1 2

# 2 3

# dtype: int32

`min()`メソッドなど

min、max、sumなどのメソッドの考え方はndarrayと同じ挙動。

`add()`などの演算メソッド

DataFrameには演算子による演算の代替となるメソッドがある（add、sub、mul、div、mod、pow）。addメソッドを例にとると、以下のように引数を指定。

add(array_like, axis=0/1)

axis=0はarray_likeを列とみなして、DataFrameオブジェクトの各列の要素との和を計算する
axis=1はarray_likeを行とみなして、DataFrameオブジェクトの各行の要素との和を計算する

print(df.add(sr, axis=0))

#     0   1   2
# 0   2   3   4
# 1   6   7   8
# 2  10  11  12

print(df.add(sr, axis=1))

#    0   1   2
# 0  2   4   6
# 1  5   7   9
# 2  8  10  12

print(df.add(sr, axis=0))

# 0 1 2

# 0 2 3 4

# 1 6 7 8

# 2 10 11 12

print(df.add(sr, axis=1))

# 0 1 2

# 0 2 4 6

# 1 5 7 9

# 2 8 10 12

`apply()`メソッド

applyメソッドは、行または列を指定した関数に渡す。

axis=0はDataFrameオブジェクトの各列を指定した関数に渡す
axis=1はDataFrameオブジェクトの各行を指定した関数に渡す

print(df.apply(lambda x: sum(x) / 3, axis=0))

# 0    4.0
# 1    5.0
# 2    6.0
# dtype: float64

print(df.apply(lambda x: sum(x) / 3, axis=1))

# 0    2.0
# 1    5.0
# 2    8.0
# dtype: float64

print(df.apply(lambda x: sum(x) / 3, axis=0))

# 0 4.0

# 1 5.0

# 2 6.0

# dtype: float64

print(df.apply(lambda x: sum(x) / 3, axis=1))

# 0 2.0

# 1 5.0

# 2 8.0

# dtype: float64

演算メソッドの図による理解

演算メソッドは少し挙動が違うので図で整理しておく。1次元のarray_likeオブジェクトがaxisの指定によって列／行としてみなされる点に注意。

補足

1次元配列の場合

1次元配列に対してaxis引数を使う場合、行ベクトルとしてaxis=1に反応しそうだが、実際にはaxis=0で各要素に対する処理が行われる。axis=1を指定すると、たとえば以下のようなエラーになる。

ary.min(axis=1)->numpy.AxisError: axis 1 is out of bounds for array of dimension 1

1	ary.min(axis=1)->numpy.AxisError: axis 1 is out of bounds for array of dimension 1

元々多次元配列を意図した引数なので、1次元配列に使うのはナンセンスだろう（axis=0を行単位の処理にしておけば自然ではあったかもしれないが）。

1行の2次元配列の場合

1行の配列（1つの1次元配列を要素に持つ2次元配列：行ベクトル）に対するaxisの効果を、sumメソッドで見てみる。

ary = np.arange(1, 4).reshape(1, -1)
print(ary)
# [[1 2 3]]

print(ary.sum(axis=0))
# [1 2 3]

print(ary.sum(axis=1))
# [6]

ary = np.arange(1, 4).reshape(1, -1)

print(ary)

# [[1 2 3]]

print(ary.sum(axis=0))

# [1 2 3]

print(ary.sum(axis=1))

# [6]

axis=0の場合は各要素が1要素の列ベクトルとみなされ、3つの列（要素）ごとに処理される。その結果は3つの要素を持つ1次元配列（行ベクトル）となる。

axis=1の場合は行ベクトル全体が1つの行とみなされ、それらの要素に対して処理がされる。その結果は1つの数値となるが、1つの要素を持つ1次元配列で返される。

1列の2次元配列の場合

1列の配列（列ベクトル）に対するaxisの効果を、sumメソッドで見てみる。

ary = np.arange(1, 4).reshape(-1, 1)
print(ary)
# [[1]
#  [2]
#  [3]]

print(ary.sum(axis=0))
# [6]

print(ary.sum(axis=1))
# [1 2 3]

ary = np.arange(1, 4).reshape(-1, 1)

print(ary)

# [[1]

# [2]

# [3]]

print(ary.sum(axis=0))

# [6]

print(ary.sum(axis=1))

# [1 2 3]

axis=0の場合は3つの要素を持つ1つの列に対して処理される。その結果は1つの数値となるが、1つの要素を持つ1次元配列で返される。

axis=1の場合は列の各要素が1要素の行とみなされ、3つの行（要素）ごとに処理される。その結果は3つの要素を持つ列ベクトルだが、3つの要素を持つ1次元配列（行ベクトル）で返される。

概要

手順

主成分と寄与率

可視化

3次元

2次元

まとめ

概要

計算の手順

主成分と寄与率

可視化

3次元

2次元

主成分分析の特徴

概要

最大化すべき分散の導出

分散の最大化

特徴量が2つの場合

第2主成分以降

主成分の意味

概要

基本

属性データが複数列の場合

属性が数値表現の場合

属性名の指定

概要

DataFrameの規模・形状

データの先頭部分と末尾部分

info()～各列の基本情報の表示

dscribe()～基本的な統計量

属性変数のカウント

1次元配列の場合

平均

分散・標準偏差

標本分散・標本標準偏差

不偏推定量

2次元配列の場合

平均

分散・標準偏差

標本分散・標本標準偏差

不偏推定量

共分散

概要

指標

PPV・NPVの計算式の導出

元データの構成

PPV・NPVの計算式

パラメーターに応じたPPV・NPVの変化

PPV

NPV

PPVとNPVの関係

シミュレーションによる挙動確認

概要

精度が高いのに性能が悪い？

誤判定の度合い

confusion matrixの活用

Confusion matrixの構成

3クラス分類の例

2クラス分類での一般化

2クラス分類の例

疾患検査

犯人特定

ヒット商品予測

4つの象限の意味・結果

指標

全体に対する率

Accuracy（正解率・正確度）

Error Rate（不正解率）

正解・事実に対する率

Sensitivity/Recall/TPR （感度・再現率・検出率・真陽性率）

Specificity/TNR（特異度・真陰性率）

FNR（偽陰性率）

FPR（偽陽性率）

予測・判定結果に対する率

Precision/PPV（適合度・精度・陽性的中率）

NPV（陰性的中率）

FDR（陽性誤り率？）

FOR（陰性誤り率？）

指標間のトレードオフに対する疑問

その他の指標

`info()`～各列の基本情報の表示

`dscribe()`～基本的な統計量

Sensitivity/Recall/TPR
（感度・再現率・検出率・真陽性率）

`ndarray`の場合

`DataFrame`の場合

合計ベクトルが`ndarray`の場合

合計ベクトルが`Series`の場合

`Series`を`ndarray`で取り出して計算する方法

`Series`のままで計算する方法

`ndarray`の場合

`max()`メソッド

`sum()`メソッド

`repeat()`メソッド

`DataFrame`の場合

`min()`メソッドなど

`add()`などの演算メソッド

`apply()`メソッド