Pythonを使って変数選択！

機械学習はデータが命です。データが精度を左右するので、精度を上げるためにデータを増やし、変数をどんどん追加してくという方向になりがちです。しかし、変数の数を多くすると、計算時間の増加をまねいたり、特定のクラスの一部のデータの影響で過学習したりなどの問題が起こります。

意味のある変数だけを抽出できたり、次元を減らすように要約できたりすれば、重要な要因がわかりますし、計算時間も減らせます。見たい番組が多すぎて色々ザッピングした結果、何も記憶に残っていないみたいなことがなくなります。今回は、このような変数の削減方法において見ていきます。

先に実装結果を示すと、各手法によって選択される変数が異なるため、同一の機械学習アルゴリズムで同一パラメータにおいてもテストデータの正答率が異なっています。今回は変数増加法の正答率が一番高く、もとの変数の1/2以下になっています。

今回取り上げる手法

• 変数選択
  単変量選択、変数減少法、変数増加法、ロジスティック回帰のL１正則化による選択、ランダムフォレストの変数重要度による選択

• 次元削減
  主成分分析

変数選択と次元削減

変数を減らすためには２種類の方法があります。それは、変数選択と次元削減です。変数選択は全変数群から変数を抜き取る方法であり、各変数の値は変わりません。それに対して次元削減は、次元を集約して新たな次元を作ることで変数を減らします。例えば、何十人の臣下を抱える王様だと、良いことを言う臣下を数人選択するのが変数選択であり、臣下の言っていることをいくつかのテーマでまとめるのが次元削減です。

実装のための準備

データの前処理

今回はUCI Machin Learning Repositoryから『Bank Marketing Data Set』のデータを拝借し、分類問題において変数削減によって正答率の向上を目指します。このデータは、ポルトガルの銀行で電話によるダイレクトマーケティングした際に、デポジットの契約をしたか否かと、その時の情報を16個の特徴量で表しています。今回使用したデータはここからダウンロードできます。リンク先のbank.zipの中に「bank.csv」、「bank-full.csv」があります。bank.csvとbank-full.csvは変数は同じですが、レコード数が異なります。今回はbank-full.csvのデータを使います。

図を見てわかる通り、特徴量の中には9個のカテゴリ変数があり、このカテゴリ変数を入力データとして使うにはダミー変数に変換しなければいけません。計算機での計算は数字で行うため、急によくわからない文字列が与えられても計算できないからです。今回のデータでいえば、2つ目の変数の「job」には9種類のカテゴリがあるため、以下のように9つのダミー変数に替えます。

またカテゴリ変数以外の変数の値は、変数によって上限がまちまちであり、変数間の比較ができません。人間とキリンの平均身長を比較しても意識高い系以外は興味もないでしょう。上限の値が大きい変数が少し違うだけで予測誤差が大きくなるので、予測誤差に対してその変数の影響が強くなります。そこで各変数を平均0、標準偏差１に正規化し変数間で比較できるようにします。
被説明変数は、契約したか否かの2択であり、これは0－1に変換します。以上のことをまとめると次のようなコードになります。

from pathlib import Path
import pandas as pd
from sklearn.preprocessing import StandardScaler
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.svm import SVC


output_dir_path = Path('./output')
if not output_dir_path.exists():
    output_dir_path.mkdir(parents=True)

# データを読み込んで、入力データと出力データに分割
datas = pd.read_csv('bank-full.csv', sep=';', encoding='utf-8')
x_datas = datas.drop('y', axis=1)
y_datas = pd.DataFrame(datas['y'], columns=['y'])

# ダミー変数に変換するカテゴリ変数の指定
categorical_variable_names = ['job', 'marital', 'education',
                               'default', 'housing', 'loan', 'contact',
                               'month','poutcome']
# カテゴリ変数に一括変換
x_dummy = pd.get_dummies(x_datas[categorical_variable_names])

x_datas.drop(categorical_variable_names, axis=1, inplace=True)

# それ以外の変数は正規化
x_datas_std = StandardScaler().fit_transform(x_datas)
x_datas_std = pd.DataFrame(x_datas_std, columns=x_datas.columns)

# 入力データの完成
x_datas = pd.concat([x_datas_std, x_dummy], axis=1)

# 出力データは0-1に変換し、１次元に整形
class_mapping = {label:idx for idx, label in enumerate(np.unique(y_datas))}
y_datas['y'] = y_datas['y'].map(class_mapping)
y_datas = y_datas.values.reshape(-1)

# データのうち３割をテストデータに分割
n_train_rate = 0.7
train_x_data, test_x_data, train_y_data, test_y_data = train_test_split(x_datas,
                                                                        y_datas,
                                                                        test_size=1-n_train_rate)

変数選択

単変量選択　単純なランキングこそ一番

単変量選択とは、ある指標のもとに変数を順位付けし、その順位に従って指定した数だけの変数を取得する方法です。オリコンランキングのCD売上げTOP10などです。順位ができてしまえば、上位から順番に選択していくだけです。
下の例では、「SelectKBest」を使っています。第一引数に順位をつける関数を入れ、引数kで選択する変数の数を指定します。分類と回帰では、第一引数に入れる関数が異なり、デフォルトでは分類用の「f_classcif」が入っています。これはF検定によって順位をつけています。分類、回帰それぞれで使える関数は以下の通りです。

• 分類：「chi2」、「f_classif」、「mutual_info_classif」

• 回帰：「f_regression」、「mutual_info_regression」

from sklearn.feature_selection import SelectKBest

col_names = train_x_data.columns

# 今回はSVMで学習（どの変数選択手法でも同じ）
model = SVC(kernel='rbf', C=1, gamma=0.1)

selected_col_names = []
scores = []
max_score = 0
selected_k = 0

# 変数の数
ks = np.arange(1, len(col_names)+1)
for k in ks:

    # ｋ個の変数を選択する(評価関数はデフォルトのf_classif)
    skb = SelectKBest(k=k)
    skb.fit(train_x_data, train_y_data)

    # skbに基づいて入力データからK個変数を選択
    train_x_data_ = skb.transform(train_x_data)
    test_x_data_ = skb.transform(test_x_data)

    model.fit(train_x_data_, train_y_data)

    # skb.get_supprt()により選択された変数のインデックスを取得
    col_names_ = np.array(col_names)[skb.get_support()]
    score = model.score(test_x_data_, test_y_data)

    if score > max_score:
        selected_k = k
        max_score = score

    scores.append(score)
    selected_col_names.append(','.join(col_names_))

    print('k {:3} : score {:5.2f}%'.format(k, score*100))

print('')
print('[MAX SCORE] : k {} : score {}%'.format(selected_k, int(max_score*100)))
print(selected_col_names[int(np.argmax(scores))])


ks = [len(cols.split(',')) for cols in selected_col_names]
results = pd.DataFrame([ks,
                        scores,
                        selected_col_names]).T
results.columns = ['n_variable', 'score', 'selected_col_name']

results.to_csv(output_dir_path.joinpath('ufs.tsv'),
               sep='\t',
               encoding='utf-8',
               index=False)

変数減少法　怪しいやつから消していこう

評価値においてモデルが改善する変数を１つずつ除去するという方法を、入力データの全変数から改善が終わるまで繰り返していきます。映画『インシテミル』(*1)のようですね。

col_names = train_x_data.columns.tolist()

model = SVC(kernel='rbf', C=1, gamma=0.1)

# 列名を減少させていくことで、入力データをヘンス指定で取得する
target_col_names = col_names.copy()
target_col_names = target_col_names[::-1]

base_score = -1       # 前ループの正答率
loop_max_score = 0    # 現ループの最大正答率
del_col_name = ''     # 削除変数

selected_col_names = []
scores = []

selected_k = 0
selected_score = 0

# 現ループの最大正答率が前ループより大きいかつ、選択できる変数があるならば続行
while (base_score < loop_max_score) and len(target_col_names) > 0:

    base_score = loop_max_score

    for target_col_name in target_col_names:
        col_names_ = target_col_names.copy()
        col_names_.remove(target_col_name)

        # 入力データから変数の選択
        train_x_data_ = train_x_data.loc[:, col_names_]
        test_x_data_ = test_x_data.loc[:, col_names_]

        model.fit(train_x_data_, train_y_data)
        score = model.score(test_x_data_, test_y_data)

        if score > loop_max_score:
            loop_max_score = score
            del_col_name = target_col_name

    if loop_max_score <= base_score:
        break

    if del_col_name in target_col_names:
        # 正答率を最も改善する変数を削除
        target_col_names.remove(del_col_name)
        selected_col_names.append(','.join(target_col_names))
        scores.append(loop_max_score)
        print('k {:3} : score {:5.2f}%'.format(len(target_col_names), loop_max_score * 100))

print('')
print('[MAX SCORE] : k {} : score {}%'.format(len(selected_col_names[int(np.argmax(scores))].split(',')),
                                              loop_max_score))
print(selected_col_names[int(np.argmax(scores))])


ks = [len(cols.split(',')) for cols in selected_col_names]
results = pd.DataFrame([ks,
                        scores,
                        selected_col_names]).T
results.columns = ['n_variable', 'score', 'selected_col_name']

results.to_csv(output_dir_path.joinpath('backward_elimination.tsv'),
               sep='\t',
               encoding='utf-8',
               index=False)

変数増加法　信頼できるやつで固めてしまおう

変数増加法は変数減少法と名前が似ている通り、変数減少法とは反対の方法です。この方法では、モデルを改善する変数を１つずつ追加していくのを改善が終わるまで繰り返します。

col_names = train_x_data.columns.tolist()

model = SVC(kernel='rbf', C=1, gamma=0.1)

target_col_names = []              # 選択変数の格納
loop_col_names = col_names.copy()  # 未選択変数を格納

base_score = -1
loop_max_score = 0
additional_col_name = ''

selected_col_names = []
scores = []

selected_k = 0
selected_score = 0

# 現ループが前ループより正答率が高くかつ、未選択変数が残されているならば続行
while (base_score < loop_max_score) and len(loop_col_names) > 0:

    base_score = loop_max_score

    for target_col_name in loop_col_names:
        col_names_ = target_col_names.copy()
        col_names_.append(target_col_name)

        # 入力データから変数の選択
        train_x_data_ = train_x_data.loc[:, col_names_]
        test_x_data_ = test_x_data.loc[:, col_names_]

        model.fit(train_x_data_, train_y_data)
        score = model.score(test_x_data_, test_y_data)

        if score > loop_max_score:
            loop_max_score = score
            additional_col_name = target_col_name

    if base_score >= loop_max_score:
        break

    # 未選択変数群から選択分数を削除
    loop_col_names.remove(additional_col_name)
    target_col_names.append(additional_col_name)
    selected_col_names.append(','.join(target_col_names))
    scores.append(loop_max_score)
    print('k {:3} : score {:5.2f}%'.format(len(target_col_names), loop_max_score * 100))

print('')
print('[MAX SCORE] : k {} : score {}%'.format(len(selected_col_names[int(np.argmax(scores))].split(',')),
                                              loop_max_score))
print(selected_col_names[int(np.argmax(scores))])


ks = [len(cols.split(',')) for cols in selected_col_names]
results = pd.DataFrame([ks,
                        scores,
                        selected_col_names]).T
results.columns = ['n_variable', 'score', 'selected_col_name']

results.to_csv(output_dir_path.joinpath('forward_selection.tsv'),
               sep='\t',
               encoding='utf-8',
               index=False)

ロジスティック回帰のL1正則化

機械学習において過学習を防ぐために、目的関数にペナルティ項を加えることがよくあります。このペナルティ項として、L1ノルムと呼ぼれる重みパラメータの絶対値の合計を使うのがL1正則化です。L１正則化を使うことで、スパースなパラメータを得やすくなります。パラメータが閾値よりも低いものは、学習器に与える影響が小さいということで、削ることができます。詳しいことはこちらをRでL1 / L2正則化を実践する - 渋谷駅前で働くデータサイエンティストのブログ。L1、L2正則化に関して図示されなが述べられているので、非常に参考になります。

from sklearn.svm import LinearSVC
from sklearn.feature_selection import SelectFromModel

col_names = train_x_data.columns.values

model = SVC(kernel='rbf', C=1, gamma=0.1)

selected_col_names = []
scores = []
selected_k = 0
max_score = 0

# ロジスティック回帰のハイパーパラメータ
Cs = [10**i for i in range(-3, 3, 1)]

for C in Cs:

    # L1正則化でロジスティック回帰を学習します
    lsvc = LinearSVC(C=C, penalty='l1', dual=False).fit(train_x_data, train_y_data)
    slm = SelectFromModel(lsvc, prefit=True)

    # 選択変数の抽出
    train_x_data_ = train_x_data.loc[:, slm.get_support()]
    test_x_data_ = test_x_data.loc[:, slm.get_support()]

    model.fit(train_x_data_, train_y_data)
    score = model.score(test_x_data_, test_y_data)

    col_names_ = col_names[slm.get_support()]
    if score > max_score:
        max_score = score
        selected_k = len(col_names_)

    selected_col_names.append(','.join(col_names_))
    scores.append(score)

    print('C {:5} : k {:3} : score {:5.2f}%'.format(C, len(col_names_), max_score * 100))

print('')
print('[MAX SCORE] : k {} : score {}%'.format(selected_k, max_score))
print(selected_col_names[int(np.argmax(scores))])

ks = [len(cols.split(',')) for cols in selected_col_names]
results = pd.DataFrame([Cs,
                        ks,
                        scores,
                        selected_col_names]).T
results.columns = ['C', 'n_variable', 'score', 'selected_col_name']

results.to_csv(output_dir_path.joinpath('lsvc.tsv'),
               sep='\t',
               encoding='utf-8',
               index=False)

ランダムフォレストの変数重要度による変数選択　インフルエンサに気をつけろ

詳細に述べると非常に長くなりそうなので、概念的な部分だけを抑えたいと思います。詳細はRandom Forestで計算できる特徴量の重要度 - なにメモで説明されており、今回調べる上で非常に勉強になりました。ランダムフォレストでは、複数の決定木を作成し、その決定木の多数決などで分類しますが、このとき学習された複数の決定木において、変数の値をぐちゃぐちゃにしたときにどのくらい精度が悪くなるかで変数の重要度を測定しています。ぐちゃぐちゃしたときに精度が落ちる変数ほど、推測するうえでキーになっているということです。日頃でたらめなことを言っている人が何を言っても影響は少ないですが、誰もが信頼している人が嘘をついたときの影響は大きいといった感じでしょうか。

from sklearn.ensemble import RandomForestClassifier

col_names = train_x_data.columns.values

model = SVC(kernel='rbf', C=1, gamma=0.1)

selected_col_names = []
scores = []
selected_k = 0
max_score = 0

rf = RandomForestClassifier()
rf.fit(train_x_data, train_y_data)

# ランダムフォレストによる変数重要度順に並び替え
col_names_ = col_names[np.argsort(rf.feature_importances_)[::-1]]

for n_variable in range(1, len(col_names)+1):

    # 指定した変数分を重要度が高いものから選択
    train_x_data_ = train_x_data.loc[:, col_names_[:n_variable]]
    test_x_data_ = test_x_data.loc[:, col_names_[:n_variable]]

    model.fit(train_x_data_, train_y_data)
    score = model.score(test_x_data_, test_y_data)

    if score > max_score:
        max_score = score
        selected_k = n_variable

    selected_col_names.append(','.join(col_names_[:n_variable]))
    scores.append(score)

    print('k {:3} : score {:5.2f}%'.format(n_variable, max_score * 100))

print('')
print('[MAX SCORE] : k {} : score {}%'.format(selected_k, max_score))
print(selected_col_names[int(np.argmax(scores))])

ks = [len(cols.split(',')) for cols in selected_col_names]
results = pd.DataFrame([ks,
                        scores,
                        selected_col_names]).T
results.columns = ['n_variable', 'score', 'selected_col_name']

results.to_csv(output_dir_path.joinpath('randomforest_importance.tsv'),
               sep='\t',
               encoding='utf-8',
               index=False)

次元削減

主成分分析　ギャップ萌えこそ良い軸

言わずわずとしれた主成分分析です。分散が大きいものほど情報を持っているということで、分散が大きくなる軸を作り、その軸を使ってデータを変換します。その際に、何番目までの軸を使うかで次元を決めていきます。詳しくはこちらを主成分分析の考え方 | Logics of Blue

from sklearn.decomposition import PCA

model = SVC(kernel='rbf', C=1, gamma=0.1)

selected_col_names = []
scores = []
selected_k = 0
max_score = 0

for n_variable in range(1, train_x_data.shape[1]+1):

    pca = PCA(n_components=n_variable)
    pca.fit(x_datas)
    train_x_data_ = pca.transform(train_x_data)
    test_x_data_ = pca.transform(test_x_data)

    model.fit(train_x_data_, train_y_data)
    score = model.score(test_x_data_, test_y_data)

    if score > max_score:
        max_score = score
        selected_k = n_variable

    scores.append(score)
    print('k {:3} : score {:5.2f}%'.format(n_variable, max_score * 100))

print('')
print('[MAX SCORE] : k {} : score {}%'.format(selected_k, max_score))

ks = [i for i in range(1, train_x_data.shape[1]+1)]
results = pd.DataFrame([ks,
                        scores,
                        selected_col_names]).T
results.columns = ['n_variable', 'score', 'selected_col_name']

results.to_csv(output_dir_path.joinpath('pca.tsv'),
               sep='\t',
               encoding='utf-8',
               index=False)

一章丸々次元削減に割り当てられており、実装コードも載せられているので、イメージだけでなく実感として理解をすることができます。また、データの前処理に関しても詳しく書かれており、機械学習をする上で大事なデータの準備も併せて学べることができます。

［第2版］Python機械学習プログラミング 達人データサイエンティストによる理論と実践 impress top gearシリーズ

参考資料

(*1)インシテミル-wikipedia