In [1]:
%matplotlib inline

import numpy as np
import pylab as pl
import math 
from sympy import *
import matplotlib.pyplot as plt
import matplotlib.animation as animation
from mpl_toolkits.mplot3d import Axes3D

In [2]:
from NN import NN

コブ・ダクラス型生産関数と課題文で例に出された関数を用いる。

いずれも定義域は0≤x≤1である。

コブ・ダグラス型生産関数は以下の通りである。

z = x_1**0.5*x_2*0.5


In [3]:
def example1(x_1, x_2):
    z = x_1**0.5*x_2*0.5
    return z

fig = pl.figure()
ax = Axes3D(fig)

X = np.arange(0, 1, 0.1)
Y = np.arange(0, 1, 0.1)

X, Y = np.meshgrid(X, Y)
Z = example1(X, Y)

ax.plot_surface(X, Y, Z, rstride=1, cstride=1)

pl.show()


/Users/nigg/anaconda/lib/python2.7/site-packages/matplotlib/collections.py:590: FutureWarning: elementwise comparison failed; returning scalar instead, but in the future will perform elementwise comparison
  if self._edgecolors == str('face'):

課題の例で使われた関数は以下の通りである。

z = (1+np.sin(4*math.pi*x_1))*x_2*1/2


In [4]:
def example2(x_1, x_2):
    z = (1+np.sin(4*math.pi*x_1))*x_2*1/2
    return z

fig = pl.figure()
ax = Axes3D(fig)

X = np.arange(0, 1, 0.1)
Y = np.arange(0, 1, 0.1)

X, Y = np.meshgrid(X, Y)
Z = example2(X, Y)

ax.plot_surface(X, Y, Z, rstride=1, cstride=1)
ax.set_zlim(-1, 1)

pl.show()


NNのクラスはすでにNN.pyからimportしてある。


In [5]:
nn = NN()

以下に使い方を説明する。

初めに、このコブ・ダグラス型生産関数を用いる。


In [6]:
x_1 = Symbol('x_1')
x_2 = Symbol('x_2')
f = x_1**0.5*x_2*0.5

入力層、中間層、出力層を作る関数を実行する。引数には層の数を用いる。


In [7]:
nn.set_input_layer(2)


Out[7]:
array([ 0.,  0.])

In [8]:
nn.set_hidden_layer(2)


Out[8]:
(array([ 0.,  0.]), array([ 0.,  0.]))

In [9]:
nn.set_output_layer(2)


Out[9]:
(array([ 0.,  0.]), array([ 0.,  0.]), array([ 0.,  0.]))

nn.set_hidden_layer()は同時にシグモイド関数で変換する前の中間層も作る。

set_output_layer()は同時にシグモイド関数で変換する前の出力層、さらに教師データを入れる配列も作る。

nn.setup()で入力層ー中間層、中間層ー出力層間の重みを入れる配列を作成する。

nn.initialize()で重みを初期化する。重みは-1/√d ≤ w ≤ 1/√d (dは入力層及び中間層の数)の範囲で一様分布から決定される。


In [10]:
nn.setup()
nn.initialize()

nn.supervised_function(f, idata)は教師データを作成する。引数は関数とサンプルデータをとる。


In [11]:
idata = [1, 2]

In [12]:
nn.supervised_function(f, idata)

nn.simulate(N, eta)は引数に更新回数と学習率をとる。普通はN=1で行うべきかもしれないが、工夫として作成してみた。N回学習した後に出力層を返す。


In [13]:
nn.simulate(1, 0.1)


Out[13]:
array([ 0.53414091,  0.4764263 ])

nn.calculation()は学習せずに入力層から出力層の計算を行う。nn.simulate()内にも用いられている。

次に実際に学習を行う。サンプルデータは、


In [14]:
X = np.arange(0, 1, 0.2)
Y = np.arange(0, 1, 0.2)
print X, Y


[ 0.   0.2  0.4  0.6  0.8] [ 0.   0.2  0.4  0.6  0.8]

の組み合わせである。


In [27]:
X = np.arange(0, 1, 0.2)
Y = np.arange(0, 1, 0.2)

a = np.array([])
b = np.array([])
c = np.array([])

nn = NN()
nn.set_network()

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)
        
for i in range(100):
        l = np.random.choice([i for i in range(len(a))])
        m = nn.main(1, f, [a[l], b[l]], 0.5)

for x in X:
    for y in Y:
        idata = [x, y]
        c = np.append(c, nn.realize(f, idata))

In [28]:
a


Out[28]:
array([ 0. ,  0. ,  0. ,  0. ,  0. ,  0.2,  0.2,  0.2,  0.2,  0.2,  0.4,
        0.4,  0.4,  0.4,  0.4,  0.6,  0.6,  0.6,  0.6,  0.6,  0.8,  0.8,
        0.8,  0.8,  0.8])

In [29]:
b


Out[29]:
array([ 0. ,  0.2,  0.4,  0.6,  0.8,  0. ,  0.2,  0.4,  0.6,  0.8,  0. ,
        0.2,  0.4,  0.6,  0.8,  0. ,  0.2,  0.4,  0.6,  0.8,  0. ,  0.2,
        0.4,  0.6,  0.8])

In [30]:
c


Out[30]:
array([ 0.62110062,  0.63653367,  0.65111218,  0.66436166,  0.67598904,
        0.63172135,  0.64662299,  0.66032343,  0.67247157,  0.68290415,
        0.64201944,  0.65613827,  0.66879249,  0.67976226,  0.68900194,
        0.65181452,  0.66495474,  0.67645812,  0.68622916,  0.69431931,
        0.66096283,  0.67299149,  0.68329938,  0.69189792,  0.69891149])

例えば(0, 0)を入力すると0.52328635を返している(つまりa[0]とb[0]を入力して、c[0]の値を返している)。

ここでは交差検定は用いていない。


In [31]:
fig = pl.figure()
ax = Axes3D(fig)

ax.scatter(a, b, c)

pl.show()


確率的勾配降下法を100回繰り返したが見た感じから近づいている。回数を10000回に増やしてみる。


In [32]:
X = np.arange(0, 1, 0.2)
Y = np.arange(0, 1, 0.2)

a = np.array([])
b = np.array([])
c = np.array([])

nn = NN()
nn.set_network()

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)
        
for i in range(10000):
        l = np.random.choice([i for i in range(len(a))])
        m = nn.main(1, f, [a[l], b[l]], 0.5)

for x in X:
    for y in Y:
        idata = [x, y]
        c = np.append(c, nn.realize(f, idata))

In [33]:
fig = pl.figure()
ax = Axes3D(fig)

ax.scatter(a, b, c)

pl.show()


見た感じ随分近づいているように見える。

同様のことを課題の例で使われた関数でも試してみる。


In [34]:
x_1 = Symbol('x_1')
x_2 = Symbol('x_2')
f = (1+sin(4*math.pi*x_1))*x_2*1/2

In [35]:
X = np.arange(0, 1, 0.2)
Y = np.arange(0, 1, 0.2)

a = np.array([])
b = np.array([])
c = np.array([])

nn = NN()
nn.set_network()

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)
        
for i in range(1000):
        l = np.random.choice([i for i in range(len(a))])
        m = nn.main(1, f, [a[l], b[l]], 0.5)

for x in X:
    for y in Y:
        idata = [x, y]
        c = np.append(c, nn.realize(f, idata))

In [36]:
fig = pl.figure()
ax = Axes3D(fig)

ax.scatter(a, b, c)

pl.show()


上手く近似できないので中間層の数を変えてみる。5層にしてみる。


In [37]:
X = np.arange(0, 1, 0.2)
Y = np.arange(0, 1, 0.2)

a = np.array([])
b = np.array([])
c = np.array([])

nn = NN()
nn.set_network(h=5)

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)
        
for i in range(1000):
        l = np.random.choice([i for i in range(len(a))])
        m = nn.main(1, f, [a[l], b[l]], 0.5)

for x in X:
    for y in Y:
        idata = [x, y]
        c = np.append(c, nn.realize(f, idata))

In [38]:
fig = pl.figure()
ax = Axes3D(fig)

ax.scatter(a, b, c)

pl.show()


目標と比べると大きく異なる。

サンプル数を、


In [ ]:
X = np.arange(-1, 1, 0.1)
Y = np.arange(-1, 1, 0.1)
print X, Y

で取り、学習の際にランダムに一個選ばれたサンプルを何十回も学習させてみた。


In [45]:
fig = pl.figure()
ax = Axes3D(fig)

f = (1+sin(4*math.pi*x_1))*x_2*1/2

X = np.arange(-1, 1, 0.1)
Y = np.arange(-1, 1, 0.1)

a = np.array([])
b = np.array([])
c = np.array([])

fig = plt.figure()

nn = NN()
nn.set_network()

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)
        c = np.append(c, nn.main2(50, f, [x, y], 0.8))

for i in range(50):
    l = np.random.choice([i for i in range(len(a))])
    m = nn.main2(20, f, [a[l], b[l]], 0.5)
    c[l] = m

a = np.array([])
b = np.array([])
c = np.array([])

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)
        c = np.append(c, nn.realize(f, [x, y]))
        
ax.scatter(a, b, c)
ax.set_zlim(0, 1)

pl.show()


<matplotlib.figure.Figure at 0x108854810>

本来ならば下のような形になるべきであるので上手くいっているとは言い難い。


In [46]:
def example2(x_1, x_2):
    z = (1+np.sin(4*math.pi*x_1))*x_2*1/2
    return z

fig = pl.figure()
ax = Axes3D(fig)

X = np.arange(-1, 1, 0.1)
Y = np.arange(-1, 1, 0.1)

X, Y = np.meshgrid(X, Y)
Z = example2(X, Y)

ax.plot_surface(X, Y, Z, rstride=1, cstride=1)
ax.set_zlim(-1, 1)

pl.show()


同じ方法でコブ・ダグラス型生産関数を学習させた様子をアニメーションにしてみた。この方法が何の意味を持つかは分からないが学習はまあまよくできていた。

mp4ファイルで添付した。

[考察]ここにはコードを書けなかったが、サンプル数が多くなった時は、ミニバッチ式の学習を試すべきだと思った。

最後に交差検定を行う。

初めに学習回数が極めて少ないNNである。


In [92]:
X = np.arange(0, 1, 0.2)
Y = np.arange(0, 1, 0.2)

a = np.array([])
b = np.array([])
c = np.array([])

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)

evl = np.array([])

for i in range(len(a)):
    nn = NN()
    nn.set_network()
    for j in range(1):
        l = np.random.choice([i for i in range(len(a))])
        if l != i:
            m = nn.main(1, f, [a[l], b[l]], 0.5)
    idata = [a[i], b[i]]
    est = nn.realize(f, idata)
    evl = np.append(evl, math.fabs(est - nn.supervised_data))

In [93]:
np.average(evl)


Out[93]:
0.31083460645411465

次に十分大きく(1000回に)してみる。


In [102]:
X = np.arange(0, 1, 0.2)
Y = np.arange(0, 1, 0.2)

a = np.array([])
b = np.array([])
c = np.array([])

nn = NN()
nn.set_network(h=7)

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)

evl = np.array([])

for i in range(len(a)):
    for j in range(10000):
        nn = NN()
        nn.set_network()
        l = np.random.choice([i for i in range(len(a))])
        if l != i:
            m = nn.main(1, f, [a[l], b[l]], 0.5)
    idata = [a[i], b[i]]
    evl = np.append(evl, math.fabs(nn.realize(f, idata) - nn.supervised_data))


---------------------------------------------------------------------------
KeyboardInterrupt                         Traceback (most recent call last)
<ipython-input-102-9daa2d1fe646> in <module>()
     22         l = np.random.choice([i for i in range(len(a))])
     23         if l != i:
---> 24             m = nn.main(1, f, [a[l], b[l]], 0.5)
     25     idata = [a[i], b[i]]
     26     evl = np.append(evl, math.fabs(nn.realize(f, idata) - nn.supervised_data))

/Users/nigg/MachineLearning/NN.py in main(self, N, f, idata, eta, i, h, o)
    125         self.initialize()
    126         self.supervised_function(f, idata)
--> 127         self.simulate(N, eta)
    128         return self.output_layer[0]
    129 

/Users/nigg/MachineLearning/NN.py in simulate(self, N, eta)
    112         self.eta = eta
    113         for i in range(N):
--> 114             self.calculation()
    115             self.output_ad()
    116             self.calculation()

/Users/nigg/MachineLearning/NN.py in calculation(self)
     67 
     68     def calculation(self):
---> 69         u = Symbol("u")
     70         hfunction = self.hfunction
     71         diff_hf = self.diff_hf

/Users/nigg/anaconda/lib/python2.7/site-packages/sympy/core/symbol.pyc in __new__(cls, name, **assumptions)
     99 
    100         """
--> 101         cls._sanitize(assumptions, cls)
    102         return Symbol.__xnew_cached_(cls, name, **assumptions)
    103 

/Users/nigg/anaconda/lib/python2.7/site-packages/sympy/core/symbol.pyc in _sanitize(assumptions, obj)
     60 
     61         # be strict about commutativity
---> 62         is_commutative = fuzzy_bool(assumptions.get('commutative', True))
     63         if is_commutative is None:
     64             whose = '%s ' % obj.__name__ if obj else ''

/Users/nigg/anaconda/lib/python2.7/site-packages/sympy/core/logic.pyc in fuzzy_bool(x)
     69     if x is None:
     70         return None
---> 71     return bool(x)
     72 
     73 

KeyboardInterrupt: 

In [106]:
evl


Out[106]:
array([ 0.17422959,  0.37774617,  0.29845184,  0.4247869 ,  0.42139651,
        0.28095943,  0.48263009,  0.12094186,  0.35090278,  0.07039276,
        0.60270475,  0.38029955,  0.47076667,  0.27789737,  0.26027828,
        0.61710335,  0.51095281,  0.29832648,  0.51007664,  0.53242457,
        0.50885681,  0.48114824,  0.53944263,  0.36336567,  0.37925241])

In [104]:
np.average(evl)


Out[104]:
0.28091957975738507

誤差の平均であるので小さい方よい。 学習回数を大幅に増やした結果、精度が上がった。

次に回数を100回にして、中間層の数を2と5で精度を比較する。


In [105]:
X = np.arange(0, 1, 0.2)
Y = np.arange(0, 1, 0.2)

a = np.array([])
b = np.array([])
c = np.array([])

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)

evl = np.array([])

for i in range(len(a)):
    for j in range(100):
        nn = NN()
        nn.set_network()
        l = np.random.choice([i for i in range(len(a))])
        if l != i:
            m = nn.main(1, f, [a[l], b[l]], 0.5)
    idata = [a[i], b[i]]
    est = nn.realize(f, idata)
    evl = np.append(evl, math.fabs(est - nn.supervised_data))

In [96]:
np.average(evl)


Out[96]:
0.3090473511199906

In [97]:
X = np.arange(0, 1, 0.2)
Y = np.arange(0, 1, 0.2)

a = np.array([])
b = np.array([])
c = np.array([])

for x in X:
    for y in Y:
        a = np.append(a, x)
        b = np.append(b, y)

evl = np.array([])

for i in range(len(a)):
    for j in range(100):
        nn = NN()
        nn.set_network(h=5)
        l = np.random.choice([i for i in range(len(a))])
        if l != i:
            m = nn.main(1, f, [a[l], b[l]], 0.5)
    idata = [a[i], b[i]]
    est = nn.realize(f, idata)
    evl = np.append(evl, math.fabs(est - nn.supervised_data))

In [98]:
np.average(evl)


Out[98]:
0.28516243998794422

コブ・ダグラス型の関数が比較的単純であるから、意味はないかと思ったが、精度はかなり良くなった。


In [ ]: