Pythonを用いて回帰分析をする。

機械学習のエッセンス」*1を勉強しているので、拙いながらもメモのため今回のように記事を残しておく。今回は、各都道府県の離婚率や収入、人口のデータ(cv.csv)を用いて前処理から線形回帰を行うことを目的とする。記載している回帰分析のコードは、機械学習のエッセンスを参考にして書いた。

まずは最初に以下のコードをimportして、準備をしておく。

import numpy as np
import matplotlib.pyplot as plt
import csv
import os
import pandas as pd
import linearreg
from mpl_toolkits.mplot3d import axes3d
%matplotlib inline

osモジュールは、作業ディレクトリの確認と変更のためにimportしておく。csvをimportする際に、必要に応じて作業ディレクトリを変更するので、その時に使えば良い。pandasは、後々データ整形がしやすくなるのでやっとく。linearregモジュールは、線形回帰のモジュールで解説は以下で行う。なお、linearregモジュールは、に記載のコードを参考にしている。

まず、データをpandasを使って整形し、離婚率divorceと収入incomeの項目以外を削除したものをdataとする。

csv_file = pd.read_csv("cv.csv")
data = csv_file.drop(["pop","car", "nm_m_30",
                                    "nm_f_30" ,"w_pop" , "conv", "gsm", "sm", "gpp",
                                   "Unnamed: 0"], axis=1)

次に、dataのincomeおよびdivorceを散布図として表す。

plt.scatter(data[data_col[0]],data[data_col[1]],marker="o")

plt.xlabel(data[data_col[0]].name)
plt.ylabel(data[data_col[1]].name)

plt.show()

出力結果が以下の画像。

f:id:bislogyaruka:20190126160337p:plain

次に、回帰分析のパラメータを返す関数reg1dimを以下のように定義する。

def reg1dim(x, y):
    n = len(data)
    a = ((np.dot(x, y) - y.sum() * x.sum() / n ) /
         ((x**2).sum() - x.sum()**2 / n ))
    b = (y.sum() - a * x.sum()) / n
    return a, b

パラメータa, b が上のようなコードで表されるのは、直線による回帰分析の一般論から導出されるのでここでは省略する。xをincome、yをdivorceとする。

x = data["income"]
y = data["divorce"]

x_income = np.array(x)
y_divorce = np.array(y)

a,bを求めるため、回帰分析を行う。

a, b = reg1dim(x_income, y_divorce)

最後に、求めたa,bを用いて回帰直線と散布図を描き、可視化する。scatterするのに、配列に変換する必要があるので、xをx_incomeのようにarray配列に置き換えた。また、x軸とy軸にラベルを貼ったほうがグラフが見やすくなるので、plt.xlabel(data[data_col[0]])でラベルの名前をつける。また、プロット図自体の範囲をどこからどこまでにするか決める必要がある。今回は、x軸を[0, xmax]、y軸を[b, a * xmax + b]の範囲で描いている。本来は0から始めないと行けないかもしれない。その場合は、[0, a*xmax +b]のように定義域を変えれば良い。

plt.scatter(x_income, y_divorce, color="k")
plt.xlabel(data[data_col[0]].name)
plt.ylabel(data[data_col[1]].name)
xmax = x_income.max()
plt.plot([0, xmax], [b, a * xmax + b], color="k")
plt.show()

出力結果が以下の図になる。出力結果を見ると、なんとなく収入が上がれば離婚率が下がる傾向にあることがわかるが、この分析からは因果関係までは読み取れない。因果推論については、現時点で私にはわからないので、パスさせてもらう。

f:id:bislogyaruka:20190126161707p:plain

ちなみに、aとbの値は以下になる。

f:id:bislogyaruka:20190126161959p:plain
aとbの値

・感想
今回は簡単な回帰分析をする際のコードの書き方を記事にしたが、思いの外データの前処理に時間がかかった。具体的には、pandasをimportして離婚率と収入の項目のみ抽出する部分と配列に直すのに時間がかかった。Python機械学習のコードを書きつつ、同時並行的に入門的なPythonの参考書で勉強する方がコードを書くスキルを効率的に身につけられると思った。わからないことがあれば、参考書の最初に戻ってみる、ググるなどをこなしていく。