Pythonで度数分布表を一発で自動生成する

# はじめに
数学・統計の分野で階級、階級値、度数、累積度数、相対度数、累積相対度数がセットの表を見ることがあると思います。どういうものかというと、こういうものです。

|階級|階級値|度数|累積度数|相対度数|累積相対度数|
|:--------:|-----:|---:|-------:|-------:|-----------:|
|0以上3未満 |   1.5|   1|       1| 0.07143|      0.0714|
|3以上6未満 |   4.5|   6|       7| 0.42857|      0.5000|
|6以上9未満 |   7.5|   2|       9| 0.14286|      0.6429|
|9以上12未満|  10.5|   2|      11| 0.14286|      0.7857|
|12以上15未満|  13.5|   3|      14| 0.21429|      1.0000|
|合計       |  -   |  14|      -|  1.00000|           -|

これをPythonで一発で出してくれる関数が意外と見つからないなということで作ってみました。

# 既存の便利関数
全部まとまった表を作成する関数はありませんが、以下の便利な関数で部分的に必要な情報を取り出せたりはします。それに加えて多少の計算を行うことで必要な値は全て揃います。

```py
# numpyのcumsum()で累積度数を取得する
data.cumsum()

# pandasのvalue_counts()で各値の出現頻度を数える
pd.Series(data).value_counts()
```

# 自動化のための工夫
#### 階級の数と階級の幅の決定
階級の数や階級の幅を決めるのに明確なルールはありません。しかし目安を知るための[スタージェスの公式](https://bellcurve.jp/statistics/glossary/7445.html)というものがあるのでそれを使います。
> **スタージェスの公式**<br>度数分布表やヒストグラムを作成するときに階級の数を決定する目安を得られる公式。Nをサンプルサイズ、kを階級数とすると、次のように計算することができる。階級の幅は、データの最小値から最大値をkで割り求める。
>
```math
k=log_2N+1
```

```py

# スタージェスの公式から階級の数を求める
class_size = 1 + np.log2(len(data))
class_size = int(round(class_size))

# 階級幅を求める
class_width = (max(data) - min(data)) / class_size # 分母は階級の数、分子は範囲。
class_width = round(class_width)
```

ただ、やはり明確なルールはなく、階級の幅を5など任意のキリがいい値にしたいこともあると思うのでそれにも対応できるようにします。
スタージェスの公式で出た値を使いたい場合は関数の第二引数に` None`を指定します。
任意の値を使いたい場合は第二引数にその任意の値を指定します。階級の数はそれに合わせて変更させます。

```py
def Frequency_Distribution(data, class_width):
    if class_width == None:
        # 階級幅を求める
        class_width = (max(data) - min(data)) / class_size #分母は階級の数、分子は範囲。
        class_width = round(class_width) # 四捨五入
    else:
        class_width = class_width
        class_size = max(x) // class_width
```

#### 階級のインデックスを動的に変更
階級は「〜以上・・・未満」というやつです。度数分布表を作成するにあたり、インデックスとして階級を設定して表に記述しておきたいわけですが、入力したデータに合わせて手動で打ち込んでいては大変です。
そこで、階級の幅と階級の数、フォーマット演算子を用いてリスト型内包表記でfor文を回すことにより、インデックスのラベルが生成できます。

```py

class_width = 5 # 階級の幅
class_size = 10 # 階級の数
['%s以上%s未満'%(w, w+class_width) for w in range(0, class_size*class_width*2, class_width)]

# ['0以上5未満','5以上10未満','10以上15未満','15以上20未満','20以上25未満','25以上30未満']
```

#### 表の作成
あとはpandasの操作で行や列の追加とカラム名やインデックス名の更新などの修正を加えるだけです。

# コード全体

```py
import pandas as pd
import numpy as np

# 度数分布表を作る
def Frequency_Distribution(data, class_width):
    # スタージェスの公式から階級の数を求める
    class_size = 1 + np.log2(len(data))
    class_size = int(round(class_size))
    if class_width == None:
        # 階級幅を求める
        class_width = (max(data) - min(data)) / class_size # 分母は階級の数、分子は範囲。
        class_width = round(class_width) # 四捨五入
    else:
        class_width = class_width
        class_size = max(x) // class_width
    # print('階級の数:', class_size)
    # print('階級幅:', class_width)
    
    # 階級に振り分ける
    # 各観測値を階級値にする
    cut_data = []
    for row in data:
        cut = row // class_width
        cut_data.append(cut)
        
    #頻度を数える
    Frequency_data = pd.Series(cut_data).value_counts()
    Frequency_data = pd.DataFrame(Frequency_data)
    #インデックスでソートし、任意の位置に行を挿入したいので一旦転置
    F_data = Frequency_data.sort_index().T
    # 度数0の階級があればデータフレームに挿入する
    for i in range(0, max(F_data.columns)):
        if (i in F_data) == False:
            F_data.insert(i, i, 0)
    F_data = F_data.T.sort_index()
    #インデックスとカラムの名前を変える
    F_data.index = ['%s以上%s未満'%(w, w + class_width) for w in range(0, class_size * class_width * 2, class_width)][:len(F_data)]
    F_data.columns = ['度数']

    F_data.insert(0, '階級値', [((w + (w + class_width)) / 2) for w in range(0, class_size * class_width * 2, class_width)][:len(F_data)])
    F_data['累積度数'] = F_data['度数'].cumsum()
    F_data['相対度数'] = F_data['度数'] / sum(F_data['度数'])
    F_data['累積相対度数'] = F_data['累積度数'] / max(F_data['累積度数'])
    F_data.loc['合計'] = [None, sum(F_data['度数']), None, sum(F_data['相対度数']), None]

    return F_data

# サンプルデータ
x = [0, 3, 3, 5, 5, 5, 5, 7, 7, 10, 11, 14, 14, 14]
Frequency_Distribution(x, None)
```
#### 結果
|階級|階級値|度数|累積度数|相対度数|累積相対度数|
|:--------:|-----:|---:|-------:|-------:|-----------:|
|0以上3未満 |   1.5|   1|       1| 0.07143|      0.0714|
|3以上6未満 |   4.5|   6|       7| 0.42857|      0.5000|
|6以上9未満 |   7.5|   2|       9| 0.14286|      0.6429|
|9以上12未満|  10.5|   2|      11| 0.14286|      0.7857|
|12以上15未満|  13.5|   3|      14| 0.21429|      1.0000|
|合計       |  -   |  14|      -|  1.00000|           -|

# 補足(2020/10/4)
@nkayさんからコメントでいただいたこちらのコードの方が非常にスマートに書けるので推奨です。

```py

def Frequency_Distribution(data, class_width=None):
    data = np.asarray(data)
    if class_width is None:
        class_size = int(np.log2(data.size).round()) + 1
        class_width = round((data.max() - data.min()) / class_size)

    bins = np.arange(0, data.max()+class_width+1, class_width)
    hist = np.histogram(data, bins)[0]
    cumsum = hist.cumsum()

    return pd.DataFrame({'階級値': (bins[1:] + bins[:-1]) / 2,
                         '度数': hist,
                         '累積度数': cumsum,
                         '相対度数': hist / cumsum[-1],
                         '累積相対度数': cumsum / cumsum[-1]},
                        index=pd.Index([f'{bins[i]}以上{bins[i+1]}未満'
                                        for i in range(hist.size)],
                                       name='階級'))


x = [0, 3, 3, 5, 5, 5, 5, 7, 7, 10, 11, 14, 14, 14]
Frequency_Distribution(x)
```

# 参考
上記のコード作成にあたり、主に以下のサイトを参考にさせていただきました
[いっかくのデータサイエンティストをいく](http://midsum-datasience.hatenablog.com/entry/2018/06/08/214710)
[統計用語集](https://bellcurve.jp/statistics/glossary/7445.html)