AOKI's copy&paste archive

高専から駅弁大学から東工大を経て大企業へ 浅く広い趣味とかキャリアの日記を

Statistics Test As Data Science

# -*- coding: utf-8 -*-
"""
Spyder Editor

This is a temporary script file.
"""
#import modules
import numpy as np
import pandas as pd
import datetime
import scipy.stats as stats
import matplotlib.pyplot as plt
from sklearn.tree import DecisionTreeClassifier
from sklearn.tree import export_graphviz
from sklearn.cluster import KMeans
from io import StringIO
import pydotplus
from IPython import display

#input master data from csv
df = pd.read_csv(r"C:\Users\*.csv", encoding = "utf-8")

#%%data cleaning
df['生年月日'] = pd.to_datetime(df['生年月日'])
df['提出日'] = pd.to_datetime(df['提出日'])

#%%
#2020/1/1からの日付差分で再格納
df['生年月日'] = (df['生年月日'] - datetime.datetime(2020, 1, 1)).dt.days
df['提出日'] = (df['提出日'] - datetime.datetime(2020, 1, 1)).dt.days

#0でないと結果がまとまらない
df = df.fillna(0)#???

#%%string to dummy
df_pref = pd.get_dummies(df['都道府県'])
df_business = pd.get_dummies(df['業種'])

#objective variables to dummy
df['Memo'] = np.where(df['提出日'] == 0, 0, 1)

#concat
df = pd.concat([df, df_pref, df_business], axis = 1)

#%%chi2(squre) test
chi2test = pd.DataFrame()
for column_name in df_pref:
    if column_name != 0:#avoid error
        stat_master = df[df['Memo'] == 1]
        pref = stat_master[stat_master['都道府県'] == column_name]['Miss']
        if len(stat_master) > 0:
            crossed = pd.crosstab(stat_master['Miss'], stat_master[column_name])
            x2, p, dof, expected = stats.chi2_contingency(crossed, correction = False)
            n = crossed.sum().sum()
            cramerV = np.sqrt(x2 / (n * (np.min(stat_master.shape) - 1)))
            chi2test = chi2test.append([['*', column_name, len(pref), sum(pref), 
            sum(pref) / len(pref) * 100, p, x2, cramerV]])


chi2test.columns = ['種類', '都道府県', '全件数', 'Miss', '*率(%)', 'P値', 'カイ自乗値', 'クラメールの連関係数']
chi2test.to_csv('chi2test_pref.csv')


#%%T test
ttest = pd.DataFrame()
for column_name in df_pref:
     stat_master = df[df['*'] == 1]
     construct = stat_master[stat_master['都道府県'] == column_name]['Miss']
     if len(construct) > 0:
         others = stat_master[stat_master['都道府県'] != column_name]['Miss']
         t, p = stats.ttest_ind(construct, others, alternative='less')
         ttest = ttest.append([['*',column_name, len(pref), sum(pref),
                                sum(pref) / len(pref) * 100, t, p]])
    
ttest.columns = ['*', '都道府県', '全件数', 'Miss', '割合', 'T', 'P']
ttest.to_csv('ttest_pref.csv')

#%%dummy source data drop
df = df.drop(columns=['ID', '都道府県', '業種', '提出日'])

#%%cal corr
df_corr = df.corr()
df_corr.to_csv('correlation.csv')

#%%
def decisiontree(x, y):
    if len(x) > 10000:
        depth = 4
    else:
        depth = 3
    tree_clf = DecisionTreeClassifier(max_depth = depth)
    tree_clf.fit(x, y)
    predicted = tree_clf.predict(x)
    print(len(y), ': ', int(sum(predicted == y) / len(y) * 100), '%')
    
    g = StringIO()
    export_graphviz(tree_clf, out_file = g, 
        feature_names = x.columns, 
        class_names = x.columns, 
        rounded=True, filled=True)
    graph = pydotplus.graph_from_dot_data(g.getvalue())
    graph.set_fontname('MS UI Gothic')
    #Setting font for Node 
    for node in graph.get_nodes():
        node.set_fontname('MS UI Gothic')
    
    #Setting font for Edges 
    for e in graph.get_edges((sad)
        e.set_fontname('MS UI Gothic')
        
    graph.progs = {'dot': u"C:\\Program Files\\Graphviz\\bin\\dot.exe"}
    display.display(display.Image(graph.create_png()))
    

 for column_name in df_pref:
     construct_decisiontree = df[df['Memo'] == 1]
     construct_decisiontree = construct_decisiontree[construct_decisiontree['都道府県'] == column_name]
     x = construct_decisiontree['Miss']
     construct_decisiontree = construct_decisiontree.drop(columns = ['Miss', 'Memo'])
     if len(x) > 1000:
         print(column_name)
         decisiontree(construct_decisiontree, x)

kmeans = KMeans(n_clusters=3, random_state=0)
clusters = kmeans.fit(df)
df['cluster'] = clusters.labels_
print(df['cluster'].value_counts())

df.to_csv('Reformed.csv')

print('fin')

最近の作業内容を内容がわからないよう多くを削除して、残った部分も一部書き換えて公開する。
コードならGitHub使えって話だし、エンジニア的キャリアを見通してポートフォリオを作るにはそのほうがいいのだが惰性。

コードの内容、目的としては、統計検定や決定木分析。
流行りっぽい言い方をすれば、AI分析だが、その実、大したことはないし、内部的な数学は大学初学年レベルといったところか。

私自身、統計は実はそれほど得意でないこともあり、試行錯誤による手戻りも多かった。
まずコード中にT検定(Ttest)があるが、これは定量データ同士でしかできない。
在学中は揃ったデータセットに恵まれたが、世の中それほど便利でもない。

統計学入門−第5章


そこで厳密な統計検定には、χ二乗検定が必要で、それがChi2testだ。
マーケのABテストの効果の評価も、同様の理論が用いられる。
私の分析はToC分析ではないが、ベースとなる数学基礎は共有できるし、経験としてマーケへの転換にも活かせるかもしれない。

「ノイズ」を避けるために、マーケターが持つべき「統計的思考」と「判断の構造化」 | AdverTimes(アドタイ) by 宣伝会議

ちなみにここでの苦労として、T検定は条件で片側ずつ検定できたが、χ二乗検定は平均との差で左右、すなわち増減が加味できない。
平均と差があるか、ないかしか、判別できない。
実際にはどっちに有意差があるかが、人間としては重要だ。
そこの判別が難しかったが、幸いにしてT検定とP値は類似傾向だったので、どっちで出たかで分類した。
あるいは平均との差の正負で。

しかしガチガチの統計検定は、有意差の有無の仮説棄却しかできない。
応用性に難がある。
そのためFor文でループして、上記の例文では都道府県で回している。
ちなみにこの前処理として、名義尺度でテキストになっているデータを、47都道府県でそれぞれダミー変数変換している。

しかしこれでもP値でランク付けするのが限界であるし、棄却域に対する真偽値以外の意味は薄い。
そこで尤度、信頼性は相対的に低下するが、決定木でデータ構造の分析も試みた。
お硬い統計では、結果自体の外部(経済)性が期待できないというのもある。
世の中、AIを評価しても、その中を読める人は限られる。
事実上、このレベルですらブラックボックスと言える。
将来的にはPytorchも自習がてら活用したいが。。
実際には上記のダミー変数などなどでカラム数は、かなりの量であるため、記載の木の深度がいい塩梅だった。

ただし上記の分析のそもそものデータソースに関して、Corrに示した相関表、相関係数自体が低い。
そもそもデータにおいて、因果関係があるかもやや曖昧ではある。。

将来的にここで”優位”だった変数をまとめて、ロジットモデルを構成して、予想式を立式したいところだ。
ちなみにこのままロジットモデルにデータフレームの元データを入れると、解が収束しない。。

これもダミー変数を主とするデータ構造自体に課題がある。
世の中、計量尺度を集積するのは非常に難しいし、名義尺度ですら正しくキレイに集積するのも難しい。
表記ゆれなど。あるいはValidationがないと。。

こう全体的な部分を半ば自前で経験すると、統計局の仕事に非常に頭が下がる。
かつてあった捏造は許せないが、その苦労もかなりのものだったろう。