クーの自由研究

かえるのクーが素朴な疑問を実験して報告します。

自由研究の準備(その2)実験コードの初版です

おおまかな内容

最初のプログラムソースコードを覚えておきます。実験の途中でいろいろ変えていきますが、初心にかえるためのものです。

なぜやるのか・どうなると思うか・どこまでやるか

CPUだけでもそこそこ早くまともに動作するものにしたいです。

自由自在に変更実験ができるように、フレームワークは使用せずに、基本的なモジュールだけで動作できるものにしたいです。

 

準備の状況

1パターンで1日中実行とかありえないので、どうにか調整してCPUだけで目的の実験が実用的にできるようにしたかったのです。GPUを買うつもりだったお金はすべて「本」になってしまいました。そしてこれからしばらく「本」か続く見込みなのです。

 

プログラムのソース

 「はてなブログの書き方勉強中」はてな表記にしたのにメニューがでてこない。。何かが違う。

>|Python|

#!/usr/bin/env python
# -*- coding: utf-8 -*-

from __future__ import print_function

# jupyter で行うときは以下の#を外しましょう。
#%matplotlib inline

import matplotlib.pyplot as plt
import numpy as np
import argparse
import time
import math
import cPickle as pickle
from sklearn.datasets import fetch_mldata

np.seterr(all='ignore')
plt.style.use('ggplot')

# テスト用なので、設定価を直接記載 やる気を出すためにパラメータだけ準備します。
par_hidden = 1000 # 隠れ層の数です。ref.青イルカ本2.3 行列Wとbの行数
par_noise = 0.2 # ノイズの割合です。0:0.999です。ref.青イルカ本5.7.2
par_untied = False # 重み共有「しない」かどうかです。Trueで共有しない。ref.青イルカ本5.2.2
par_sub_scale = 100 # 100とすれば1/100のデータ量でepuchをおわります。GPUを買ったら1にします。独自実装:CPUだけでも(精度は無視して)それなりに動く魔法(いんちき)
par_batch_size = 10 # このPGでは1バッチ内に0:9を均等に埋めるので必ず10の倍数にしてください。さもなくば、学習がかたよります。20とすれば、0:9を2コづつバッチにします。ref.青イルカ本3.3
par_epoch = 128 # Epochの回数です。2**n回毎にPeriod(Milestone)をいれてるので64とか128の数字がいいかもです。
par_display_volune = 30 # 中間結果をサンプル表示させる数
par_dropout = 0.4 # Todo 未実装です。そのうちどこからから切り貼りします。ref.青イルカ本3.5.3
par_beta = 0.0 # Todo 青イルカ本5.4のβです。未実装です。
par_train_shuffle = False # Todo 学習時にデータ取得をシャッフルするかどうか ref.青イルカ本3.6.7

# 以下 Todo & 妄想
# 強化学習, 畳みこみ, プーリング, スパース正則化, 階層化・深層化
# 蒸留(distillation)用学習済先生パラメータの定義・取り込み 蒸留って量子力学用語?
# 量子化係数, Huffman Coding, 圧縮手法, Low-Rank Regulation, LearningUsingPrevillageInformation,
# DataDependentInitialization, GeneralizedDistillation=Distillation+LUPI,
# VAEのZ軸フィッティングのためのパラメータ M1+M2
# ニューロ自体の動的追加・削除・統合・分離、分散学習と統合、
# 学習時ニューロのノイズもしくは生物学的要素の再導入試行(1/f揺らぎ)

# ここから本体
# MNIST のデータを準備します。
print('load MNIST dataset')
print('---Start---')
mnist = fetch_mldata('MNIST original')

# mnist.data : 70,000件の784次元ベクトルデータ
mnist.data = mnist.data.astype(np.float32)
mnist.data /= 255 # 0-1のデータに変換

# mnist.target : 正解データ(教師データ)AutoEncoderでは使いません。
mnist.target = mnist.target.astype(np.int32)

N = 60000
y_train, y_test = np.split(mnist.data.copy(), [N])
N_test = y_test.shape[0]

x_train, x_test = np.split(mnist.data, [N])


# AutoEncoderを定義します。
# --------------------------------------------------------------------------------------------
class KooAutoEncoder(object):
# 定数
k_mnist_start_index = [0, 5923, 12665, 18623, 24754, 30596, 36017, 41935, 48200, 54051, 60000] # MNISTの0:9範囲
# 以下の定数は妄想です。
# Todo 白色化は未実装です。こうやって書いておけばいつかは実装されるかもしれません。 ref.青イルカ本5.5
k_NO_whitening = 'NO' # 白色化しない
k_ZCA_whitening = 'ZCA' # ZCA白色化
k_PCA_whitening = 'PCA' # PCA白色化
# Todo そのうち活性化関数を切り替えれるようにするためのなにか。ref.青イルカ本2.2
k_activation_sigmoid = 'sigmoid'
k_activation_step = 'step'
k_activation_relu = 'relu'
k_activation_softsign = 'softsign'
k_activation_softplus = 'softplus'
k_activation_maxout = 'maxout'
k_activation_tanx = 'tanx'
k_activation_prelu = 'prelu'
k_activation_leakyrelu = 'leakyrelu'
k_activation_elus = 'elus'
# hybrid
k_activation_tranctanx = 'tranctanx'
k_activation_leakytanx = 'leakytanx'
k_activation_steprelu = 'steprelu'
# Todo そのうちオプティマイザを切り替えれるようにするためのなにか。ref.青イルカ本3,4
k_optimize_SGD = "SGD"
k_optimize_AdaGrad = "AdaGrad"
k_optimize_Adam = "Ada"
k_optimize_RMSprop = "RMSprop"
k_optimize_AdaDelta = "AdaDelta"
k_optimize_Adamax = "Adamax"

def __init__(self, n_visible = 784, n_hidden = 784,
W1 = None, W2 = None, b1 =None, b2 = None,
noise = 0.0, untied = False, dropout = 0.0, sub_scale = 10):

self.rng = np.random.RandomState(1)

# 中間結果を後で確認するための記録用。メモリ使用量がはんぱなくなるので、とりあえずPeriodで覚えてみます。
self.W1rec = []
self.W2rec = []
self.b1rec = []
self.b2rec = []

initRange = np.sqrt(6. / (n_hidden + n_visible + 1)) # 6.は調査中。いまのところ魔法の数字

# Wの初期化 初期化方法は実験のしがいがありそうです。ref.青イルカ本3.6.6
if W1 == None:
self.W1 = self.random_init(initRange, (n_hidden, n_visible))

if W2 == None:
if untied:
W2 = self.random_init(initRange, (n_visible, n_hidden))
else:
W2 = self.W1.T

self.W2 = W2

# bの初期化 お約束のゼロスタートです。ref.青イルカ本3.6.6
if b1 == None:
self.b1 = np.zeros(n_hidden)
if b2 == None:
self.b2 = np.zeros(n_visible)

# 活性化関数を設定します。
self.func = Sigmoid()

# パラメータを覚えておきます。
self.n_visible = n_visible
self.n_hidden = n_hidden
self.alpha = 0.1
self.noise = noise
self.untied = untied
self.dropout = dropout
self.sub_scale = sub_scale
self.entropyrec = []

def random_init(self, r, size):
return np.array(self.rng.uniform(low = -r, high = r, size=size))

def corrupt(self, x, noise):
return self.rng.binomial(size = x.shape, n = 1, p = 1.0 - noise) * x

def dropW(self, p):
# Todo dropoutは未実装
dWone = (self.rng.uniform(low=-(p * 100), high=(1 - p) * 100, size=self.n_hidden) > 0) * 1
dW = np.vstack((dWone, dWone) * (self.n_visible / 2))
return dW

def encode(self, x, dropW):
# ref.青イルカ本5.1
# Todo dropoutは未実装
#return self.func.activate_func(np.dot(dropW.T * self.W1, x) + (dropW[0] * self.b1))
return self.func.activate_func(np.dot(self.W1, x) + (self.b1))

def encode_by_snap(self, W, b, x):
# あとから中間結果を確認するためのencode
return self.func.activate_func(np.dot(W, x) + b)

def decode(self, y, dropW):
# ref.青イルカ本5.1
# Todo dropoutは未実装
#return self.func.activate_func(np.dot(dropW * self.W2, y) + (self.b2))
return self.func.activate_func(np.dot(self.W2, y) + (self.b2))

def decode_by_snap(self, W, b, y):
# あとから中間結果を確認するためのdecode
return self.func.activate_func(np.dot(W, y) + b)

def get_entropy(self, x, z):
eps = 1e-10
# ref.青イルカ本2.4.3式(2.8) もしくは 5.2.1
return - np.sum((x * np.log(z + eps) + (1.-x) * np.log(1.-z + eps)))

def get_entropy_and_grad(self, x_batch, dnum):
entropy = 0.
dropW = [] #Dummy
grad_W1 = np.zeros(self.W1.shape)
grad_W2 = np.zeros(self.W2.shape)
grad_b1 = np.zeros(self.b1.shape)
grad_b2 = np.zeros(self.b2.shape)

for x in x_batch:
tilde_x = self.corrupt(x, self.noise)
# Todo DropOutは苦戦中(現在機能なし)
#dropW = self.dropW(self.dropout)
p = self.encode(tilde_x, dropW)
y = self.decode(p, dropW)

entropy += self.get_entropy(x, y)

delta1 = - (x - y)

if self.untied:
grad_W2 += np.outer(delta1, p)
else:
grad_W1 += np.outer(delta1, p).T

grad_b2 += delta1

delta2 = np.dot(self.W2.T, delta1) * self.func.differential_func(p)
grad_W1 += np.outer(delta2, tilde_x)
grad_b1 += delta2

entropy /= len(x_batch)
grad_W1 /= len(x_batch)
grad_W2 /= len(x_batch)
grad_b1 /= len(x_batch)
grad_b2 /= len(x_batch)

return entropy, grad_W1, grad_W2, grad_b1, grad_b2


def train(self, X, epochs = 15, batch_size = 20):
batch_num = len(X) / (batch_size * self.sub_scale)
# epochは表示は1から、内部は0からにしています。
for epoch in range(epochs):
total_entropy = 0.0
start_time = time.clock()
for i in range(batch_num):
# CPUなので実際の1/Subscaleで回す。(実際には必要のない余分なロジック)
offset = epoch * batch_num * (batch_size / 10) + batch_size / 10 * i
subbatch = batch_size / 10
batch = X[self.k_mnist_start_index[0] + offset:self.k_mnist_start_index[0] + offset + subbatch]
# 数字が均等にバッチに入るようにしている。
for j in range(10)[1:]:
batch = np.concatenate([batch, X[self.k_mnist_start_index[j] + offset:self.k_mnist_start_index[j] + offset + subbatch]], axis = 0)
# GPUならまともに流しましょう。
#batch = X[i*batch_size : (i+1)*batch_size]
# CPUなので実際の1/Subscaleで回す。(実際には必要のない余分なロジック)ここまで

entropy, gradW1, gradW2, gradb1, gradb2 = self.get_entropy_and_grad(batch, len(X))

total_entropy += entropy
self.W1 -= self.alpha * gradW1
self.W2 -= self.alpha * gradW2
self.b1 -= self.alpha * gradb1
self.b2 -= self.alpha * gradb2

grad_sum = gradW1.sum() + gradW2.sum() + gradb1.sum() + gradb2.sum()
#print ("\x1b[1Aage:batch/Total = %d/%d" % (i + 1, batch_num))

end_time = time.clock()
entropy = (1. / (batch_num * self.sub_scale)) * total_entropy
print ("epoch = %d: entropy = %f: \ttime = %.3f" % (epoch + 1, entropy, (end_time - start_time)))
self.entropyrec.append(entropy)
# period 2**n - 1 で途中結果を記録
if epoch in (2 ** np.array(range(20)) - 1):
print ("period = %d ----------" % (np.int(np.log2(epoch + 1)) + 1))
self.W1rec.append(self.W1.copy())
self.W2rec.append(self.W2.copy())
self.b1rec.append(self.b1.copy())
self.b2rec.append(self.b2.copy())


def dump_weights(self, save_path):
with open(save_path, 'w') as f:
d = {
"W1" : self.W1,
"W2" : self.W2,
"b1" : self.b1,
"b2" : self.b2,
}

pickle.dump(d, f)

# 中間結果取得用
def get_W1(self, i):
return self.W1rec[i]

def get_W2(self, i):
return self.W2rec[i]

def get_b1(self, i):
return self.b1rec[i]

def get_b2(self, i):
return self.b2rec[i]

def get_entropyrec(self):
return self.entropyrec

# MNIST開始INDEX取得(固定値)
def get_mnist_start_index(self, i):
return self.k_mnist_start_index[i]
# --------------------------------------------------------------------------------------------

# 活性化関数の元クラスです。
# --------------------------------------------------------------------------------------------
class ActivateFunction(object):
# __metaclass__ = ActivateFunctionMeta
def __init__(self):
self.name = self.__class__.__name__

@classmethod
#@abstractmethod
def activate_func(cls, x):
raise NotImplementedError()

@classmethod
#@abstractmethod
def differential_func(cls, x):
raise NotImplementedError()
# --------------------------------------------------------------------------------------------

# 活性化関数を定義していきます。
# --------------------------------------------------------------------------------------------
class Sigmoid(ActivateFunction):
# ref.青イルカ本2.2
def activate_func(cls, x):
return 1. / (1. + np.exp(-x))

def differential_func(cls, x):
return x * (1. - x)
# --------------------------------------------------------------------------------------------

class ReLU(ActivateFunction):
# ref.青イルカ本2.2
def activate_func(cls, x):
return x * (x > 0)

def differential_func(cls, x):
return x * (1. - x)
# --------------------------------------------------------------------------------------------

# オートエンコーダを準備して学習実行します。
ae = KooAutoEncoder(n_hidden = par_hidden, noise = par_noise, untied = par_untied, sub_scale = par_sub_scale, dropout = par_dropout)
try:
# 学習の実行
ae.train(x_train, epochs = par_epoch, batch_size = par_batch_size)
except KeyboardInterrupt:
exit()
pass

# 学習は以上で終了です。
# 中間結果を表示させてみます。

# draw digit images
def draw_digit(data, i, row, colmn):
size = 28
plt.subplot(row + 1, colmn, i)
Z = data.reshape(size, size) # convert from vector to 28x28 matrix
Z = Z[::-1, :] # flip vertical
plt.xlim(0, size)
plt.ylim(0, size)
plt.pcolor(Z)
plt.title("%d" % i, size=9)
plt.gray()
plt.tick_params(labelbottom="off")
plt.tick_params(labelleft="off")

for period in range(np.int(np.log2(par_epoch) + 1)):
plt.figure(figsize=(10, 10))
targetW = ae.get_W1(period)
# たくさん表示すると遅いので高々par_display_volune個表示にします。
displayVolume = min((len(targetW), par_display_volune))
for i in range(displayVolume):
draw_digit(targetW[i], i + 1, 5, 10)
print ("W1 Fig. Period %d (Epoch %d)" % (period + 1, 2 ** period))
plt.show()

# これ以降をvar_offsetをかえてjupyterなんかでいろいろ試せば感じがわかるはずです。
# var_offsetはMNISTデータへのインデックスのオフセットに使用します。
var_offset = 0
for period in range(np.int(np.log2(par_epoch) + 1)):
plt.figure(figsize=(10, 10))
for sample in range(10):
y_x = ae.encode_by_snap(ae.get_W1(period), ae.get_b1(period), x_train[ae.get_mnist_start_index(sample) + var_offset])
y = ae.decode_by_snap(ae.get_W2(period), ae.get_b2(period), y_x)
draw_digit(y, sample + 1, 5, 10)
print ("x_hat Fig. by Period %d (Epoch %d)" % (period + 1, 2 ** period))
plt.show()
# Original Mnist Image
print (" Fig. Original MNIST Data")
plt.figure(figsize=(10, 10))
for sample in range(10):
draw_digit(x_train[ae.get_mnist_start_index(sample) + var_offset], sample + 1, 5, 10)
plt.show()


plt.style.use('ggplot')
plt.figure(figsize=(10,7))
plt.plot(ae.get_entropyrec(), lw=1)
plt.title("")
plt.ylabel("entropy")
plt.xlabel("epoch")
plt.show()
# --------------------------------------------------------------------------------------------

||<

何かが違う。。。

テスト用なので、main形式にしていないです。

まとめ(わかったことわからなかったこと)

長くてすみません。

あと、ソースの貼り方がわからなかったのですが、それなりに読めるので今回はこのようにしておきます。(わからないので断念)

 あ、そうそう、DropOutをいれようとおもいましたが、うまくいかなかったので宿題にしました。時間はかかるは、収束せず発散していくわ、ダメダメでした。勉強しなおします。

感想・思ったこと・考えたこと

なんといってもGPUをもっていないことにつきます。

ですが、実験したいことろはGPUでぶんぶん回す「ところ」でないことも多いので、当面CPUだけでがんばります。

本で読んだのと(切り貼りとはいえ)プログラムしてみたのでは大違いでした。デバックのやり方をなんとかマスターしないと。。。

テストはまだほとんどできていないので、パラメータを変えればいろいろ問題があるかもしれません。実験プログラムなので、随時改定していきます。

そして、改めてネットの諸先生に感謝です。

 

参考にした資料

pythonの本いろいろ。

資料ではないですが、Udemyでがぜんpythonというよりnumpy, pandas勉強中です。。。

 

「蛇の足」は「セルピエンテ・タコーン」というらしい

 禅の心はpythonの心。合掌。

気になる実行結果は自由研究の準備(その3)をご期待ください。

こわごわ、機械学習の「機械学習を崇拝する人々の集い 」へ参加してみました。「崇拝する」に惹かれ、「人々」にしりごみしました。「かえる」はNGなら抜けようと思います。

(追記)グループの他の方のレベルの高さとあまり差があるので、かえるのクーはひとりに「かえる」ことにしました。自由研究というスタイル(そしてページのスタイル)は思いのほか「浮いている」ことを自覚した次第です。でも初志貫徹します!!