numpy – bincount

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

概要

numpy.bincount()関数の仕様

整数型の配列を引数にとる
配列中、同じ値の要素の個数をカウントする
0～要素の最大値を要素とし、各要素番号に対応する値の個数を要素とする配列を返す
元のデータの要素ごとの重みを指定することができる

使い方

基本形

引数で与えた整数型配列中の同じ値をカウントして、各値ごとの個数を要素とする配列を返す。

import numpy as np

a = np.array([0, 1, 1, 2, 2, 2, 3, 3])
print(np.bincount(a))

# [1 2 3 2]

import numpy as np

a = np.array([0, 1, 1, 2, 2, 2, 3, 3])

print(np.bincount(a))

# [1 2 3 2]

上の結果の意味は、0が1個、1が2個、2が3個、3が2個。

値が飛んでいる場合

引数の配列中、0～最大値までの整数値に対する数をカウントする。値が存在しない場合の個数は0。

b = np.array([1, 3, 3, 5, 5, 5])
print(np.bincount(b))

# [0 1 0 2 0 3]

b = np.array([1, 3, 3, 5, 5, 5])

print(np.bincount(b))

# [0 1 0 2 0 3]

上の例では、0～5までの個数がカウントされ、0, 2, 4は配列中に存在しないので0となっている。

順不同

引数の配列中の要素は昇順である必要はない。

c = np.array([2, 2, 3, 3, 3, 1])
print(np.bincount(c))

# [0 1 2 3]

c = np.array([2, 2, 3, 3, 3, 1])

print(np.bincount(c))

# [0 1 2 3]

`weights`の意味

引数にweightsを指定する場合。

データの配列と同じ要素数のweightsの配列を与える。
要素をカウントの場合に1ずつ足すのではなく、各要素の位置に対応した重みが加算されていく

w = np.array([0.1, 0.2, 0.3, 0.4, 0.5, 0.6])
d = np.array([1, 2, 4, 2, 4, 4])
print(np.bincount(d, weights=w))

# [0.  0.1 0.6 0.  1.4]

w = np.array([0.1, 0.2, 0.3, 0.4, 0.5, 0.6])

d = np.array([1, 2, 4, 2, 4, 4])

print(np.bincount(d, weights=w))

# [0. 0.1 0.6 0. 1.4]

上の例では以下のように動作している。

0は存在しないので0
1は0番目に1つだけ存在し、その位置のweightsの値は0.1
2は1番目と3番目に存在するので、weightsの第1要素0.2と第3要素0.4を加えて0.6
3は存在しないので0
4は2番目、4番目、5番目にあるので、weightsの第2要素0.3、第4要素0.5、第5要素0.6を加えて1.4

pyplot.imshow – 画像表示

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

概要

matplotlib.pyplot.imshow()は画像表示用のメソッドで、表示対象として、画像ファイルや画像情報を格納した配列を指定する。

pyplotやsubplotで直接実行するほか、Axesオブジェクトのメソッドとしても実行できる。

ピクセルデータのレンジのデフォルト設定と与えるデータのレンジによって予期しない結果になることもあり、vmin、vmaxを明示的に指定した方がよい。

画像ファイルの表示

以下のコードは、JPEGファイルを読み込んで表示する。

ここではpyplot.subplotのメソッドとしてimshow()を実行している。画像が1つの場合、pyplot.imshow()でもよい。

1つは画像ファイルをそのまま引数にし、もう1つは画像ファイルを配列の形にしてから引数に渡している。画像の配列の形については後述。

import numpy as np
import matplotlib.pyplot as plt
from PIL import Image

image = Image.open("./coffeemill.jpg")

image_array = np.asarray(image)

plt.subplot(121).imshow(image)
plt.subplot(122).imshow(image_array)
plt.show()

print(image_array.shape)
# (871, 653, 3)

print(image_array)
[[[171 147 111]
  [167 143 107]
  [170 147 113]
  ...
  [ 93  82  86]
  [ 92  81  87]
  [ 93  83  91]]

 [[170 146 110]
  [167 143 107]
  [170 147 113]
  ...
  [ 92  81  87]
  [ 85  76  81]
  [ 77  67  75]]

  .....

 [[ 36  23  15]
  [ 37  24  16]
  [ 41  26  21]
  ...
  [248 196 123]
  [249 197 124]
  [247 195 122]]]

import numpy as np

import matplotlib.pyplot as plt

from PIL import Image

image = Image.open("./coffeemill.jpg")

image_array = np.asarray(image)

plt.subplot(121).imshow(image)

plt.subplot(122).imshow(image_array)

plt.show()

print(image_array.shape)

# (871, 653, 3)

print(image_array)

[[[171 147 111]

[167 143 107]

[170 147 113]

...

[ 93 82 86]

[ 92 81 87]

[ 93 83 91]]

[[170 146 110]

[167 143 107]

[170 147 113]

...

[ 92 81 87]

[ 85 76 81]

[ 77 67 75]]

.....

[[ 36 23 15]

[ 37 24 16]

[ 41 26 21]

...

[248 196 123]

[249 197 124]

[247 195 122]]]

配列の画像表示

基本形

imshow()は配列を引数にとることができる。

以下の例では、カラーマップを指定して2×2=4要素の2次元配列を表示している。最小値0がカラーマップbwrの青に、最大値255が赤に対応し、その間の数値の大きさに応じたカラーマップ上の色が選択されている（デフォルトのcmapはvirいdis）。

なお、この例ではpyplotから直接imshow()を実行している。

import numpy as np
import matplotlib.pyplot as plt

im_array = np.array([
    [0, 85],
    [170, 255]
])

plt.imshow(im_array, cmap='bwr')

plt.show()

import numpy as np

import matplotlib.pyplot as plt

im_array = np.array([

[0, 85],

[170, 255]

])

plt.imshow(im_array, cmap='bwr')

plt.show()

レンジ

imshow()に配列を渡して描画させるとき、数値のレンジに留意する必要がある。

デフォルトでは、imshow()は渡された配列の中の最小値と最大値をカラーマップの下限値と上限値に対応させ、線形にマッピングする。

なお、この例ではarray-likeとして2次元のリストを渡していて、Axesからimshow()を呼び出している。

import numpy as np
import matplotlib.pyplot as plt

images = []
images.append([[0, 0.5, 1]])
images.append([[0, 127, 255]])
images.append([[-1, 0, 1]])
images.append([[1000, 1000.5, 1001]])
print(images)

fig, axs = plt.subplots(2, 2)

for image, ax in zip(images, axs.flatten()):
    ax.imshow(image, cmap='bwr')

plt.show()

import numpy as np

import matplotlib.pyplot as plt

images = []

images.append([[0, 0.5, 1]])

images.append([[0, 127, 255]])

images.append([[-1, 0, 1]])

images.append([[1000, 1000.5, 1001]])

print(images)

fig, axs = plt.subplots(2, 2)

for image, ax in zip(images, axs.flatten()):

ax.imshow(image, cmap='bwr')

plt.show()

4つの配列はそれぞれ最小値と最大値が異なり、かつその中央の値を持つ。値は異なるが全て最小値がカラーマップbwrの下限値に対応する青、最大値が上限値に対応する赤、中央値は白となっている（特段フランス国旗を意図したものではない）。

viminとvmax

imshow()の引数でvminとvmaxを設定すると、配列の値に関わらず、vminとvmiaxをカラーマップの下限値と上限値に対応させる。

以下の例では最小値0、最大値1の2要素の配列を、vminとvmaxを変えてカラーマップbwrで描画させている。

import numpy as np
import matplotlib.pyplot as plt

image = [[0, 1]]

fig, axs = plt.subplots(2, 2)

axs[0, 0].imshow(image, cmap='bwr')
axs[0, 0].set_title("Default")
axs[0, 1].imshow(image, cmap='bwr', vmin=0, vmax=2)
axs[0, 1].set_title("vmin=0, vmax=2")
axs[1, 0].imshow(image, cmap='bwr', vmin=-1, vmax=2)
axs[1, 0].set_title("vmin=-1, vmax=2")
axs[1, 1].imshow(image, cmap='bwr', vmin=0.2, vmax=0.8)
axs[1, 1].set_title("vmin=0.2, vmax=0.8")

plt.show()

import numpy as np

import matplotlib.pyplot as plt

image = [[0, 1]]

fig, axs = plt.subplots(2, 2)

axs[0, 0].imshow(image, cmap='bwr')

axs[0, 0].set_title("Default")

axs[0, 1].imshow(image, cmap='bwr', vmin=0, vmax=2)

axs[0, 1].set_title("vmin=0, vmax=2")

axs[1, 0].imshow(image, cmap='bwr', vmin=-1, vmax=2)

axs[1, 0].set_title("vmin=-1, vmax=2")

axs[1, 1].imshow(image, cmap='bwr', vmin=0.2, vmax=0.8)

axs[1, 1].set_title("vmin=0.2, vmax=0.8")

plt.show()

左上はデフォルトなので、最小値0がカラーマップ下限値に対応した青に、最大値1が上限値に対応した赤になっている。

右上はvmin=0で配列の最小値0と同じだが、vmax=2としている。このため配列の0はカラーマップ下限の青で、配列の1はカラーマップ中央の白になっている。

左下はvmin=-1も設定されているので、配列の0、1はカラーマップの左から1/3、2/3に相当する色となっている。

右下はvminとvmaxが配列の最小値と最大値の範囲より内側にある。このため、配列の最小値・最大値はそれぞれカラーマップの下限・上限に対応する青・赤となっている。

RGB

array-likeの次元が3次元になると、RGB/RGBA形式だと認識される。

[rows, cols, 3]: 3次元目のサイズが3の時はRGB表現と認識される。1次元目と2次元目はそれぞれ画像の行数と列数とみなされ、3次元目は3つの列がR, G, Bの値に対応する。
[rows, cols, 4]: 3次元目のサイズが4の時はRGBA表現と認識される。1次元目と2次元目はそれぞれ画像の行数と列数とみなされ、3次元目は3つの列がR, G, Bの値に対応し、4つ目の列が透明度に対応する。

R, G, B, Aの値は、配列のdtypeがint形式の時には0～255、floatの時には0～1の範囲が想定される。

以下の例の内容。

画像サイズを2行×4列として、R, G, Bごとに画像のピクセルデータを設定→shape=(3, 2, 4)
ピクセル並び替え後の配列を4つ準備
forループでピクセル並び替え
画像表示とデータ内容の表示

3次元配列のピクセルの並び替えは、泥臭くforループで回しているが、もっとエレガントな方法があるかもしれない（もとから(3, rows, cols)の形にしてくれればよかったのに）。

# --- --B -G- -GB
# R-- R-B RG- RGB

import numpy as np
import matplotlib.pyplot as plt

image_org = np.array([
    [
        [  0,   0,   0,   0],
        [255, 255, 255, 255]
    ],
    [
        [  0,   0, 255, 255],
        [  0,   0, 255, 255]
    ],
    [
        [  0, 255,   0, 255],
        [  0, 255,   0, 255]
    ]
])

row = image_org.shape[1]
col = image_org.shape[2]

image_int1 = np.empty((row, col, 3), dtype=int)
image_int2 = np.empty((row, col, 3), dtype=int)
image_float1 = np.empty((row, col, 3), dtype=float)
image_float2 = np.empty((row, col, 3), dtype=float)

for m in range(row):
    for n in range(col):
        for rgb in range(3):
            image_int1[m, n, rgb] = image_org[rgb, m, n]
            image_int2[m, n, rgb] = image_org[rgb, m, n] /255
            #-> all the pixels become BLACK
            image_float1[m, n, rgb] = image_org[rgb, m, n]
            image_float2[m, n, rgb] = image_org[rgb, m, n] / 255
            #- >works, but
            # "Clipping input data to the valid range for imshow with RGB data ([0..1] for floats or [0..255] for integers)."

fig, axs = plt.subplots(2, 2)

axs[0, 0].imshow(image_int1)
axs[0, 0].set_title("dtype=int, 0~255")
axs[0, 1].imshow(image_int2)
axs[0, 1].set_title("dtype=int, 0.0~1.0")
axs[1, 0].imshow(image_float1)
axs[1, 0].set_title("dtype=float, 0.0~255.0")
axs[1, 1].imshow(image_float2)
axs[1, 1].set_title("dtype=float, 0.0~1.0")

plt.show()

print("original image shape  :{}".format(image_org.shape))
print("rearranged image shape:{}".format(image_int1.shape))
print("rearranged image:\n{}".format(image_int1))

# Clipping input data to the valid range for imshow with RGB data ([0..1] for floats or [0..255] for integers).
# original image shape  :(3, 2, 4)
# rearranged image shape:(2, 4, 3)
# rearranged image:
# [[[  0   0   0]
#   [  0   0 255]
#   [  0 255   0]
#   [  0 255 255]]
# 
#  [[255   0   0]
#   [255   0 255]
#   [255 255   0]
#   [255 255 255]]]

# --- --B -G- -GB

# R-- R-B RG- RGB

import numpy as np

import matplotlib.pyplot as plt

image_org = np.array([

[

[ 0, 0, 0, 0],

[255, 255, 255, 255]

[

[ 0, 0, 255, 255],

[ 0, 0, 255, 255]

[

[ 0, 255, 0, 255],

[ 0, 255, 0, 255]

]

])

row = image_org.shape[1]

col = image_org.shape[2]

image_int1 = np.empty((row, col, 3), dtype=int)

image_int2 = np.empty((row, col, 3), dtype=int)

image_float1 = np.empty((row, col, 3), dtype=float)

image_float2 = np.empty((row, col, 3), dtype=float)

for m in range(row):

for n in range(col):

for rgb in range(3):

image_int1[m, n, rgb] = image_org[rgb, m, n]

image_int2[m, n, rgb] = image_org[rgb, m, n] /255

#-> all the pixels become BLACK

image_float1[m, n, rgb] = image_org[rgb, m, n]

image_float2[m, n, rgb] = image_org[rgb, m, n] / 255

#- >works, but

# "Clipping input data to the valid range for imshow with RGB data ([0..1] for floats or [0..255] for integers)."

fig, axs = plt.subplots(2, 2)

axs[0, 0].imshow(image_int1)

axs[0, 0].set_title("dtype=int, 0~255")

axs[0, 1].imshow(image_int2)

axs[0, 1].set_title("dtype=int, 0.0~1.0")

axs[1, 0].imshow(image_float1)

axs[1, 0].set_title("dtype=float, 0.0~255.0")

axs[1, 1].imshow(image_float2)

axs[1, 1].set_title("dtype=float, 0.0~1.0")

plt.show()

print("original image shape :{}".format(image_org.shape))

print("rearranged image shape:{}".format(image_int1.shape))

print("rearranged image:\n{}".format(image_int1))

# Clipping input data to the valid range for imshow with RGB data ([0..1] for floats or [0..255] for integers).

# original image shape :(3, 2, 4)

# rearranged image shape:(2, 4, 3)

# rearranged image:

# [[[ 0 0 0]

# [ 0 0 255]

# [ 0 255 0]

# [ 0 255 255]]

# [[255 0 0]

# [255 0 255]

# [255 255 0]

# [255 255 255]]]

imshow()に渡す配列のdtypeがint型の時は、ピクセルデータのレンジが0～255になる。

左上は元の配列のままR, G,Bが0か255なので、想定した組み合わせの色となっている
右上は想定されているレンジに対して0.0～1.0の値を与えていることから、どのピクセルともR, G, Bが0か1（ほぼゼロ）となり黒くなっている（そのまま実行され、特にメッセ維持は出ない）

配列のdtypeがfloatの時は、ピクセルデータの想定レンジは0.0～1.0になる。

左下は最小値0と最大255を与えているが、結果は左上と同じで、imshow()のデフォルトのレンジ0～255に変更されているようである（特にメッセージは出ない）
右下は与えるピクセルデータを0.0～1.0としたところ、”入力データをクリップしている”というメッセージが出たが、レンジが修正されたらしく結果は意図通り

並べ替えた後の配列は、直感的にはわかりにくい形になっている。

グレースケール

グレースケールの場合は、cmap='gray'を指定する。vmin、vmaxは省略しても同じ結果となるが念のため。

import numpy as np
import matplotlib.pyplot as plt

row_L = np.arange(0, 256, 8, dtype=int)
row_R = np.arange(256, 0, -8)
row = np.append(row_L, row_R)
n_size = 256 // 8

image_array = np.empty((n_size, n_size), dtype=int)

for r in range(n_size):
    image_array[r, :] = row[:n_size]
    row = np.roll(row, -1)

plt.imshow(image_array, vmin=0, vmax=255, cmap='gray')
plt.show()

import numpy as np

import matplotlib.pyplot as plt

row_L = np.arange(0, 256, 8, dtype=int)

row_R = np.arange(256, 0, -8)

row = np.append(row_L, row_R)

n_size = 256 // 8

image_array = np.empty((n_size, n_size), dtype=int)

for r in range(n_size):

image_array[r, :] = row[:n_size]

row = np.roll(row, -1)

plt.imshow(image_array, vmin=0, vmax=255, cmap='gray')

plt.show()

matplotlib – カラーマップ

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

matplotlibのカラーマップの一覧。cmap引数の指定のために。ソースはこちら。

../../_images/sphx_glr_colormaps_001.png

../../_images/sphx_glr_colormaps_002.png

../../_images/sphx_glr_colormaps_003.png

../../_images/sphx_glr_colormaps_004.png

../../_images/sphx_glr_colormaps_005.png

../../_images/sphx_glr_colormaps_006.png

../../_images/sphx_glr_colormaps_007.png

PCA – Boston house-pricesデータセット

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

概要

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

主成分が適切に発見されてよい相関が得られることを期待したが、IrisデータやBreast cancerデータの場合のようなクラス分類データにおける良好な結果は得られなかった。

ただし、Boston housing pricesデータはIrisやcancerのデータよりも複雑な社会行動に関するものであり、その指標も限定されていることから、これをもってPCAが回帰系のデータに不向きとまでは言い切れない。

なお、Boston housing pricesデータの特徴量には属性データ（カテゴリーデータ、クラスデータ）が含まれることから、DataFrameのget_dummis()メソッドによるone-hot encodingを行っている。

計算の手順

必要なパッケージをインポート
Boston housing pricesデータセットを準備
データセットをスケーリング／エンコーディング
1. 属性データの列を取り出して、get_dummiesでone-hot化
2. StandardScalerで残りの特徴量データを標準化
3. 上記2つを結合して前処理済みデータとして準備
PCAモデルのインスタンスを生成
- 引数n_components=2として、2つの特徴量について計算
fit()メソッドにより、モデルにデータを学習させる
主成分やその寄与率を確認
- 主成分はPCA.comonents_を、寄与率はPCA.explained_variance_ratio_を確認
transform()メソッドによって、主成分に沿ってデータを変換
3つの主成分について3次元可視化
2つの主成分について2次元可視化

前処理

特徴量のうちの1つCHASについては、「チャールズ川に関するダミー変数（1：川沿い、0：それ以外）」～”Charles River dummy variable (= 1 if tract bounds river; 0 otherwise)“となっていて、0か1の属性変数である。この変数をDataFrameのget_dummies()メソッドでone-hot化する。

また、その他のデータについてはStandardScalerで標準化する。

CHASのデータのみone-hot化
CHASの列を除いたデータをStandardScalerで標準化
上記2つのデータをjoin()で結合

df_chas_encoded = pd.get_dummies(df['CHAS'], prefix='CHAS')
df_wo_chas = df.drop('CHAS', axis=1)
df_wo_chas_scaled = pd.DataFrame(
    StandardScaler().fit_transform(df_wo_chas), columns=df_wo_chas.columns)
df_preprocessed = df_wo_chas_scaled.join(df_chas_encoded)

df_chas_encoded = pd.get_dummies(df['CHAS'], prefix='CHAS')

df_wo_chas = df.drop('CHAS', axis=1)

df_wo_chas_scaled = pd.DataFrame(

StandardScaler().fit_transform(df_wo_chas), columns=df_wo_chas.columns)

df_preprocessed = df_wo_chas_scaled.join(df_chas_encoded)

可視化

2次元

ここではまず、2次元可視化の結果を確認する。

クラス分類の場合は2次元で2つの主成分を確認できるが、回帰データの場合はターゲットの量を確認する必要があるため、グラフの軸を1つ消費する。このため、2次元による表現では1つの主成分による説明性を確認することになる。

各点の色や大きさをターゲットの値によって変化させ、2つの軸を2つの主成分に割り当てる方法も考えられるが、直感的にとらえにくくなる。

この結果を見る限り、あまり美しい結果とはなっていない。データを俯瞰した際、各特徴量であまりいい説明ができなかったが、その中でもある程度関係がみられたMDEVやLSTATとの相関と変わらないくらい。

3次元

そこで3次元の可視化にして、2つの主成分による説明性を確認する。

これでもあまりいい結果にならない。ただしグラフを見ると、大きく2つの塊に分かれているように見える。ターゲットである住居価格とは別に、特徴量の組み合わせに隠れている、性質の違うグループがあるのかもしれない。

主成分と寄与率

2つの主成分と寄与率について表示させてみる。

Components:
     feature    comp_0    comp_1
0       CRIM  0.251012  0.399530
1         ZN -0.256281  0.436197
2      INDUS  0.346630 -0.109264
3        NOX  0.342782 -0.169987
4         RM -0.189333  0.063747
5        AGE  0.313601 -0.318114
6        DIS -0.321462  0.331591
7        RAD  0.319816  0.381156
8        TAX  0.338513  0.318517
9    PTRATIO  0.205059  0.182159
10         B -0.203023 -0.333020
11     LSTAT  0.309817 -0.027564
12  CHAS_0.0 -0.001091  0.035239
13  CHAS_1.0  0.001091 -0.035239

Explained :[0.50514047 0.1109301 ]

Components:

feature comp_0 comp_1

0 CRIM 0.251012 0.399530

1 ZN -0.256281 0.436197

2 INDUS 0.346630 -0.109264

3 NOX 0.342782 -0.169987

4 RM -0.189333 0.063747

5 AGE 0.313601 -0.318114

6 DIS -0.321462 0.331591

7 RAD 0.319816 0.381156

8 TAX 0.338513 0.318517

9 PTRATIO 0.205059 0.182159

10 B -0.203023 -0.333020

11 LSTAT 0.309817 -0.027564

12 CHAS_0.0 -0.001091 0.035239

13 CHAS_1.0 0.001091 -0.035239

Explained :[0.50514047 0.1109301 ]

寄与率は第1主成分が50%程度で低いため高い相関が出ないといえるかもしれない。だが、Breast cancerデータセットのクラス分類では第1主成分の寄与率が40%台だが、明確なクラス分類ができていた。やはり回帰系の問題にはPCAは不向きなのかもしれない。

主成分の要素について、先の散布図が第1主成分と負の相関があることから、第1主成分の各特徴量は価格低下に寄与するものはプラス、価格上昇に寄与するものはマイナスとなるはずである。

たとえばZNやRMがマイナスなのは頷けるが、DISがマイナスなのは微妙。TAXやPRATIOがプラスなのも逆のような気がする。

先にも書いたが、Boston housing pricesデータセットで取りそろえられた特徴量は、住居の価格以外の何かを特徴づける傾向が強いのかもしれない。

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

概要

使い方

基本形

値が飛んでいる場合

順不同

weightsの意味

概要

画像ファイルの表示

配列の画像表示

基本形

レンジ

viminとvmax

RGB

グレースケール

概要

計算の手順

前処理

可視化

2次元

3次元

主成分と寄与率

概要

手順

主成分と寄与率

可視化

3次元

2次元

まとめ

概要

計算の手順

主成分と寄与率

可視化

3次元

2次元

主成分分析の特徴

概要

最大化すべき分散の導出

分散の最大化

特徴量が2つの場合

第2主成分以降

主成分の意味

概要

基本

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

属性が数値表現の場合

属性名の指定

概要

DataFrameの規模・形状

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

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

dscribe()～基本的な統計量

属性変数のカウント

1次元配列の場合

平均

分散・標準偏差

標本分散・標本標準偏差

不偏推定量

2次元配列の場合

平均

分散・標準偏差

標本分散・標本標準偏差

不偏推定量

共分散

`weights`の意味

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

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