Empty

2013年4月から社会人になりました。

マンガでわかる統計学を読んだ

概要

最近は就活の時期ということもあり、人と会うことが多くなり、どんなことをしたいのか? と聞かれるたびに、「データマイニング機械学習の技術を使って世の中をより良くしたい。」といったことを良く言及している。
しかし、いざ専門的な知識を持っているのかと聞かれると自信を持ってあるとは言えない自分がいる、かつ基礎的な知識についても怪しい節を感じている始末。
で、これではイカンっ!と思って、前々から読みたかったけど積んでいるままであった 「マンガでわかる統計学」シリーズを手に取って再勉強したので、そのメモなどをエントリーする。

読んだ本

今回は上の本のまとめ。回帰分析編はまた後にまとめる予定。

目次

第1章 データの種類をたしかめよう!
第2章 データ全体の雰囲気をつかもう! (数量データ編)
第3章 データ全体の雰囲気をつかもう! (カテゴリーデータ編)
第4章 基準値と偏差値
第5章 確立を求めよう!
第6章 2変数の関連を調べよう!
第7章 独立性の検定をマスターしよう!

カテゴリーデータと数量データ (第1章あたり)

カテゴリーデータ
測れないデータ
選択肢の幅が等間隔でないもの。

e.g.

  • Q. ご感想
    • ["おもしろかった", "ややおもしろかった", "ふつう", "ややつまらなかった", "とてもつまらなかった"]
  • Q. 英検の級
    • ["1級", "準1級", "2級", ... ]
      • 級間の対象語数や難易度の差は等しくないのでカテゴリーデータとして扱うのが一般的。
数量データ
測れるデータ
ある目盛りと隣の目盛りが等間隔であるデータ。標準化されているもの?(上手く定義できない...)

e.g.

  • 身長、体重、実行時間、とか
カテゴリーデータを数値データとして扱う (第1章あたり)

扱いやすいという理由から、カテゴリーデータを数量データとしてみなす場合も多い。
e.g.

とてもおもしろかった 5 points
ややおもしろかった 4 points
どちらともいえない 3 points
ややつまらなかった 2 points
とてもつまらなかった 1 points

度数分布表 & ヒストグラム (第2章~4)

用語
階級
ある範囲で区切ったときのそれぞれのカテゴリー e.g. 300円以上400円未満の商品
階級値
ある階級の値域の中央値 e.g. 350(円)
度数
ある階級に所属している数 e.g. 15(個)
相対度数
全体を1としたときの割合 e.g. 全体:100(個) ある階級の度数:15(個) その階級の相対度数: 0.15
ヒストグラム
底辺を階級の幅に、面積が度数を表すように高さを定めた長方形を密着して並べて作った統計グラフの一種。なんとなーくデータ全体の雰囲気を直観的に掴むことができる.
平均(arithmetric mean)
相加平均・算術平均を指すことが多い。全データの合計/データ数 で求まる。 (その他に、幾何平均、調和平均といった違う概念を指す場合もある。)
中央値(median)
有限個のデータを小さい順位並べたとき、中央に位置する値。データ数が偶数個の場合は、真ん中2つの算術平均をとる。
標準偏差(standard deviation)
散らばりの程度を表す。(1データあたりの平均からのズレ) 0以上。母集団の標準偏差と、標本の標準偏差を求める式は若干異なるが、ほとんどの場合は母集団全体を調べられないので、標本の標準偏差を用いる場合が多い。
階級の幅を決め方
  • スタージェスの公式を用いて幅を決めるやり方が紹介されていた。
  • そのほかの方法、[入門解析法 永田先生 p.9]にも載ってた.
コードに起こす (理解度チェック)

本の内容ではないけど、理解しているならコードに起こせるよね?という考えで組み込み関数を極力使わずに、
本の内容を計算するためのスクリプトを書いた。

#!/usr/bin/env python
# coding: utf-8

import math
import random
from prettytable import PrettyTable


def struges_fomula(dataset):
        """スタージェスの公式"""
        category_num = 1 + math.log(len(dataset), 10) / math.log(2, 10)
        max_v, min_v = max(dataset), min(dataset)
        return float(max_v - min_v) / category_num


def original_fomula(dataset):
        """教科書に載ってた方法
        参考: 入門統計解析法 日科技連出版社 p.9
        """
        m = 1 # 測定範囲
        n = len(dataset) # データ数
        max_v, min_v = max(dataset), min(dataset)
        R = max_v - min_v # 範囲
        h = math.sqrt(n) # 仮の区間数
        return int( (R/h) / m) * m # 区間の幅 測定範囲の整数倍に丸め込む


class statistics():

        def __init__(self, n=100, max=100):
                self.make_dataset(n, max)

        def make_dataset(self, n, max):
                """無作為にデータセットを作る"""
                self.dataset =  [random.random() * 100 for i in range(n)]

        def set_dataset(self, dataset):
                """指定のデータセットを用いる"""
                self.dataset = dataset

        def calculate(self):
                """統計量を算出する"""
                self._calculate_average()
                self._calculate_standard_deviation()
                self._calculate_T_score()

        def _calculate_average(self):
                """平均を算出する"""
                sum = 0.0
                for data in self.dataset:
                        sum += data
                self.average = sum / len(self.dataset)

        def _calculate_standard_deviation(self):
                """標準偏差を算出する"""
                sum = 0.0
                for data in self.dataset:
                        sum += pow( data - self.average, 2)
                self.standard_deviation = math.sqrt(sum / len(self.dataset))

        def _calculate_T_score(self):
                """偏差値を算出する"""
                self.T_scores = []
                for data in self.dataset:
                        standard_value = float(data - self.average) / self.standard_deviation
                        self.T_scores.append(standard_value * 10 + 50.0)

        def dump_data(self):
                """算出した統計量を表示する"""
                # display each score and each T-score
                fields = ["id", "score", "T-score"]
                pt = PrettyTable(fields=fields)
                for i, data in enumerate(self.dataset):
                        pt.add_row([i+1, data, "%2.2f" % self.T_scores[i]])
                pt.printt()

                # display average and standard deviation
                static_fields = ["average", "standard_deviation"]
                static_pt = PrettyTable(fields=static_fields)
                static_pt.add_row(["%2.2f" % self.average, "%2.2f" % self.standard_deviation])
                static_pt.printt()


        def histogram(self, func=struges_fomula):
                """ヒストグラムを作成する"""
                category_distance = func(self.dataset)
                sorted_dataset = sorted(self.dataset)
                category_start_v = min(self.dataset)
                category_end_v = category_start_v + category_distance
                category_freq = []
                _count = 0
                for val in sorted_dataset:
                        while True:
                                _key = "%s ~ %s" % tuple(map(int,[category_start_v, category_end_v]))
                                if category_start_v <= val < category_end_v:
                                        _count += 1
                                        break
                                else:
                                        category_start_v = category_end_v
                                        category_end_v += category_distance
                                        category_freq.append({"key":_key, "count":_count})
                                        _count = 0
                # display histogram
                fields = ["range", "visualization", "freqquency", "conposition ratio (%)"]

                pt = PrettyTable(fields=fields)
                pt.set_field_align("visualization" , "l")
                for info in category_freq:
                        pt.add_row([info["key"], "#" * info["count"], info["count"], "%3.2f" % (float(info["count"]) * 100 / len(self.dataset))])
                pt.add_row(["sum", "", len(self.dataset), "%3.2f" % 100])
                pt.printt()


if __name__ == "__main__":
        # 日本史 得点
        history_dataset = [73, 61, 14, 41, 49, 87, 69, 65, 36, 7, 53, 100, 57, 45, 56, 34, 37, 70]
        # 生物 得点
        biology_dataset = [59, 73, 47, 38, 63, 56, 15, 53, 80, 50, 41, 62, 44, 26, 91, 35, 53, 68]
        # p.78
        practice_dataset = [16.3, 22.4, 18.5, 18.7, 20.1]


        sta = statistics()

        # 基本統計量 算出 about 日本史
        print("----- 日本史 -----")
        sta.set_dataset(history_dataset)
        sta.calculate()
        sta.dump_data()
        sta.histogram()

        # 基本統計量 算出 about 生物
        print("----- 生物 -----")
        sta.set_dataset(biology_dataset)
        sta.calculate()
        sta.dump_data()
        sta.histogram()
        #sta.histogram(original_fomula)

        sta.set_dataset(practice_dataset)
        sta.calculate()
結果(一例)
----- 日本史 -----
+----+-------+---------+
| id | score | T-score |
+----+-------+---------+
| 1  |   73  |  58.79  |
| 2  |   61  |  53.52  |
| 3  |   14  |  32.85  |
| 4  |   41  |  44.72  |
| 5  |   49  |  48.24  |
| 6  |   87  |  64.95  |
| 7  |   69  |  57.04  |
| 8  |   65  |  55.28  |
| 9  |   36  |  42.53  |
| 10 |   7   |  29.77  |
| 11 |   53  |  50.00  |
| 12 |  100  |  70.67  |
| 13 |   57  |  51.76  |
| 14 |   45  |  46.48  |
| 15 |   56  |  51.32  |
| 16 |   34  |  41.65  |
| 17 |   37  |  42.96  |
| 18 |   70  |  57.47  |
+----+-------+---------+
+---------+--------------------+
| average | standard_deviation |
+---------+--------------------+
|  53.00  |       22.74        |
+---------+--------------------+
+---------+---------------+------------+-----------------------+
|  range  | visualization | freqquency | conposition ratio (%) |
+---------+---------------+------------+-----------------------+
|  7 ~ 24 | ##            |     2      |         11.11         |
| 24 ~ 42 | ####          |     4      |         22.22         |
| 42 ~ 60 | #####         |     5      |         27.78         |
| 60 ~ 78 | #####         |     5      |         27.78         |
| 78 ~ 96 | #             |     1      |          5.56         |
|   sum   |               |     18     |         100.00        |
+---------+---------------+------------+-----------------------+
...

5章 確立密度関数

ヒストグラムにおける階級の幅を極限まで狭めた曲式が確率密度関数であるという話で、
主要な確立分布を紹介していた。

  • 正規分布
  • 標準正規分布
  • カイ2乗分布
  • t分布
  • F分布

分布表の見方などを説明していた。x軸と確率密度関数の囲む面積は1である。とかの話。

2変数の関連を調べよう! (第6章)

数量データ と 数量データ
単相関係数
数量データ と カテゴリーデータ
相関比
カテゴリーデータ と カテゴリーデータ
クラメールの連関係数

負の相関, 無相関, 正の相関
○○以上なら相関があるといったことが言える統計学的な基準はない.

これらの指標はあくまでも、直線的な相関が見られる場合のみの判断ということに注意。

コードに起こす (理解度チェック)

また、少しコードを書いて理解を深めた。

  • 数量データ と 数量データ: 単相関係数
# numerable_data * numerable_data
dataset = [
        (3000, 7000),
        (5000, 8000),
        (12000, 25000),
        (2000, 5000),
        (7000, 12000),
        (15000, 30000),
        (5000, 10000),
        (6000,15000),
        (8000, 20000),
        (10000, 18000),]

def cal_cc(dataset):
        """calculate simple correlation coefficient"""
        N = len(dataset)
        sum_x, sum_y = 0.0, 0.0
        Sxx, Syy, Sxy = 0.0, 0.0, 0.0
        for x, y in dataset:
                sum_x += x
                sum_y += y
        ave_x = sum_x / N
        ave_y = sum_y / N

        for x, y in dataset:
                Sxx += pow(x - ave_x, 2)
                Syy += pow(y - ave_y, 2)
                Sxy += (x - ave_x) * (y - ave_y)

        return Sxy / math.sqrt(Sxx * Syy)
  • 数量データ と カテゴリーデータ: 相関比
# numerable_data * category_data
dataset2 = [
        (27, "Hermes"),
        (33, "CHANEL"),
        (16, "BURBERRY"),
        (29, "BURBERRY"),
        (32, "CHANEL"),
        (23, "Hermes"),
        (25, "CHANEL"),
        (28, "Hermes"),
        (22, "BURBERRY"),
        (18, "BURBERRY"),
        (26, "CHANEL"),
        (26, "Hermes"),
        (15, "BURBERRY"),
        (29, "CHANEL"),
        (26, "BURBERRY"),]

def cal_cr(dataset):
        """calculate correlation ratio"""
        N = len(dataset)

        label_count = {}
        sum_all = 0.0
        sum_each_label = {}
        sum_squre = {}
        # 級間変動 aboung class variation
        among_class_variation = 0.0
        # 級内変動 within class variation
        within_class_variation = 0.0

        for age, label in dataset:
                sum_each_label.setdefault(label, 0.0)
                label_count.setdefault(label, 0)
                sum_each_label[label] += age
                sum_all += age
                label_count[label] += 1

        for age, label in dataset:
                sum_squre.setdefault(label, 0.0)
                sum_squre[label] += pow(age - sum_each_label[label]/label_count[label], 2)

        within_class_variation = sum(sum_squre.values())

        for label in sum_each_label.keys():
                among_class_variation += label_count[label] * \
                                pow( (sum_each_label[label] / label_count[label]) - (sum_all / N), 2)
        return among_class_variation / (within_class_variation + among_class_variation)
  • カテゴリーデータ と カテゴリーデータ: クラメールの連関係数
# category_data * category_data
dataset3 = {
        "WOMAN":{"PHONE":34, "MAIL":61, "FACE":53},
        "MAN": {"PHONE":38, "MAIL":40, "FACE":74},}

# category_data * category_data
practice_set = {
        "Japanese":{"COFFEE":43, "TEA":33},
        "Westerm":{"COFFEE":51, "TEA":53},
        "Chinese":{"COFFEE":29, "TEA":41}}

def cal_cramer(dataset):
        """calculate Cramer's coefficient of association"""
        sum_each_sex = {}
        sum_each_method = {}
        row_num = len(dataset.keys())
        column_num = len(dataset[dataset.keys()[0]].keys())

        for sex in dataset.keys():
                sum_each_sex[sex] = sum(dataset[sex].values())
                for method in dataset[sex].keys():
                        sum_each_method.setdefault(method, 0.0)
                        sum_each_method[method] += dataset[sex][method]

        N = sum(sum_each_sex.values())

        # カイ二乗統計量
        chi_squared_statistic = 0.0

        for sex in dataset.keys():
                for method in dataset[sex].keys():
                        # 期待度数の計算
                        expected_frequency = float(sum_each_sex[sex] * sum_each_method[method]) / N
                        chi_squared_statistic += pow(dataset[sex][method] - expected_frequency, 2) / expected_frequency

        return math.sqrt(chi_squared_statistic / \
                                (N * ( min([row_num, column_num]) -1 )) )

独立性の検定をマスターしよう! (第7章)

「検定」とは母集団について分析者が立てた仮説が正しいかどうかを標本のデータから推測する分析手法のこと
正確には「統計的仮説検定」という。

独立性の検定, 相関比の検定, 無相関の検定, 母平均の差の検定, 母比率の差の検定など検定の種類はいろいろとあるが、分析手順は基本的に同じ。

  1. 母集団の定義する.
  2. 帰無仮説と対立仮説を立てる.
  3. どの「検定」を行うかを選択する.
  4. 有意水準を決定する.
  5. 標本のデータから検定統計量の値を求める.
  6. 求めた検定統計量の値が棄却域に入っているかどうかを調べる.
  7. 棄却域に入っていたならば、「対立仮説は正しいとの結論」そうでなければ、「帰無仮説は誤っているとはいえない」との結論を下す.
疑問

理解はできたつもりだが、実務で「検定」をどのようにして利用すればいいのかが思いつかなかった。
無相関であることを示したり、この因子は目的変数に効果を与えていないことを示すなど例だけでは、現実世界のどこで使われているのかがイメージできなかった。
また、教科書などで検定をする際に統計量を求めるときには、予め、「自由度2のカイ2乗分布に従う」、「平均が○○で標準偏差が××の正規分布とみなして」...といった前提が与えられているから、問題は解けるが、実際のデータでは、どの分布に従っているかなんてわからんし、そこをどうすればいいのかな~なんてモヤモヤしていた。

まだまだ勉強なり続けます。近いうちに回帰分析編の内容もまとめられたらいいな。