Python
python2.7

Pythonでデータの分析を出来るようになりたい(その1)

More than 1 year has passed since last update.

Pythonでデータの分析を出来るようになりたい

ちなみに Python 初心者です。
Pythonは文法がすっきりしていて、ライブラリも豊富でパワフルだという評判を聞いてたんで、年末休みを利用してちょっとかじってみることにした。

題材

何事も、勉強するなら具体的な題材・問題があった方がいい。
そこで最近読み始めた統計に関する本の内容をベースにすることにした。

統計学が最強の学問である[実践編]データ分析のための思想と方法 - 西内啓 著

第1章 統計学の実践は基本の見直しから始まる
 05 なぜ、平均値は真実を捉えることができるのか?

この章に「裏が出る確率が2/3で表が出る確率が1/3、というコイン」を使った例がいくつかある。そこで Python を使って実際に「コイン」を何回も投げて(もちろんシミュレーションだけど)例と同じ結果になるか確認してみることにした。
(詳しい内容は本を買うか図書館で借りて調べるかしてください。P59~P64の部分です)

Pythonのバージョン

使ったPython のバージョンは2.7。夏に Python のセミナーみたいなのを受けたときに言われるがままに入れた Spyderってのがインストールされてたんで、それを使った。

コインを2回投げた時の確率

まずは、コインを2回投げるときの全ての組み合わせの確率についての例をPythonで実験してみた。

可能な組み合わせは:
裏・裏 (表0枚)
裏・表 (表1枚)
表・裏 (表1枚)
表・表 (表2枚)

Python コード

from random import randint
from decimal import Decimal
from prettytable import PrettyTable
import numpy as np

def tossBiasedCoin():
    """ Returns 0 or 1 with 0 having 2/3 chance """
    return randint(0,2) % 2

# Make a 2x2 array
counts = [[0 for j in range(2)] for i in range(2)]

# Toss a coin many times to get counts
sampleCount = 500000
for num in range(sampleCount):
    first = tossBiasedCoin()
    second = tossBiasedCoin()
    counts[first][second] += 1

# Conert all counts to perentage
TWOPLACES = Decimal(10) ** -2 
for i in range(2):
    for j in range(2):
        value = counts[i][j]        
        counts[i][j] = (100 * Decimal(counts[i][j])/Decimal(sampleCount)).quantize(TWOPLACES)
        print("Converted the value {} to percentage {}".format(value, counts[i][j]))

# Make summaries of number of heads.
keys = np.arange(3)
values = [counts[0][0], 
          counts[0][1]+counts[1][0],
          counts[1][1]]

# Add row descriptions
counts[0].insert(0, '1st tail')
counts[1].insert(0, '1st head')

# Create table with column descriptions, add rows, then show it.
table = PrettyTable(["", "2nd tail", "2nd head"])
table.padding_width = 1
table.add_row(counts[0])
table.add_row(counts[1])
print table

# Draw a bar chart
import matplotlib.pyplot as plt
fig, ax = plt.subplots()
rects = plt.bar(keys,
                 values, 
                 0.5,
                 alpha=0.4,
                 align="center", 
                 color='b')

plt.xlabel('Number of heads')
plt.ylabel('Probability (%)')
plt.title('Probabilities heads with a biased coin')
plt.xticks(keys, np.arange(3))

plt.tight_layout()
plt.show()

実行結果

まずはP59の図表1-17とほぼ同じ表が表示される。4通りの組み合わせと、それぞれの確率。
(この表はPrettyTableを使った。後述)
image

「表・表」の組み合わせが出にくいのは、そもそも表が出る確率が1/3というコインを使ってるから。

つぎにp60の図表1-18とほぼ同じ棒グラフが表示される。

image

これは4通りのパターンを「表が何枚出るか」で3つにまとめたもの。「表が1枚」がパターンが(表・裏)と(裏・表)の2通りあったから、それは一つにまとめられている。

コードの内容

ライブラリのインポート

使用する関数などを実行環境に import する。
便利なライブラリは Python 標準ライブラリ にあるかもしれないが、膨大な数の third-party ライブラリが Python パッケージインデックス というレポジトリにある。

from random import randint
from decimal import Decimal
from prettytable import PrettyTable
import numpy as np
  • 乱数生成にrandomモジュールの randint を使う。
  • float型を下2桁に整えるのにdecimalモジュールの Decimal を使う。
  • 表を作るにはprettytableモジュールの PrettyTable を使う。
  • arange という関数が便利なのでnumpyモジュールをインポートするが、特にヘビーな使い方はしてない。

Pythonでの関数定義

def tossBiasedCoin():
    """ Returns 0 or 1 with 0 having 2/3 chance """
    return randint(0,2) % 2

関数にするほどの内容じゃないけど、関数の定義の練習としてコインの表・裏(1か0)を返す関数を作った。0, 1, 2 のいずれかの値を乱数で生成して、値が偶数なら0を返し、それ以外なら1を return する。三つの値のうちふたつが偶数だから2/3の確率になる。

2x2の変数を準備

出現頻度を記録するために2x2の sequence を使う。この場合は list 型になる。

# Make a 2x2 array
counts = [[0 for j in range(2)] for i in range(2)]

それぞれの変数は0に初期化しておく。
for 文 を、組み込み関数の一つである range関数 で作った即席の list の上でループ実行させている。2x2にするために、list を要素として持つ list を作る。

コインを投げて結果を記録する

ここでは50万回投げるけど、そんなに投げなくてもいいかな(笑)

# Toss a coin many times to get counts
sampleCount = 500000
for num in range(sampleCount):
    first = tossBiasedCoin()
    second = tossBiasedCoin()
    counts[first][second] += 1

結果は0か1なので、ちょうど2x2構造に対してのインデックスとして使用できる。インデックスしたマスの数を1づつ増やしていく。

ちなみにPythonではC言語でおなじみの ++ オペレータは無いもよう。Pythonのオペレータのリストはこちら

パーセントに変換する

頻度を投げた総数で割って、割合を出す。

# Convert all counts to perentage
TWOPLACES = Decimal(10) ** -2
for i in range(2):
    for j in range(2):
        value = counts[i][j]
        counts[i][j] = (100 * Decimal(counts[i][j])/Decimal(sampleCount)).quantize(TWOPLACES)
        print("Converted the value {} to percentage {}".format(value, counts[i][j]))

2x2 構造のそれぞれのマスを訪れることになるんで、iとjの二つのインデックスを使ってループ実行し、[i][j]の形で任意のマスの値にアクセスする。値はin-placeで置き換えるけど、デバッグ用に変換前と変換後の値を表示している。

Decimal.quantize 関数に0.01を渡すことで値を下2桁に丸める。

棒グラフ用にデータを準備

棒グラフの棒は3つ。表が0枚、表が1枚、そして表が2枚。

# Make summaries of number of heads.
keys = np.arange(3)
values = [counts[0][0],
          counts[0][1]+counts[1][0],
          counts[1][1]]

表が1枚の頻度だけは、表がコイン1投目で出たか2投目で出たかに関係ないのでその両方の値を合わせている。

表用にデータを準備

表の行となる list の左側に list.insert を使って「1投目が裏」と「1投目が表」という欄を加える。

# Add row descriptions
counts[0].insert(0, '1st tail')
counts[1].insert(0, '1st head')

表を作る

PrettyTableというサード・パーティ・ライブラリ関数を使う。

# Create table with column descriptions, add rows, then show it.
table = PrettyTable(["", "2nd tail", "2nd head"])
table.padding_width = 1
table.add_row(counts[0])
table.add_row(counts[1])
print table

棒グラフを作る

サード・パーティー製ライブラリ、matplotlibを使う。この matplotlib はそれだけで本が一冊書けそうなぐらい充実しているみたい(ギャラリー参照)
matplotlib の中の pyplot棒グラフを書く。

# Draw a bar chart
import matplotlib.pyplot as plt
fig, ax = plt.subplots()
rects = plt.bar(keys,
                 values,
                 0.5,
                 alpha=0.4,
                 align="center",
                 color='b')

plt.xlabel('Number of heads')
plt.ylabel('Probability (%)')
plt.title('Probabilities heads with a biased coin')
plt.xticks(keys, np.arange(3))

plt.tight_layout()
plt.show()

(その2)につづく