Create a gist now

Instantly share code, notes, and snippets.

PyMC3でDNNをADVIで学習
# -*- coding: utf-8 -*-
import pymc3 as pm
from pymc3 import Model, Normal
import pandas as pd
import theano.tensor as T
import theano
import sklearn
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
sns.set_style('white')
##データの前処理
#ファイルの読み込み
train_df = pd.read_csv("outlier_removed_high_feature_train.csv", index_col=0, header=0)
test_df = pd.read_csv("high_feature_test.csv", header=0, index_col=0, names=['is_Higgs','m_jj', 'm_jjj', 'm_lv', 'm_jlv', 'm_bb', 'm_wbb', 'm_wwbb'])
train_df = train_df.reindex(np.random.permutation(train_df.index)).reset_index(drop=True)
#トレーニングデータの目的変数をまず-1を0に変換してから抽出
train_df.loc[(train_df.is_Higgs == -1 ),"is_Higgs"] = 0
Y_train = train_df.values[0::,0]
#トレーニングデータの説明変数の正規化
train_df_means = train_df.mean().values
train_df_max = train_df.max().values
train_df_min = train_df.min().values
X_train = train_df.apply(lambda x: (x-x.mean())/(x.max()-x.min()), axis=0).fillna(0).values[0::,1::]
#テストデータのの目的変数もまず-1を0に変換してから抽出
test_df.loc[(test_df.is_Higgs == -1 ),"is_Higgs"] = 0
Y_test = test_df.values[0::,0]
#テストデータの説明変数を上の正規化定数を使って正規化
X_test = test_df.values[0::,1::]
for i in range(X_test.shape[1]):
X_test[0::,i] = (X_test[0::,i] - train_df_means[i+1])/(train_df_max[i+1]-train_df_min[i+1] )
total_size = len(Y_train) #トレーニングデータの総数
n_feature = X_train.shape[1] #特徴量の数
n_first_layer_hidden = 20 #第一隠れ層のニューロンの数
n_second_layer_hidden = 15 #第二隠れ層のニューロンの数
mini_batch_size = 10000 #ミニバッチサイズ
#ミニバッチの作成
def create_minibatch(data, mini_batch_size=50,seed=0):
rng = np.random.RandomState(seed)
while True:
ixs = rng.randint(len(data), size=mini_batch_size)
yield data[ixs]
minibatches = zip(
create_minibatch(X_train, mini_batch_size),
create_minibatch(Y_train, mini_batch_size),
)
##sharedテンソルの作成
ann_input = theano.shared(X_train)
ann_output = theano.shared(Y_train)
with pm.Model() as neural_network:
#事前分布 on weights (標準正規分布で正則化)
weights_in_1 = Normal('w_in_1', mu=0, sd=1, shape=(n_feature, n_first_layer_hidden) )
weights_1_2 = Normal('w_1_2', mu=0, sd=1, shape=(n_first_layer_hidden, n_second_layer_hidden))
weights_2_out = Normal('w_2_out', mu=0, sd=1, shape=(n_second_layer_hidden,))
#事前分布 on bias (標準正規分布で正則化)
bias_on_1 = Normal("b_on_1", mu=0, sd=1,shape=(n_first_layer_hidden,))
bias_on_2 = Normal("b_on_2", mu=0, sd=1,shape=(n_second_layer_hidden,))
# Build neural-network using tanh activation function
act_1 = T.tanh(T.dot(ann_input, weights_in_1) + bias_on_1)
act_2 = T.tanh(T.dot(act_1, weights_1_2) + bias_on_2)
act_out = T.nnet.sigmoid(T.dot(act_2, weights_2_out))
# Binary classification -> Bernoulli likelihood
out = pm.Bernoulli('out', act_out, observed=ann_output)
# Tensors and RV that will be using mini-batches
minibatch_tensors = [ann_input, ann_output]
minibatch_RVs = [out]
#minibatch_ADVIで推定
v_params = pm.variational.advi_minibatch(n=50000,
minibatch_tensors=minibatch_tensors,
minibatch_RVs=minibatch_RVs,
minibatches=minibatches,
total_size=total_size,
learning_rate=1e-2, epsilon=1.0)
trace = pm.variational.sample_vp(v_params, draws=5000)
# Replace shared variables with testing set
ann_input.set_value(X_test)
ann_output.set_value(Y_test)
#(パラメータの)サンプルを事後分布から生成
ppc = pm.sample_ppc(trace, model=neural_network, samples=50)
#アウトプットのサンプルを取り出し、サンプルの平均を取る。
#そして0.5より大きければTrueを取るテストデータサイズのベクトルをpredに格納。
pred = ppc['out'].mean(axis=0) > 0.5
print('Accuracy = {}%'.format((Y_test == pred).mean() * 100))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment