Pythonによる統計分析入門

平均・分散・相関係数・グラフ描画（ヒストグラム、散布図）・ファイル入出力

簡単な統計量とグラフ描画のプログラムをPythonで作ってみた。１変量および２変量の記述統計量、グラフ描画、ファイル入出力の３セッションに分けて説明する。

なお、Pythonによるデータ分析の解説書として、岡本安晴（著）「いまさら聞けないPythonでデータ分析―多変量解析，ベイズ統計分析（PyStan，PyMC）―」丸善出版を用意した。

本ウェブサイトのプログラム例は圧縮ファイルscriptfiles.zipにまとめたが、ファイル名scriptfiles.zipを右ボンクリックによりダウンロード・展開して、自由に御利用頂きたい

記述統計量

１変量の場合

表１の５０個の値の平均値と分散などを求めてみる。

表１　５０個のスコア

66, 59, 62, 64, 63,68, 65, 59, 68, 64,

65, 51, 67, 64, 83,59, 61, 62, 57, 72,

65, 64, 54, 60, 53,65, 67, 60, 53, 79,

74, 53, 61, 68, 75,50, 57, 55, 66, 56,

55, 61, 70, 71, 49,69, 70, 80, 73, 72

プログラムStatOneVar.pyを以下のように作成した。

＝＝＝＝＝＝＝＝＝＝＝＝ StatOneVar.py ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

import numpy as np

data = np.array( # Numpyの配列型ndarray型として用意する

[66, 59, 62, 64, 63,

68, 65, 59, 68, 64,

65, 51, 67, 64, 83,

59, 61, 62, 57, 72,

65, 64, 54, 60, 53,

65, 67, 60, 53, 79,

74, 53, 61, 68, 75,

50, 57, 55, 66, 56,

55, 61, 70, 71, 49,

69, 70, 80, 73, 72])

n = len(data) # データ数

print("n =%d"%n)

print("Data...")

i = 1

for x in data:

print(f"{i}: {x}")

i += 1

v_sum = np.sum(data)

mean = v_sum / n

print("Mean = %f"%mean)

v_ssum = np.sum((data - mean)**2)

var = v_ssum / n

sd = var**0.5

print(f"Variance = {var:.3f} sd = {sd:.3f}")

uvar = v_ssum / (n - 1.0)

sd_uvar = uvar**0.5

print(f"Unbiased Var. = {uvar:.3f} sd(uvar) = {sd_uvar:.3f}")

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

表１のデータをリストdataとして設定（assignment記号 = の働きについては、本節末に補足説明を用意した。）している。平均値、分散、不偏分散、および分散あるいは不偏分散の平方根としての標準偏差が算出されている。それぞれは、次式

、、、

、

で算出される。

プログラムStatOneVar.pyを実行すると、図１のような結果になる。

図１

２変量の場合

体重と身長のように２つの変量がペアになっている場合について考える。

表２　［体重, 身長］の５０人分のデータ

[70, 168], [68, 165], [62, 159], [56, 156], [70, 170],

[64, 162], [59, 162], [70, 169], [68, 168], [69, 170],

[67, 168], [63, 166], [58, 163], [74, 170], [65, 168],

[60, 160], [65, 165], [72, 175], [61, 164], [77, 174],

[68, 164], [64, 161], [70, 166], [47, 147], [66, 165],

[68, 170], [73, 174], [74, 172], [63, 164], [83, 182],

[58, 159], [64, 166], [67, 167], [75, 174], [67, 163],

[68, 166], [68, 167], [66, 165], [63, 165], [55, 154],

[61, 163], [60, 160], [62, 164], [62, 160], [59, 160],

[51, 155], [59, 163], [68, 170], [56, 157], [79, 180]

プログラムStatTwoVars.pyでは、表２のデータがnumpyの配列ndarray型pairsとして設定されている。

＝＝＝＝＝＝＝＝＝＝＝＝＝ StatTwoVars.py ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

import numpy as np

pairs = np.array(

[[70, 168], [68, 165], [62, 159], [56, 156], [70, 170],

[64, 162], [59, 162], [70, 169], [68, 168], [69, 170],

[67, 168], [63, 166], [58, 163], [74, 170], [65, 168],

[60, 160], [65, 165], [72, 175], [61, 164], [77, 174],

[68, 164], [64, 161], [70, 166], [47, 147], [66, 165],

[68, 170], [73, 174], [74, 172], [63, 164], [83, 182],

[58, 159], [64, 166], [67, 167], [75, 174], [67, 163],

[68, 166], [68, 167], [66, 165], [63, 165], [55, 154],

[61, 163], [60, 160], [62, 164], [62, 160], [59, 160],

[51, 155], [59, 163], [68, 170], [56, 157], [79, 180]])

n = len(pairs)

print("n = %d"%n)

sum_x = np.sum(pairs.T[0])

sum_y = np.sum(pairs.T[1])

mean_x = sum_x / n

mean_y = sum_y / n

print(f"Mean_x = {mean_x:.3f}")

print(f"Mean_y = {mean_y:.3f}")

ssum_x = np.sum((pairs.T[0] - mean_x)**2)

ssum_y = np.sum((pairs.T[1] - mean_y)**2)

sum_xy = np.sum((pairs.T[0] - mean_x) * (pairs.T[1] - mean_y))

var_x = ssum_x / n

var_y = ssum_y / n

print(f"Var_x = {var_x:.3f}")

print(f"Var_y = {var_y:.3f}")

sd_x = var_x ** 0.5

sd_y = var_y ** 0.5

print(f"sd_x = {sd_x:.3f}")

print(f"sd_y = {sd_y:.3f}")

cor = sum_xy / ((ssum_x * ssum_y)**.5)

print("r = %f"%cor)

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

２変量の関係を表すものとして相関係数（correlation coefficient）が求められている。これは、次式

によって算出されている。

プログラムStatTwoVars.pyを実行すると図２のような結果になる。

図２

補足説明：Assignment = について

Assignment = については、以下のコードに示されるような注意が必要である。

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

a = [1, 2]

b = a # Check-11

b[0] = 3 # Check-12

print("Check-1")

print(a)

print(b)

a1 = [1, 2]

b1 = a1 # Check-21

a1[0] = 5 # Check-22

print("Check-2")

print(b1)

print(a1)

a2 = [1, 2]

b2 = a2 # Check-31

b2 = [6, 7] # Check-32

print("Check-3")

print(a2)

print(b2)

a3 = 10

b3 = a3

b3 = 15

print("Check-4")

print(a3)

print(b3)

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

上のコードを実行すると以下の結果を得る。

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

Check-1

[3, 2]

Check-2

[5, 2]

Check-3

[1, 2]

[6, 7]

Check-4

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

「Check-1」に続く出力では、同じ出力「[3, 2]」が得られている。これは、Check-11における代入文によって、変数bが変数aと同じものを指し示すものとなっていることを意味している。

「Check-2」に続く出力も、代入がb1ではなくa1（b1 = a1の右辺）の要素に行われているが、2つの変数b1とa1が同じものを表していることを意味していることは変数aとbの場合と同様であることを意味するものである。

「Check-3」に続く出力では、変数b2の値の変更は変数a2には及んでいない。これは、変数の値全体を変更する場合は、変数は新しい変更値を指し示すものとして設定され、それまでの同じものを指し示していた変数a2から切り離されたことを意味している。

「Check-4」に続く出力は整数の場合を扱っているが、整数はそれ１つで変数の表す値全体であるので、変数b3に値を設定すると、それまで同じものを表していた変数a3とは別のものを指し示すことになる。

グラフ描画

ヒストグラム

上の表１のデータに対して、ヒストグラムを描くプログラムHistogram.pyを作成してみた。描画プログラムではパッケージmatplotlibを用いている（importしている）が、これはPythonに追加インストールするものである。

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝ Histogram.py ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

import matplotlib.pyplot as plt

data = [66, 59, 62, 64, 63,

68, 65, 59, 68, 64,

65, 51, 67, 64, 83,

59, 61, 62, 57, 72,

65, 64, 54, 60, 53,

65, 67, 60, 53, 79,

74, 53, 61, 68, 75,

50, 57, 55, 66, 56,

55, 61, 70, 71, 49,

69, 70, 80, 73, 72]

n = len(data) # データ数

print("n =", n)

print("Data...")

for i, x in enumerate(data):

print(f"{i+1}: {x}")

plt.title("Histogram") # ヒストグラムのタイトル

plt.xlabel("Score") # X軸のラベル

plt.ylabel("Frequency") # Y軸のラベル

results = plt.hist(data, bins = 10) # ヒストグラムの作成、bins（級数）を10に指定

print(results)

plt.show() # ヒストグラムの表示

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

上のスクリプトを実行すると、図３のヒストグラムが描かれる。

図３

画面下のフロッピーディスクのアイコンのボタンをクリックすると、表示画像を保存することができる。右上のＸ印のクリックで、表示画面は閉じられ、show()の実行終了となる。

端末には、次のスクリプト

results = plt.hist(data, bins = 10) # ヒストグラムの作成、bins（級数）を10に指定

print(results)

の出力が、図４のように表示されている。

図４

１番目の配列の要素は、ヒストグラムの各bin（級）内の度数を表し、２番目の配列は級の境界値を表している。

境界値を切りの良い値とするため、binsを次のように設定する。

bins = np.arange(45, 85.1, 5)

np.arangeは、

np.arange(a, b, c)

の形で呼び出すと、aから始まり、間隔cでbより小さい範囲で数値の配列が設定される。

図３のヒストグラムは、binの関係を調整した方が良さそうである。この調整を行ったスクリプトHistogram_bins.pyを以下のように作成してみた、

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝ Histogram_bins.py ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

import matplotlib.pyplot as plt

import numpy as np

data = [66, 59, 62, 64, 63,

68, 65, 59, 68, 64,

65, 51, 67, 64, 83,

59, 61, 62, 57, 72,

65, 64, 54, 60, 53,

65, 67, 60, 53, 79,

74, 53, 61, 68, 75,

50, 57, 55, 66, 56,

55, 61, 70, 71, 49,

69, 70, 80, 73, 72]

n = len(data) # データ数

print("n =", n)

print("Data...")

for i, x in enumerate(data):

print(f"{i+1}: {x}")

plt.title("Histogram") # ヒストグラムのタイトル

plt.xlabel("Score") # X軸のラベル

plt.ylabel("Frequency") # Y軸のラベル

results = plt.hist(data, bins = np.arange(45, 85.1, 5))

print(results)

plt.show() # ヒストグラムの表示

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

上のスクリプトを実行すると、図５のヒストグラムを得る。

図５

級境界値が切りの良い値になっていることが分かる。

散布図

散布図を描くプログラムScatter.pyを以下のように作成してみた。グラフで日本語表記を行うためにjapanize_matplotlibを用いている。(注：現在2025.06.15、Python 3.13.5ではjapanize_matplotlibはエラーが出る。Python3.13.2/miniconda/Windowsで実行した。)

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝ Scatter.py ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

import numpy as np

import scipy.stats as ss

import matplotlib.pyplot as plt

import japanize_matplotlib

# Data

pairs = np.array(

[[70, 168], [68, 165], [62, 159], [56, 156], [70, 170],

[64, 162], [59, 162], [70, 169], [68, 168], [69, 170],

[67, 168], [63, 166], [58, 163], [74, 170], [65, 168],

[60, 160], [65, 165], [72, 175], [61, 164], [77, 174],

[68, 164], [64, 161], [70, 166], [47, 147], [66, 165],

[68, 170], [73, 174], [74, 172], [63, 164], [83, 182],

[58, 159], [64, 166], [67, 167], [75, 174], [67, 163],

[68, 166], [68, 167], [66, 165], [63, 165], [55, 154],

[61, 163], [60, 160], [62, 164], [62, 160], [59, 160],

[51, 155], [59, 163], [68, 170], [56, 157], [79, 180]])

n = len(pairs) # The number of pairs (xi, yi)

print("n = %d"%n)

cor = ss.pearsonr(pairs.T[0], pairs.T[1]).statistic

print(f"r = {cor:.3f}")

plt.title(f"散布図...r = {cor:.3f}", fontsize=16)

plt.xlabel("体重", fontsize=14)

plt.ylabel("身長", fontsize=14)

plt.plot(pairs.T[0], pairs.T[1], 'bo')

plt.show()

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

プログラムScatter.pyを実行すると図６の散布図が表示される。

図６

ファイル入出力

テキストファイルの入出力について説明する。他のファイル、CSVファイルやExcelファイルについては、別のウェブサイトで扱っている。

入力データは、プログラム中に書き込むよりは別ファイルとして用意する方が便利である。また、出力もコンソールではなく、ファイルに出力した方が後で使いやすい。

１変量データの場合

表１のデータを、DataOneVar.txtとして以下のようなテキストファイルDataOneVar.txtで用意した。このテキストファイルは圧縮ファイルscriptfiles.zipに含まれている。

＝＝＝＝＝＝＝＝＝＝＝＝ DataOneVar.txt ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

・

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

データは/（半角スラッシュ）で始まる行の次の行から、１行に１つずつ書かれている。

最後のデータの次の行は、/（半角スラッシュ）で始まる行としている。

ファイル入出力に対応したプログラムStatOneVarF.pyを、以下のように作成した。

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝ StatOneVarF.py ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

import numpy as np

import matplotlib.pyplot as plt

# 入力データファイルの準備

s_in = input("Input data file = ")

f_in = open(s_in, "r")

# ファイル全体を読み込み、ファイルの各行を文字列で表してリストとする

data_f = f_in.readlines()

f_in.close()

# 出力用テキストファイルの準備

s_out = input("Output data file = ")

f_out = open(s_out, "w")

f_out.write("Input data file = " + s_in + "\n\n")

print('data_f =\n', data_f)

pos = 0

while data_f[pos][0] != '/': # 先頭が/である文字列まで進む

pos += 1

data = []

pos += 1

while data_f[pos][0] != '/': # 先頭が/である文字まで処理を進める

data.append(data_f[pos])

pos += 1

print('data =\n', data)

data = np.array(data, dtype=float) # 文字列の表す数値に変換してnumpyの配列型とする

print('data =\n', data)

d_mean = np.mean(data)

d_std = np.std(data)

print('mean =', d_mean, ' sd =', d_std)

for v in data:

f_out.write(f'{v}\n')

f_out.write(f'mean = {d_mean:.3f} sd = {d_std:.3f}\n')

plt.title("Histogram")

plt.xlabel("Scores")

plt.ylabel("Frequencies")

plt.hist(data)

plt.show()

f_out.close()

print(f'\n出力ファイル{s_out}が保存されました。\n')

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

次のコードによって、コンソールから入力した名前のファイルを入力ファイルとしている。

s_in = input("Input data file = ")

f_in = open(s_in, "r")

出力ファイルは以下のコードによって設定している。

s_out = input("Output data file = ")

f_out = open(s_out, "w")

入力データファイルの内容は、コード

data_f = f_in.readlines()

によってリストdata_fに設定されている。入力データファイルの構成（これはリストdata_fに文字列のリストとして取り込まれている）は、データの書かれている部分がスラッシュ / で始まる行で挟まれているので、以下のコードによってデータ（data_fの内容）をNumpyの配列型dataに設定している。

pos = 0

while data_f[pos][0] != '/': # 先頭が/である文字列まで進む

pos += 1

data = []

pos += 1

while data_f[pos][0] != '/': # 先頭が/である文字まで処理を進める

data.append(data_f[pos])

pos += 1

print('data =\n', data)

data = np.array(data, dtype=float) # 文字列の表す数値に変換してnumpyの配列型とする

プログラムを実行すると、まず、入力データファイル名を聞いてくる（図９）。

図７

入力データファイル名の設定後、Enterキーを押すと、次に出力ファイル名を聞いてくる（図８）。

図８

出力ファイル名として適当なテキストファイル名を設定してEnterキーを押すと計算が始まり、図９の画面が表示される。

図９

プログラムの実行終了後、出力ファイルを開くと、入力データファイルDataOneStat.txtの場合は、出力ファイルは以下のようになっている。

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝ ResultsOneVar.txt ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

Input data file = DataOneVar.txt

66.0

59.0

62.0

　・

80.0

73.0

72.0

mean = 63.680 sd = 7.760

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

２変量データの場合

表２の２変量（体重、身長）のデータを、以下のようにテキストファイルDataTwoVars.txtとして作成した。このファイルも圧縮ファイルscriptfiles.zipに含まれている。

＝＝＝＝＝＝＝＝＝＝＝＝＝ DataTwoVars.txt ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

70 168

68 165

62 159

56 156

70 170

64 162

59 162

・

62 160

59 160

51 155

59 163

68 170

56 157

79 180

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

データがスラッシュ / で始まる行で挟まれているので、読み込み方は１変量のプログラムStatOneVarF.pyの場合と基本的に同じである。２変量用のプログラムScatterF.pyを以下のように作成した。

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝ ScatterF.py ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

import numpy as np

import scipy.stats as ss

import matplotlib.pyplot as plt

# 入力データファイルの準備

s_in = input("Input data file = ")

f_in = open(s_in, "r")

# ファイル全体を読み込み、ファイルの各行をそれぞれ文字列としてリストに格納

data_f = f_in.readlines()

f_in.close()

# 出力テキストファイルの「準備

s_out = input("Output data file = ")

f_out = open(s_out, "w")

f_out.write("Input data file = " + s_in + "\n\n")

pos = 0

while data_f[pos][0] != '/': # 先頭が/で始まる行まで進む

pos += 1

data = []

while data_f[pos][0] != '/': # 先頭が/で始まる行まで処理を進める

data.append(data_f[pos])

pos += 1

v_data = []

for v in data:

v_data.append(v.split()) # 空白文字で区切られた文字列を空白文字で分解してリストとして

# 　リストv_datanoの要素として追加

data = np.array(v_data, dtype=float) # v_dataの要素をその表す数値に変換してnumpyの配列型とする

for v in data:

print(v)

f_out.write(f'{v}\n')

cor = ss.pearsonr(data.T[0], data.T[1]).statistic

f_out.write(f'\nr = {cor:.5f}')

plt.title(f"Scattergram...r = {cor:.3f}", fontsize=16)

plt.xlabel("Weight")

plt.ylabel("Height")

plt.plot(data.T[0], data.T[1], 'g^') # 'g^' -> green triangle

plt.show()

f_out.close()

print(f'\n出力テキストファイル{s_out}が保存されました\n')

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

１行に１組のデータ、ＸとＹ、が含まれているので、それを分離してNumpyの配列型datasに格納するため、以下のスクリプトが用いられている。

v_data = []

for v in data:

v_data.append(v.split()) # 空白文字で区切られた文字列を空白文字で分解してリストとして

# 　リストv_datanoの要素として追加

data = np.array(v_data, dtype=float) # v_dataの要素をその表す数値に変換してnumpyの配列型とする

また、上のプログラムScatterF.pyでは、散布図の点が緑色の三角形で描かれるように、第３パラメータの値が「'g^’」となっている。

plt.plot(data.T[0], data.T[1], 'g^')

このとき、散布図は図10のように描かれる。

図10

出力ファイを図11に示されている名前ResultsTwoVars.txtで設定すると、

図11

その内容は以下に示すようになっている。

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝ ResultsTwoVars.txt ＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

Input data file = DataTwoVars.txt

[ 70. 168.]

[ 68. 165.]

[ 62. 159.]

　　・

[ 68. 170.]

[ 56. 157.]

[ 79. 180.]

r = 0.94281

＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝

Home