実装関連事項

各種プログラミング言語の基本的な書き方やソフトウェア等の使用方法について.

ロジスティック回帰 (logistic regression) は非線形多重回帰法のひとつ.回帰とついているが,普通,分類のために用いられることが多い.サポートベクターマシンをはじめとするアドホックな解析をしなければ変数の重要度を出せないような手法,すなわち,ブラックボックス的な学習法とは異なり,計算の過程において重要な特徴量が自然に明らかにされるホワイトボックス的な手法である.

学習させるデータの構造

以下のようなデータセットを考える.このデータセットでは,69次元 (アトリビュート) のインプットベクトルに対して0または1のスカラーの値が対応している.インプットベクトルの各アトリビュートはコンマ (,) によって分割されており,インプットベクトルとターゲットベクトル (スカラー) はタブ (\t) によって分割されている.

64,-1,50,6,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,3,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0	0
60,-1,50,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0	1
66,-1,100,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0	0
64,1,25,12,0,0,0,0,0,0,0,1,0,0,0,0,0,0,3,1,0,2,1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1	0
89,-1,25,6,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,0,0,1,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0	0
76,-1,25,2,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0	0
40,1,25,4,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0	0
66,-1,50,4,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0	1
77,-1,50,7,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,1,0,1,0,0,0,0,0,1,0,0,0,0,0,1,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0	0
72,-1,50,5,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,1,0,0,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0	0
.
.
.

データセット

このような構造のインスタンスが173個含まれるデータをラーニングデータセット (classification_02_learning.txt) とする.以下からダウンロードできる.

classification_02_learning.txt

学習による予測器の生成および新規データにおけるテスト

このラーニングデータセットを学習して予測器を作るには以下のようにする.Theanoのロジスティック回帰の公式ページの例を改変.ロジスティック回帰は,少し拡張するだけでニューラルネットワークへと変化させることができるため以下ではニューラルネットワークへと発展させやすいようにコードを書く.

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
import sys,re
import numpy as np
import theano
import theano.tensor as T
np.random.seed(0)
np.set_printoptions(threshold=np.inf,linewidth=np.inf,suppress=True,precision=12)

def main():
	# 1. setting hyper parameters
	EPOCH,TN,VN,TBATCHSIZE,VBATCHSIZE=1000,123,50,10,4
	TBATCHN,VBATCHN=TN//TBATCHSIZE,VN//VBATCHSIZE
	LEARNINGFILE="classification_02_learning.txt"
	
	# 2. reading data
	litrainx,litraint,livalidx,livalidt=[],[],[],[]
	fin=open(LEARNINGFILE,"r")
	for i,line in enumerate(fin):
		litmp=re.split("\t",line.rstrip())
		lix=[float(i) for i in litmp[0].split(",")]
		lit=[int(litmp[1])]
		if i<TN:
			litrainx.append(lix)
			litraint.append(lit)
		elif i<TN+VN:
			livalidx.append(lix)
			livalidt.append(lit)
	fin.close()
	litrainx,litraint,livalidx,livalidt=np.asarray(litrainx,dtype=np.float32),np.asarray(litraint,dtype=np.int32),np.asarray(livalidx,dtype=np.float32),np.asarray(livalidt,dtype=np.int32)
	
	# 3. designing network
	x=T.fmatrix("x")
	t=T.imatrix("t")
	layer=Layer(69,1,T.nnet.sigmoid)
	liparameter=layer.liparameter
	y=layer.forward(x)
	cost=-T.mean(t*T.log(y)+(1-t)*T.log(1-y))+layer.regularization
	ligradient=T.grad(cost=cost,wrt=liparameter)
	liupdate=adam(liparameter,ligradient)
	train=theano.function(inputs=[x,t],outputs=[cost,y],updates=liupdate)
	valid=theano.function(inputs=[x,t],outputs=[cost,y])
	
	# 4. learning data
	for epoch in range(EPOCH):
		index=np.random.permutation(TN)
		traincost,trainacc=0,0
		for i in range(TBATCHN):
			start=i*TBATCHSIZE
			end=start+TBATCHSIZE
			currentcost,currentprediction=train(litrainx[index[start:end]],litraint[index[start:end]])
			traincost+=currentcost
			trainacc+=accuracy(litraint[index[start:end]],map(round,currentprediction))
		index=np.random.permutation(VN)
		validcost,validacc=0,0
		for i in range(VBATCHN):
			start=i*VBATCHSIZE
			end=start+VBATCHSIZE
			currentcost,currenctprediction=valid(livalidx[index[start:end]],livalidt[index[start:end]])
			validcost+=currentcost
			validacc+=accuracy(livalidt[index[start:end]],map(round,currenctprediction))
		print('Epoch {0:6d}: Training cost={1:7.4f}, acc={2:7.4f}, Validation cost={3:7.4f}, acc={4:7.4f}'.format(epoch+1,float(traincost/TBATCHN),float(trainacc/TBATCHN),float(validcost/VBATCHN),float(validacc/VBATCHN)))
	
	# 5. printing importances of the predictor
	print("\nConnection weight: ",np.dot(layer.W.eval()[:,0],np.identity(69,dtype=np.float32)),sep="")
	
	# 6. prediction for new data
	predict=theano.function(inputs=[x],outputs=y)
	newdata=np.asarray(np.random.randn(1,69),dtype=theano.config.floatX)
	print("\nPrediction result: ",round(predict(newdata)),sep="")

class Layer(object):
	def __init__(self,din,dout,activation):
		self.W=theano.shared(value=np.asarray(np.random.randn(din,dout),dtype=theano.config.floatX),name="W")
		self.b=theano.shared(value=np.asarray(np.random.randint(size=dout,low=0,high=10),dtype=theano.config.floatX),name="b")
		self.activation=activation
		self.liparameter=[self.W,self.b]
		self.regularization=0.01*(self.W**2).sum()
	def forward(self,x):
		u=T.dot(x,self.W)+self.b
		self.z=self.activation(u)
		return self.z

def accuracy(lix,liy):
	tp=0
	for x,y in zip(lix,liy):
		if x==y: tp+=1
	return tp/len(lix)

def adam(liparameter,ligradient,a=0.001,b1=0.9,b2=0.999,e=1e-8):
	liupdate=[]
	t=theano.shared(value=np.float32(1),name="t")
	liupdate.append((t,t+1))
	for parameter,gradient in zip(liparameter,ligradient):
		m=theano.shared(value=np.zeros(parameter.get_value().shape,dtype=theano.config.floatX),name='m')
		v=theano.shared(value=np.zeros(parameter.get_value().shape,dtype=theano.config.floatX),name='v')
		m_=b1*m+(1-b1)*gradient
		v_=b2*v+(1-b2)*gradient**2
		m_hat=m_/(1-b1**t)
		v_hat=v_/(1-b2**t)
		parameter_=parameter-(a*m_hat)/(T.sqrt(v_hat)+e)
		liupdate.append((m,m_))
		liupdate.append((v,v_))
		liupdate.append((parameter,parameter_))
	return liupdate

def relu(x):
	return x*(x>0)

def round(x):
	return int(x+0.5)

if __name__ == '__main__':
	main()

以上のプログラムにおいて学習は以下の6ステップからなる.

  1. ハイパーパラメーターの設定.
  2. データの読み込み.
  3. ネットワークの設計.
  4. ミニバッチ法による学習.
  5. フィーチャーの重要度の出力.
  6. 新たなデータに対する予測.

最初に学習のハイパーパラメーターを決める.この場合,学習の最大エポックは1000回,トレーニングデータセットのサイズは123,バリデーションデータセットのサイズは50,ミニバッチのサイズはトレーニングとバリデーションでそれぞれ10と4となる.

	# 1. setting hyper parameters
	EPOCH,TN,VN,TBATCHSIZE,VBATCHSIZE=1000,123,50,10,4
	TBATCHN,VBATCHN=TN//TBATCHSIZE,VN//VBATCHSIZE
	LEARNINGFILE="classification_02_learning.txt"

32行目からのネットワークのデザインの部分では最初に Layer(69,1,T.nnet.sigmoid) で69次元の入力に対して,1次元の出力をし,このとき活性化関数にシグモイド関数を使う指定をする.liparameter には学習のパラメーター,すなわち重み行列 W とバイアスベクトル b が格納される.y=layer.forward(x) では順伝播計算をすることで出力層の値を出す.cost は損失関数を計算して求まる損失.L2ノルム正則化項を加える.ligradient には cost を各パラメーターで偏微分した値が入り,liupdate には最適化手法 Adam によって更新されるパラメーターが入る.それらを train および valid で表す,学習器に入れる.train に値が入力されるとパラメーターが更新されるように updates= の部分に liupdate を指定する.

	# 3. designing network
	x=T.fmatrix("x")
	t=T.imatrix("t")
	layer=Layer(69,1,T.nnet.sigmoid)
	liparameter=layer.liparameter
	y=layer.forward(x)
	cost=-T.mean(t*T.log(y)+(1-t)*T.log(1-y))+layer.regularization
	ligradient=T.grad(cost=cost,wrt=liparameter)
	liupdate=adam(liparameter,ligradient)
	train=theano.function(inputs=[x,t],outputs=[cost,y],updates=liupdate)
	valid=theano.function(inputs=[x,t],outputs=[cost,y])

64行目からの2行は変数の重要度をコネクションウエイト法で計算するための記述.ロジスティック回帰の場合,単に W を出力するだけでいいが,ニューラルネットワークへ拡張したときにも計算しやすいように以下のように書く.

	# 5. printing importances of the predictor
	print("\nConnection weight: ",np.dot(layer.W.eval()[:,0],np.identity(69,dtype=np.float32)),sep="")

90行目の関数は2016年時点で最も優れている最適化法であるAdamの実装.これまでに記憶した勾配の平均と分散を用いて現在の勾配を最適化する手法.

def adam(liparameter,ligradient,a=0.001,b1=0.9,b2=0.999,e=1e-8):
	liupdate=[]
	t=theano.shared(value=np.float32(1),name="t")
	liupdate.append((t,t+1))
	for parameter,gradient in zip(liparameter,ligradient):
		m=theano.shared(value=np.zeros(parameter.get_value().shape,dtype=theano.config.floatX),name='m')
		v=theano.shared(value=np.zeros(parameter.get_value().shape,dtype=theano.config.floatX),name='v')
		m_=b1*m+(1-b1)*gradient
		v_=b2*v+(1-b2)*gradient**2
		m_hat=m_/(1-b1**t)
		v_hat=v_/(1-b2**t)
		parameter_=parameter-(a*m_hat)/(T.sqrt(v_hat)+e)
		liupdate.append((m,m_))
		liupdate.append((v,v_))
		liupdate.append((parameter,parameter_))
	return liupdate

以上のプログラムの名前を lr.py とすると,これは以下のようなコマンドで実行することができる.

$ THEANO_FLAGS=mode=FAST_RUN,device=cpu,floatX=float32 lr.py

結果は以下のようになる.この場合,1000エポックの段階で,トレーニングデータセットにおける正確度は約0.8,バリデーションデータセットにおける正確度は約0.7で過学習は起きていなそうなことがわかる.Connection weight: の部分には69個の各アトリビュートの重要度が出力されており,これの絶対値が大きいほどそのアトリビュートが結果に与える影響が強いことを表す.最後の新たなデータに対する結果は0である.

Epoch      1: Training cost=130.0827, acc= 0.2417, Validation cost=119.5624, acc= 0.3125
Epoch      2: Training cost=129.6084, acc= 0.2417, Validation cost=118.2159, acc= 0.3125
Epoch      3: Training cost=126.6161, acc= 0.2500, Validation cost=117.6241, acc= 0.3125
.
.
.
Epoch    998: Training cost= 0.4380, acc= 0.8167, Validation cost= 0.7641, acc= 0.7083
Epoch    999: Training cost= 0.4505, acc= 0.8083, Validation cost= 0.7414, acc= 0.7292
Epoch   1000: Training cost= 0.4478, acc= 0.8167, Validation cost= 0.7785, acc= 0.7083

Connection weight: [-0.013053309172 -0.773799180984  0.005231603514 -0.119090951979 -0.307564377785 -0.              0.             -0.              0.058809589595 -0.052471533418 -0.191868528724  0.141893491149  0.              0.305930733681  0.             -0.151497423649 -0.114544339478 -0.329956889153 -0.262738287449 -0.401144057512 -0.184271275997 -0.041341982782  0.491266459227  0.199132636189 -0.045660510659 -0.              0.087909527123 -0.470117539167 -0.099915862083 -0.022986473516 -0.518820166588 -0.             -0.              0.125468358397 -0.348094284534  0.             -0.09798681736   0.              0.             -0.008578025736  0.198688134551 -0.408965945244 -0.03944331035   0.             -0.146395936608 -0.059272922575 -0.              0.             -0.14760479331  -0.022076519206 -0.             -0.              0.311165213585  0.415008485317  0.06547921896   0.252530395985  0.             -0.             -0.101379200816  0.542908549309  0.              0.09248958528  -0.             -0.             -0.              0.             -0.130057036877 -0.09130936861  -0.063522666693]

Prediction result: 0
  1. Olden JD et al., An accurate comparison of methods for quantifying variable importance in artificial neural networks using simulated data. Ecological Modeling, 178:389-397,2004.
  2. Theanoロジスティック回帰公式ページ
Hatena Google+