ノートブック

古き良きスタイルのノートブックです。

NNの基本とガウス型関数による関数近似

f:id:primitivem:20170929131332p:plain

この記事では、ニューラルネットワークを使ったガウス型関数による任意の関数への近似を行います。 それにあたりニューラルネットワークの基本的なことを確認するため計算の粗過程を1ステップずつ手計算で確認してみたいと思います。 またニューラルネットワークを計算するのにchainerというライブラリを使用しました。

ライブラリの使用法や導入についてはこの記事のテーマではないので割愛しますが、この記事に登場するスクリプトはすべて 以下のインポート文の元おこなっています。

# -*- coding: utf-8 -*-
import random,math
import numpy as np
import matplotlib.pyplot as plt
from chainer import Function, gradient_check, Variable, FunctionSet, optimizers, serializers, utils
from chainer import Link, Chain, ChainList
import chainer.functions as F
import chainer.links as L
import chainer.computational_graph as c

さらに本記事の内容と関連してchainerでの独自関数の定義方法についての公式ドキュメントの解説と 合成関数の微分に関する連鎖律が参考になります。

以下でいろいろな構造のニューラルネットワークの例を出してそれついてその計算法を確認してゆきます。 なおニューラルネットワークや勾配法による最適化問題についてのそもそものコンセプトの説明は本記事では割愛し端的にモデル例と計算だけを示してゆきます。

線型結合のみの構造

ではまず入力と出力の間に活性化関数がなく(恒等写像のみとも言える)線型結合だけの関係性になっているシンプルな以下のような例。

f:id:primitivem:20170929131242p:plain

この構造の教師信号に対する勾配法の過程1ステップを手計算で実際に行い、それと計算機の結果が同じになることを追ってみます。

数理モデルとしてはこうなります。

Y=wX+b
L=(Y-t)^2

これは単純な線型結合モデルとそのバックプロパゲーション学習のためのロス関数(エラー関数) になっています。 そしてこのロス関数 L微分は以下のようになります。

\frac{\partial}{\partial w} L = \frac{\partial Y}{\partial w} \frac{\partial L}{\partial Y}
\frac{\partial}{\partial b} L = \frac{\partial Y}{\partial b} \frac{\partial L}{\partial Y}

すなわちパラメータ  w b を勾配法で以下のようにして最適化します。 1ステップの最適化プロセスで矢印の右側が施され、矢印の左側にパラメータが変化します。そしてこのステップを繰り返すことで教師信号に近づけてゆきます。

 w \leftarrow w - \eta \frac{\partial L}{\partial w}
 b \leftarrow b - \eta \frac{\partial L}{\partial b}

ここで  \eta は学習率

 L の解析解は以下のようになります。

\frac{\partial}{\partial w} L = 2wx^2+2bx-2xt
\frac{\partial}{\partial b} L = 2wx-2t+2b

では実際に適当な数値をこの式に代入して計算してみましょう。
例えば w=1,b=0,トレーニングデータが (x,t)=(1,2) のとき、 \eta=0.0 で, L微分

\frac{\partial}{\partial w} L = -2
\frac{\partial}{\partial b} L = -2

このように計算できて、勾配法(バックプロパゲーション)による計算の一番最初のステップは

w \leftarrow 1-0.01\times(-2) = 1.02
w \leftarrow 0-0.01\times(-2) = 0.02

このようになります。 ではこの手計算と同じことをchainerを使って計算して比較してみましょう。 スクリプトは以下のようになり

model = L.Linear(1,1)
optimizer = optimizers.SGD(lr=0.01)
optimizer.setup(model)
x = Variable(np.array([[1]], dtype=np.float32))
t = Variable(np.array([[2]], dtype=np.float32))
model.W.data = np.array([[1]], dtype=np.float32)
model.b.data = np.array([0], dtype=np.float32)
optimizer.zero_grads()
y = model(x)
loss = F.mean_squared_error(y, t)
loss.backward()
optimizer.update()
print model.W.data  #[[ 1.01999998]]と出力される
print model.b.data  #[ 0.02]と出力される

"print model.W.data" と "print model.b.data" がちょうど最初の1ステップのバックプロパゲーション過程の計算を終えたあとの結果を出力します。 そしてこの数値が先ほどの手計算と同じであることが確認できます。

線型結合に活性化関数を加えた場合

次に先ほどの恒等写像の代わりにシグモイド関数写像に使った場合を見てみましょう。 やりかたは先程とおなじでイメージはこんな感じになり、

f:id:primitivem:20170929140718p:plain

数理モデルは以下のようになり、

Y = w \sigma(X) + b

シグモイド関数は以下のようになります。

\sigma(x) = \frac{1}{1+exp(-x^2)}

先ほどと同じロス関数(Mean Squar Error)を使うと L微分は以下のようになります。

\frac{\partial}{\partial w} L = \frac{\partial Y}{\partial w} \frac{\partial L}{\partial Y} = \sigma(X) \cdot 2(Y-t)
\frac{\partial}{\partial b} L = \frac{\partial Y}{\partial b} \frac{\partial L}{\partial Y} = 1 \cdot 2(Y-t)

では適当な数値を代入して1ステップ計算してみましょう。
例えば w=1,b=0,\eta=0.01,(x,t)=(1,2) のときだと、

\sigma(1) = 0.73
\frac{\partial}{\partial w} L = -1.8542
\frac{\partial}{\partial b} L = -2.56

このような L微分への代入結果になりました。 この手計算と同じことをスクリプトで以下のように計算してみると、

model = L.Linear(1,1)
optimizer = optimizers.SGD(lr=0.01)
optimizer.setup(model)
x = Variable(np.array([[1]], dtype=np.float32))
t = Variable(np.array([[2]], dtype=np.float32))
model.W.data = np.array([[1]], dtype=np.float32)
model.b.data = np.array([0], dtype=np.float32)
optimizer.zero_grads()
h = F.sigmoid(x)
y = model(h)
loss = F.mean_squared_error(y, t)
loss.backward()
optimizer.update()
print model.W.data  #[[ 1.01855338]]
print model.b.data  #[ 0.02537883]
print model.W.grad  #[[-1.85534108]]
print model.b.grad  #[-2.5378828]
print loss.grad  #1.0
print loss.data  #1.61021232605

"print model.l1.W.data" と "print model.l2.W.data" この部分での出力結果が手計算と全く同じ結果になることが確認できます。

線型結合2層

つぎはこのような場合です。

f:id:primitivem:20170929141914p:plain

線型結合2層で繋いだ場合、数理モデルは以下のようになり、

Y = w_{2}h + b_{2}
h = w_{1}X + b_{1}

その微分は以下のようになります。

\frac{\partial}{\partial w_{2}} L = \frac{\partial Y}{\partial w_{2}} \frac{\partial L}{\partial Y} = h \cdot 2(Y-t)
\frac{\partial}{\partial w_{1}} L = \frac{\partial Y}{\partial w_{1}} \frac{\partial L}{\partial Y} = \frac{\partial h}{\partial w_{1}} \frac{\partial Y}{\partial h} \frac{\partial L}{\partial Y} = x \cdot w_{2} \cdot 2(Y-t)
\frac{\partial}{\partial b_{2}} L = \frac{\partial Y}{\partial b_{2}} \frac{\partial L}{\partial Y} = 1 \cdot 2(Y-t)
\frac{\partial}{\partial b_{1}} L = \frac{\partial h}{\partial b_{1}} \frac{\partial Y}{\partial h} \frac{\partial L}{\partial Y} = 1 \cdot w_{2} \cdot 2(Y-t)

例えば w_{1}=w_{2}=1,b_{1}=b_{2}=0,\eta=0.01,(x,t)=(1,2) の場合、

h=1 \cdot 1 +0
Y=1 \cdot 1 +0=1
\frac{\partial}{\partial w_{2}} = -2
\frac{\partial}{\partial w_{1}} = -2
\frac{\partial}{\partial b_{2}} = -2
\frac{\partial}{\partial b_{1}} = -2

このように数値が求まり、同じ計算をするスクリプトは以下で

model = Chain(l1=L.Linear(1,1),
                l2=L.Linear(1,1))
optimizer = optimizers.SGD(lr=0.01)
optimizer.setup(model)
x = Variable(np.array([[1]], dtype=np.float32))
t = Variable(np.array([[2]], dtype=np.float32))
model.l1.W.data.fill(1)
model.l1.b.data.fill(0)
model.l2.W.data.fill(1)
model.l2.b.data.fill(0)
optimizer.zero_grads()
h = model.l1(x)
y = model.l2(h)
loss = F.mean_squared_error(y, t)
loss.backward()
optimizer.update()
print model.l1.W.data
print model.l1.b.data
print model.l2.W.data
print model.l2.b.data
print model.l1.W.grad #この部分が手計算で求めた勾配
print model.l1.b.grad
print model.l2.W.grad #この部分が手計算で求めた勾配
print model.l2.b.grad
print loss.grad
print loss.data

この結果は手計算と同じ結果を出力することが確認できます。

線型結合2層にシグモイド関数

このように

f:id:primitivem:20170929142202p:plain

線型結合2層のネットワークに活性化関数としてシグモイド関数を使用した場合の数理モデル

Y = w_{2}h + b_{2}
h = \sigma( w_{1}X + b_{1} )

このようになり、その微分

{\sigma}'(x) = \left [  1-\sigma(x)  \right ]\sigma(x)
\frac{\partial}{\partial w_{2}} L = \frac{\partial Y}{\partial w_{2}} \frac{\partial L}{\partial Y} = h \cdot 2(Y-t)
z =w_{1}x +b
\frac{\partial h}{\partial w_{1}} = \frac{\partial z}{\partial w_{1}} \frac{\partial h}{\partial z} = x \cdot \left [ 1-\sigma(z) \right]\sigma(z)
\frac{\partial}{\partial w_{1}} L = \frac{\partial h}{\partial w_{1}} \frac{\partial Y}{\partial h} \frac{\partial L}{\partial Y} =x \cdot \left [ 1-\sigma(z) \right]\sigma(z) \cdot w_{2} \cdot2(Y-t)

このようになります。ちなみに面倒くさいので b_{1},b_{2} は省略しました。 例えば w_{1}=w_{2}=1,b_{1}=b_{2}=0,\eta=0.01,(x,t)=(1,2) の場合だと

z = 1 \cdot 1+0
\sigma(1)=0.73
h=0.73
Y=1 \cdot 0.73 +0=0.73
\frac{\partial}{\partial w_{2}} = -1.8542
\frac{\partial}{\partial w_{1}} = -0.500634

こういう計算結果になって、同じことを以下のスクリプトでも計算させると

model = Chain(l1=L.Linear(1,1),
                l2=L.Linear(1,1))
optimizer = optimizers.SGD(lr=0.01)
optimizer.setup(model)
x = Variable(np.array([[1]], dtype=np.float32))
t = Variable(np.array([[2]], dtype=np.float32))
model.l1.W.data.fill(1)
model.l1.b.data.fill(0)
model.l2.W.data.fill(1)
model.l2.b.data.fill(0)
optimizer.zero_grads()
h = F.sigmoid(model.l1(x))
y = model.l2(h)
loss = F.mean_squared_error(y, t)
loss.backward()
optimizer.update()
print model.l1.W.data
print model.l1.b.data
print model.l2.W.data
print model.l2.b.data
print model.l1.W.grad #この部分が手計算で求めた勾配
print model.l1.b.grad
print model.l2.W.grad #この部分が手計算で求めた勾配
print model.l2.b.grad
print loss.grad
print loss.data

手計算と同じ結果を得ることができます。

線型結合1層にガウス関数

次はこんな感じで

f:id:primitivem:20170929142515p:plain

線型結合1層にガウス関数を活性化関数として使用した場合の数理モデルは以下のようになり、

Y = w \phi(X) + b

またガウス関数は以下のかたちをしています。

\phi(x) = exp(-x^2)

この数式の微分は以下のようになります。

\frac{\partial}{\partial w} L = \frac{\partial Y}{\partial w} \frac{\partial L}{\partial Y} = 2\phi(X)(w\phi(X)+b-t)
\frac{\partial}{\partial b} L = \frac{\partial Y}{\partial b} \frac{\partial L}{\partial Y} = 2(w\phi(X)+b-t)

例えば w=1,b=0,\eta=0.01,(x,t)=(1,2) の場合だと

\phi(1) = e^{-1}
\frac{\partial}{\partial w} L = -1.200847
\frac{\partial}{\partial b} L = -3.2642

こういう結果になり、同じことをスクリプトで計算させると以下のようになります。

class Gaussian(Function):
    def forward_cpu(self, inputs):
        r = inputs[0]
        z = np.exp(-r*r)
        return z,

    def backward_cpu(self, inputs, grad_outputs):
        r = inputs[0]
        gz = grad_outputs[0]
        gr = -2.0*r*np.exp(-r*r) * gz
        return gr,

def gaussian(x):
    return Gaussian()(x)

model = L.Linear(1,1)
optimizer = optimizers.SGD(lr=0.01)
optimizer.setup(model)
x = Variable(np.array([[1]], dtype=np.float32))
t = Variable(np.array([[2]], dtype=np.float32))
model.W.data = np.array([[1]], dtype=np.float32)
model.b.data = np.array([0], dtype=np.float32)
optimizer.zero_grads()
h = gaussian(x)
y = model(h)
loss = F.mean_squared_error(y, t)
loss.backward()
optimizer.update()
print model.W.data #この部分が手計算で求めた勾配
print model.b.data
print model.W.grad #この部分が手計算で求めた勾配
print model.b.grad

これも結果は手計算と同じ結果が確かめられます。 このスクリプトではガウス関数を活性化関数として独自に定義しています。独自関数をchainerで使用するには
define your own function この公式ドキュメントを参考にしました。 またこスクリプトで書いた関数定義のクラスは以下に登場するスクリプトでは記述を省略します。

線型結合2層にガウス関数

つぎはこんなの

f:id:primitivem:20170929143253p:plain

数理モデルはこうなり

Y = w_{2} h + b_{2}
h = \phi(w_{1} X + b_{1})

微分は以下のようになります。

{\phi}'(x) = -2x \cdot e^{-x^{2}}
\frac{\partial}{\partial w_{2}} L = \frac{\partial Y}{\partial w_{2}} \frac{\partial L}{\partial Y} = h \cdot 2(Y-t)
z =w_{1}x +b
\frac{\partial h}{\partial w_{1}} = \frac{\partial z}{\partial w_{1}} \frac{\partial h}{\partial z} = x \cdot (-2ze^{-z^{2}})
\frac{\partial}{\partial w_{1}} L = \frac{\partial h}{\partial w_{1}} \frac{\partial Y}{\partial h} \frac{\partial L}{\partial Y} =-x \cdot \left [ 2ze^{-z^{2}} \right] \cdot w_{2} \cdot2(Y-t)

ここでも b_{1},b_{2} は割愛しました。 例えば w_{1}=w_{2}=1,b_{1}=b_{2}=0,\eta=0.01,(x,t)=(1,2) のとき、

z = 1 \cdot 1+0
h=e^{-1}
Y=1 \cdot e^{-1} +0=e^{-1}
\frac{\partial}{\partial w_{2}} = -1.20084
\frac{\partial}{\partial w_{1}} = 2.4016

数値はこのようになり、同じ計算のスクリプトは以下のようになり、

model = Chain(l1=L.Linear(1,1),
                l2=L.Linear(1,1))
optimizer = optimizers.SGD(lr=0.01)
optimizer.setup(model)
x = Variable(np.array([[1]], dtype=np.float32))
t = Variable(np.array([[2]], dtype=np.float32))
model.l1.W.data.fill(1)
model.l1.b.data.fill(0)
model.l2.W.data.fill(1)
model.l2.b.data.fill(0)
optimizer.zero_grads()
h = gaussian(model.l1(x))
y = model.l2(h)
loss = F.mean_squared_error(y, t)
loss.backward()
optimizer.update()
print model.l1.W.data
print model.l1.b.data
print model.l2.W.data
print model.l2.b.data
print model.l1.W.grad #この部分が手計算で求めた勾配
print model.l1.b.grad
print model.l2.W.grad #この部分が手計算で求めた勾配
print model.l2.b.grad
print loss.grad
print loss.data

手計算と同じ結果が確認できます。

ガウス関数のノードが2つ

長々といろんなパターンのネットワーク粗過程を見てきましたが、最後にこのようなモデル。

f:id:primitivem:20170929143544p:plain

中間層にノードが2つあります。数理モデル

Y = w_{21} h_{1} + w_{22} h_{2} + b_{21}
h_{1} = \phi(w_{11} X +b_{11})
h_{2} = \phi(w_{12} X +b_{12})
L=(Y-t)^2

こうなって、微分は以下のようになります。

\frac{\partial}{\partial w_{21}} L = \frac{\partial Y}{\partial w_{21}} \frac{\partial L}{\partial Y} =h_{1} \cdot 2(Y-t)
\frac{\partial}{\partial w_{22}} L = \frac{\partial Y}{\partial w_{22}} \frac{\partial L}{\partial Y} =h_{2} \cdot 2(Y-t)

z_{1}=w_{11}x+b_{11}
\frac{\partial h_{1}}{\partial w_{11}} =\frac{\partial z_{1}}{\partial w_{11}} \frac{\partial h_{1}}{\partial z_{1}} = x \cdot (-2z_{1}e^{-z_{1}^{2}})

\begin{split}\frac{\partial}{\partial w_{11}} L = \frac{\partial Y}{\partial w_{11}} \frac{\partial L}{\partial Y} \
=\frac{\partial h_{1}}{\partial w_{11}} \frac{\partial Y}{\partial h_{1}} \frac{\partial L}{\partial Y}\
=\frac{\partial h_{1}}{\partial w_{11}} \cdot w_{21} \cdot 2(Y-t) \
= -2 \cdot z_{1} e^{-z_{1}^{2}} \cdot w_{21} \cdot 2(Y-t) \end{split}

z_{2}=w_{12}x+b_{12}
 \frac{\partial h_{2}}{\partial w_{12}} =\frac{\partial z_{2}}{\partial w_{12}} \frac{\partial h_{2}}{\partial z_{2}} = x \cdot (-2z_{2}e^{-z_{2}^{2}})

\begin{split}\frac{\partial}{\partial w_{12}} L = \frac{\partial Y}{\partial w_{12}} \frac{\partial L}{\partial Y} \
=\frac{\partial h_{2}}{\partial w_{12}} \frac{\partial Y}{\partial h_{2}} \frac{\partial L}{\partial Y}\
=\frac{\partial h_{2}}{\partial w_{12}} \cdot w_{22} \cdot 2(Y-t) \
= -2 \cdot z_{2} e^{-z_{2}^{2}} \cdot w_{22} \cdot 2(Y-t) \end{split}

ここでも b_{11},b_{12},b_{21} は割愛しました。 そして例えば w_{11}=w_{12}=w_{21}=w_{22}=1,b_{11}=b_{12}=b_{21}=0,\eta=0.01,(x,t)=(1,2) のとき

z_{1} = 1 \cdot 1+0 =1
z_{2} = 1 \cdot 1+0 =1
h_{1}=\phi(z_{1})=e^{-1}
h_{2}=\phi(z_{2})=e^{-1}
Y=1 \cdot h_{1} + 1 \cdot h_{2} +0=2e^{-1}
\frac{\partial}{\partial w_{21}} L = -0.930176
\frac{\partial}{\partial w_{22}} L = -0.930176
\frac{\partial}{\partial w_{11}} L = 1.860353
\frac{\partial}{\partial w_{12}} L = 1.860353

このような手計算結果となり、同様の計算を行う以下のスクリプトでも

model = Chain(l1=L.Linear(1,2),
                l2=L.Linear(2,1))
optimizer = optimizers.SGD(lr=0.01)
optimizer.setup(model)
x = Variable(np.array([[1]], dtype=np.float32))
t = Variable(np.array([[2]], dtype=np.float32))
model.l1.W.data.fill(1)
model.l1.b.data.fill(0)
model.l2.W.data.fill(1)
model.l2.b.data.fill(0)
optimizer.zero_grads()
h = gaussian(model.l1(x))
y = model.l2(h)
loss = F.mean_squared_error(y, t)
loss.backward()
optimizer.update()
print model.l1.W.data
print model.l1.b.data
print model.l2.W.data
print model.l2.b.data
print model.l1.W.grad #この部分が手計算で求めた勾配
print model.l1.b.grad
print model.l2.W.grad #この部分が手計算で求めた勾配
print model.l2.b.grad
print loss.grad
print loss.data

手計算と同じ結果を確認することができました。

任意の関数にシグモイド関数で近似

上記のようなニューラルネットワークの最も基本的な計算の確認は終わりにして、次にシグモイド関数をつかったノードを複数個用意して適当な関数へのフィッティングを行ってみます。 スクリプトは以下のようになります。

def f(x):
    y = np.sin(x)
    return y

#教師信号の用意
buff = [[i*0.1] for i in range(-3*10,3*10)]
x_data = np.array(buff, dtype=np.float32)
buff = [[f(i*0.1)] for i in range(-3*10,3*10)]
t_data = np.array(buff, dtype=np.float32)
x = Variable(x_data)
t = Variable(t_data)
plt.plot(x.data, t.data, color="blue", linewidth=0.5, linestyle="-")

#フィッティング
model = Chain(l1=F.Linear( 1, 200),l2=F.Linear(200, 1))
optimizer = optimizers.SGD()
optimizer.setup(model)
for i in range(100000):
    if i%100 == 0:
        print i
    optimizer.zero_grads()
    h = F.sigmoid(model.l1(x))
    y = model.l2(h)
    loss = F.mean_squared_error(y, t)
    loss.backward()
    optimizer.update()

#プロット
x = np.arange(-7.0, 7.0, 0.1,dtype=np.float32)
buff = [[-7+i*0.1] for i in range(len(x))]
x_data = np.array(buff, dtype=np.float32)
x = Variable(x_data)
plt.axis([-7, 7, -2, 3])
h = F.sigmoid(model.l1(x))
y_plot = model.l2(h)
plt.plot(x.data, y_plot.data, color="red", linewidth=0.5, linestyle="-")
plt.show()

そして出力結果が以下

f:id:primitivem:20170929144033p:plain

青色の線がターゲットとなる任意の関数のを表しており、赤色の線がそれに近似した関数を表しています。 うまく関数近似ができているのがわかります。

任意の関数にガウス型関数で近似

つぎに先ほどと同じ任意の関数に対して、今度はガウス型関数で関数近似をしてみます。 これを行うにあたり、近似関数を構成する個々の各ガウス型関数の中心座標を任意の等間隔な位置に固定します。モデルの表現としては以下のようなものを考えます。

\begin{align*} y(\mathbf{x}) = \sum_{i=1}^{N} w_{i} \phi(|| \mathbf{x}-\mathbf{c}_{i} ||) \end{align*}

 \mathbf{c}_{i}ガウス型関数の中心座標の等間隔な位置を表すパラメータになっています。
しかし今回は入力ベクトルが一次元なので簡単のために

\begin{align*} y(x) = \sum_{i=1}^{N} w_{i} \phi(\left | x-c_{i} \right |) \end{align*}

このように表現した式を用いて以下のスクリプト関数近似を行ってみます。

#ガウス型関数の定義
class Gaussian(Function):
    def forward_cpu(self, inputs):
        r = inputs[0]
        z = np.exp(-r*r)
        return z,

    def backward_cpu(self, inputs, grad_outputs):
        r = inputs[0]
        gz = grad_outputs[0]
        gr = -2.0*r*np.exp(-r*r) * gz
        return gr,

def gaussian(x):
    return Gaussian()(x)

#RBF用のリンク関数の定義
class LinearFunction_RBF(Function):
    def forward(self, inputs):
        x, W = inputs
        return x.dot(W.T),

    def backward(self, inputs, grad_outputs):
        x, W = inputs
        gy, = grad_outputs

        gx = gy.dot(W)
        gW = gy.T.dot(x)
        return gx, gW

def linear_rbf(x, W):
    return LinearFunction_RBF()(x, W)

class Linear_RBF(Link):
    def __init__(self, in_size, out_size):
        super(Linear_RBF, self).__init__(W=(out_size, in_size))
        self.W.data.fill(1.0)

    def __call__(self, x):
        return linear_rbf(x, self.W)


"任意の関数にガウス型関数の中心座標を固定してフィッティング"

#任意の関数を用意
def f(x):
    y = np.sin(x)
    return y

#教師信号の用意
buff = [[i*0.1] for i in range(-3*10,3*10)]
x_data = np.array(buff, dtype=np.float32)
buff = [[f(i*0.1)] for i in range(-3*10,3*10)]
t_data = np.array(buff, dtype=np.float32)
x = Variable(x_data)
t = Variable(t_data)
plt.plot(x.data, t.data, color="blue", linewidth=0.9, linestyle="-")

#model = Chain(l1=L.Linear(10, 1))
model = Chain(l1=Linear_RBF(10, 1)) #これが重みだけを調節するリンク関数
optimizer = optimizers.SGD(lr=0.01)
optimizer.setup(model)

#隠れ層への入力データを生成する
n = 10
c = [-(n/2)*0.5+i*0.5 for i in range(n)]
#print c
def input(x,c):
    xbuff = []
    for d in x.data:
        buff = [math.fabs(d[0]-i) for i in c]
        #入力データが多次元ならばノルムをちゃんと求めないといけないが、今回は入力元なので
        #アドホックな書き方として絶対値にしている。ノルムはnp.linalg.norm()を使えい   
        xbuff.append(buff)
    return Variable(np.array(xbuff, dtype=np.float32))

x_in = input(x,c)
#print "input",x_in.data

#学習
for i in range(5000):
    optimizer.zero_grads()
    h = gaussian(x_in)
    #print h.data
    y = model.l1(h)
    #print "y.data  ",y.data
    loss = F.mean_squared_error(y, t)
    loss.backward()
    optimizer.update()

#プロット
plt.axis([-7, 7, -2, 2])

#近似関数のプロット
x_plot = np.arange(-7.0, 7.0, 0.1,dtype=np.float32)
x_ = []
for d in x_plot:
    buff = [d-i for i in c]
    x_.append(buff)
x_ = Variable(np.array(x_, dtype=np.float32))
#print x_.data

h = gaussian(x_)
y_plot = model.l1(h)
plt.plot(x_plot, y_plot.data, color="red", linewidth=0.9, linestyle="-")

#基底関数のプロット
x_ = Variable(x_plot)
w = model.l1.W.data[0]
for i in range(n):
    y = w[i]*gaussian(x_-c[i])
    plt.plot(x_plot, y.data, color="black", linewidth=0.2, linestyle="-")

#トレーニングデータをプロット
for i in range(len(x.data)):
    point_x = x.data[i][0]
    point_y = t.data[i][0]
    #plt.plot(point_x, point_y, "ro")

plt.show()

これを実行してみた結果は以下のように出力されます。

f:id:primitivem:20170929131332p:plain

青色の線がターゲットの任意の関数をプロットしていて、赤色の線が近似関数を表しています。これもうまくフィッティングできているのがわかります。 そして灰色で示された個々の基底関数が等間隔に並んでなめらかにターゲット関数にフィッティングしている様子がきれにに確認できました。

以上、ニューラルネットワークの基本的な計算過程の確認とガウス関数による関数近似についての記事でした。