デジタル・デザイン・ラボラトリーな日々

アラフィフプログラマーが数学と物理と英語を基礎からやり直す。https://qiita.com/yaju

ディープラーニング(深層学習)を理解してみる(識別問題 その２)

機械学習 tensorflow 人工知能

はじめに

3月の問題が解決しないまま２ヶ月が経ってしまい、4月分はうっかり記事を飛ばしてしまった。
このままだと5月分の記事も飛ばしてしまうので、あとで埋めるために5月分の記事を立てた。

しばしお待ち下さい。

ディープラーニング(深層学習)を理解してみる(識別問題)

機械学習 tensorflow 人工知能

はじめに

前々回の「対数logを理解してみる」と前回の「自然対数の底(ネイピア数) e を理解してみる」では、人工知能に使用する基礎的な数学知識が足りなかったのでシリーズとは脱線して書いてみました。 また、このシリーズで書いていきます。

ニューラルネットワークは、回帰問題と分類問題があります。
これまでやってきたリンゴとミカンの値段を予測するというのが回帰問題で、今回やるのは分類問題(識別問題)となります。

以前フレームワークのTensorFlowを使用してやりましたが、今回はフレームワークを使用しないでやってみます。
yaju3d.hatenablog.jp

識別問題とは

入力データから判断して分類させる問題。今回の例では買えるのか買えないのかを判断させる。 識別分類器としてロジスティック回帰を使う。

ロジスティック回帰とは

ロジスティック回帰は、発生確率を予測する手法です。 基本的な考え方は線形回帰分析と同じなのですが、予測結果が 0 から 1 の間を取るように数式やその前提に改良が加えられています。 今回のように購入有無（買える or 買えない）の2値しかとりえない値を従属変数の実績値として用い、説明変数を用いてその発生確率を説明するという構造になっています。

例題

たかしくんは八百屋へ財布を預かってお使いに行きました。
しかし、たかしくんはお金を 数えられません。

気まぐれおやじ曰く、
リンゴ２個＋ミカン３個、リンゴ０個＋ミカン１６個 なら買えるが、リンゴ３個＋ミカン１個、 リンゴ２個＋ミカン８個は買えないとのこと。

リンゴ１個＋ミカン１１個は買えますか？ → 識別問題 f:id:Yaju3D:20160421011949p:plain

式で表そうとしてみる…
f:id:Yaju3D:20160421012015p:plain f:id:Yaju3D:20160421012037p:plain

シグモイド曲線

f:id:Yaju3D:20160421012100p:plain
※シグモイド曲線とは、入力した値を0から1の間に収めてくれる関数の１つ
　多くの自然界に存在する事柄は、このようなＳ字曲線を取る。
　生物の神経細胞、細胞の生存率曲線などなど

予測式(モデル)が作れた!
f:id:Yaju3D:20160421012204p:plain

ロジスティック回帰
f:id:Yaju3D:20160421012219p:plain

プログラム

import numpy as np
import pandas as pd

def p_y_given_x(x, w, b):
    # x, w, b から y の予測値 (yhat) を計算
    def sigmoid(a):
        return 1.0 / (1.0 + np.exp(-a))
    return sigmoid(np.dot(x, w) + b)

def grad(x, y, w, b):
    # 現予測値から勾配を計算
    error = y - p_y_given_x(x, w, b)
    w_grad = -np.mean(x.T * error, axis=1)
    b_grad = -np.mean(error)
    return w_grad, b_grad
  
def gd(x, y, w, b, eta=0.1, num=1000):
    for i in range(1, num):
        # 入力をまとめて処理
        w_grad, b_grad = grad(x, y, w, b)
        w -= eta * w_grad
        b -= eta * b_grad
        e = np.mean(np.abs(y - p_y_given_x(x, w, b)))
        yield i, w, b, e

x = np.array([[2., 3.], [0., 16.], [3., 1.], [2., 8.]])
y = np.array([1., 1., 0., 0.])

# w, b の初期値を作成
w = [-10,-10]
b = 200
gen = gd(x, y, w, b)

for i, w, b, e in gen:
  if (i+1) % 100 == 0:
    print(i+1,b,w,e);

実行結果

100 196.2993226598431 [-22.375      -12.28583744] 0.6420233463035018
200 192.7195561228789 [-34.875      -12.06210203] 0.6420233463023709
300 189.1722948623358 [-47.25687076 -11.79079916] 0.48881302062549137
400 186.7724229463934 [-55.05210822 -11.5078109 ] 0.2665264782982885
500 185.57859640809747 [-58.91482296 -11.29546188] 0.02159297746827703
600 185.4754835638496 [-59.25851476 -11.22680836] 0.007859710540613895
700 185.42589243161035 [-59.4240343  -11.19266882] 0.004766951474845803
800 185.39307584996868 [-59.53360836 -11.16986379] 0.003414942999345553
900 185.36853902132037 [-59.61555193 -11.15273521] 0.002658633400582877
1000 185.34894411228476 [-59.68099871 -11.13901984] 0.002175879579138274

TensorFlowで実行した結果とほぼ同じになっています。
f:id:Yaju3D:20160421012400p:plain

予測結果

0.95848435 ≒ 0.96 f:id:Yaju3D:20160424232047p:plain

自然対数の底(ネイピア数) e を理解してみる

機械学習

はじめに

前回、対数を理解してみました。 yaju3d.hatenablog.jp

今回は機械学習を学ぶ上で出てくる自然対数の底(ネイピア数 e)とは何かを理解していきます。
間違える人は結構多いですが、e は「自然対数」ではありません。「ネイピア数」あるいは「自然対数の底」と呼ばれる定数です。

自然対数の底とは

対数の底には、10進法を使う常用対数とネイピア数と呼ばれる数学定数を使う自然対数があります。
記号として通常は e が用いられ、その値は、e = 2.71828 18284 59045 23536 02874 71352\cdots と続く無理数で超越数となっています。 とても不思議な値ですよね。
不思議な値といえば円周率 \pi = 3.14159265359\cdots も同様に無理数で超越数となっています。

qiita.com

最後に

今回はQiitaの方に書いてしまったので、リンクするだけにします。

対数logを理解してみる

機械学習

はじめに

機械学習を学んでいると対数logがでてくる。基礎的なことから対数を理解してみたい。
指数はイメージし易いが、対数は分かりにくいと思われている。指数と対数はペアの関係にあり、かけ算とわり算のように逆関係にある。 先ずは、指数の大きさを視覚的にイメージするために、アメリカ・ワシントンにある航空宇宙博物館で公開されていた9分半の映画「パワーズ・オブ・テン」を紹介する。

9分半の映画

パワーズ・オブ・テンです、10の冪(10^{n})の違いを視覚でご覧ください。


Powers of Ten with Japanese translation

対数とは

例えば、2を3回かけ算すると 2 \times 2 \times 2 = 8 になります。これを「2を3乗したら8になる」と言い、以下のように書きます。

2 ^ 3 = 8

このとき、2 の右上に乗っている 3 のことを「指数」と言います。指数は「1つの数を何回掛けるか」を表しています。
一方、「◯を何乗すれば△になるか」を表す数のことを「対数」と言います。 例えば「2を何乗すれば8になるか」を表す数は以下のように表記され、これを「2を底とする8の対数」と言います。

log_{2}8 = 3

2 を底をする 8 の対数は 3 ということになります。
もう少し身近な値としてみます。1000 は、10 の 3乗 で「10を底とする1000の対数」は以下のようになります。

log_{10}1000 = 3

計算の簡略化

対数には大きな桁数でも簡単に計算できるようになるというメリットがある。
それには指数法則と対数法則を知っておく必要がある。

指数法則①

a^{p} \times a^{q} = a^{(p+q)}

3^{4} \times 3^{5} について考える。3^{4} \times 3^{5} を指数を使わない形で表現すると(3 \times 3 \times 3 \times 3) と(3 \times 3 \times 3 \times 3 \times 3) となる。 3 を繰り返しかけ算する回数は、 4 + 5 = 9 回である。
よって、3^{4} \times 3^{5} = 3^{(4+5)} となる。

指数法則②

(a^{p})^{q} = a^{(p \times q)}

(5^{3})^{4} について考える。(5^{3})^{4} は、5^{3} を4回繰り返しかけ合わせるので、(5^{3})^{4} = 5^{3} \times 5^{3} \times 5^{3} \times 5^{3} である。指数の部分の合計は、3乗が4回、つまり (3 \times 4) によって計算できる。
よって、 (5^{3})^{4} = 5^{3 \times 4} となる。

指数法則③

(a \times b)^{p} = a^{p} \times b^{p}

(3 \times 7)^{5} について考える。(3 \times 7)^{5} は、(3 \times 7) を5回繰り返しかけ合わせるので、 (3 \times 7)^{5} = (3 \times 7) \times (3 \times 7) \times (3 \times 7) \times (3 \times 7) \times (3 \times 7)
= (3 \times 3 \times 3 \times 3 \times 3) \times (7 \times 7 \times 7 \times 7 \times 7) 、つまり = 3^{5} \times 7^{5} によって計算できる。
よって、 (3 \times 7)^{5} = 3^{5} \times 7^{5} となる。

対数法則①

かけ算を足し算に変換

log_{a}(M \times N) = log_{a}M + log_{a}N

log_{2}(8 \times 4) = log_{2}32 = 5 = 3 + 2 = log_{2}8 + log_{2}4
注目すべきは、M \times N というかけ算が、log_{a}M + log_{a}N という足し算に変換されているということである。このことが対数を利用してかけ算を足し算へ簡略化して計算する際に生かされる。

対数法則②

わり算を引き算に変換

log_{a}(M \div N) = log_{a}M - log_{a}N

log_{2}(8 \div 4) = log_{2}2 = 1 = 3 - 2 = log_{2}8 - log_{2}4
今回は、M \div N というわり算が、log_{a}M - log_{a}N という引き算の形に変換される。

対数法則③

累乗を簡単なかけ算に変換

log_{a}M^{k} = k \times log_{a}M

log_{2}8^{3} = log_{2}512 = 9 = 3 \times 3 = 3log_{2}8
例えば、2^{29} のような大きな計算も、この対数法則③を使えば一瞬にして近似値を求められる。

底の変換公式

底も自由に変えられるようになります。
\displaystyle \log_{a}b = \frac{(\log_{c}b)}{(\log_{c}a)} = 対数

\displaystyle \log_{2}8 = \frac{(\log_{10}8)}{(\log_{10}2)}
\displaystyle = \frac{0.903089987}{0.301029995} = 3

C# 数学１０ 「常用対数、ネイピア数e-基本１、自然対数-基本１」e log ln exp

対数の起源

対数を考え出したのはジョン・ネイピアです。
対数(logarithm)の名前の由来は、logos (比、神の言葉)とギリシャ語のarithmos (数) を合わせて logarithms(ロガリズム) という造語でネイピアが考案しました。

ジョン・ネイピア

1550年、宗教戦争が続くスコットランドにネイピアは生まれました。時代はヨーロッパ列強諸国が覇権を争う大航海時代のまっただ中。ネイピアはマーキストン城の城主として官の仕事の他、領民のために農業土木、軍事技術など多くの発明を行うエンジニアとして活躍しました。(中略) 遠洋航海に必要な天文学（船の測位）を支える数学が球面三角法です。数学者でも天文学者でもないネイピアは現実の切実な問題解決を目標としていました。それが天文学的計算の克服です。三角関数の計算の中に現れる大きい数の計算は天文学者を苦しめました。大航海時代は計算との闘いでもあったのです。天文学者は直面する天文学的計算を克服する手立てを見つけることができませんでした。彼らの計算を助けるために、ネイピアはついに新しい計算法を見つけ出す決心をします。時にネイピア44歳、1594年。その20年後の1614年、ついに人類は青天の霹靂として「対数」を手にします。第5回　ジョン・ネイピア対数誕生物語

また、日々私たちが使っている小数点(.)を使う表記はネイピアの発明です。小数点はネイピアが対数を生み出す過程で考え出した副産物だったのです。第6回　ジョン・ネイピア小数点「・」誕生物語
※小数の考え方はネイピア(1550〜1617)とほぼ同時期のシモン・ステヴィン(1548-1620) ですが19.178と表す小数を「19⓪1①7②8」のように表記していた。ステヴィンの本は1585年に出版、これはネイピアが対数表の計算を開始する1594年の9年前のことです。ネイピアはステヴィンの小数は使いづらいと判断し、もっと機能的で使いやすい小数の表記法を求め続けたのです。

発明を行うエンジニア

ネイピアはマーキストン城の城主であり、困った人を助けるために自分の才能を使う、優秀なエンジニアだったのです。
自分の領地の収穫を増やすために肥料や揚水機の研究をしたり、スペインの.侵攻を恐れて潜水艦や戦車などの兵器も数多く開発しています。これも領地内の人民を安心させるためだったのでしょう。
また、対数の計算を効率化するためにネイピアの骨と呼ばれる、かけ算や割り算などを簡単に行うための計算器具を発明しています。これは1617年(ネイピアが亡くなる年)に論文「小さな棒による計算術」としてエディンバラで発表されました。

対数表

対数表には、指定した数を底とする対数(例 底=2)、常用対数(底=10)、自然対数(底e=2.718281828…)があります。試験などで底を省略した対数表が出た場合は常用対数となります。

ネイピアが作成した対数表は不思議な底(0.9999999)だったため、使いやすい底を10にした常用対数表をネイピアの意思をついだブリッグスが作成しました。1617年に1000まで計算して出版し、1624年には1から20000までと90000から100000まで、小数点以下14桁まで計算した対数表を出版しました。1628年にブリッグスの対数表で抜けていた20,000～90,000の対数表をオランダのアドリアン・ブラックが完成させました。ただ、ブラックの対数表は10桁のものでしたが実際に利用するには十分な桁数でした。対数の発見

qiita.com

対数をとるメリット

対数のメリットは、極めてケタ数の大きい数の面倒なかけ算や割り算を、易しい足し算や引き算に直してくれることにあります。但し、求めた値はあくまで近似値であり、正確な値ではありません。
電卓や計算機のなかった昔の時代において非常に便利な計算方法でした。ピエール・シモン・ラプラスからは「天文学者の寿命を2倍にした(天文学者が一生にこなせる計算量を2倍にした)」と絶賛されるようになりました。

対数表を使った計算

例 131 \times 219 \times 563 \times 608 を計算する場合

10進数を底とする対数で表すと、log_{10}(131 \times 219 \times 563 \times 608) となる。
log_{10}( (1.31 \times 10^{2}) \times (2.19 \times 10^{2}) \times (5.63 \times 10^{2}) \times (6.08 \times 10^{2}) )
= log_{10}((1.31 \times 2.19 \times 5.63 \times 6.08) \times 10^{8})
ここから、かけ算が足し算になります。
= log_{10}1.31 + log_{10}2.19 + log_{10}5.63 + log_{10}6.08 + log_{10}10^{8}
= log_{10}1.31 + log_{10}2.19 + log_{10}5.63 + log_{10}6.08 + 8

ここで、 log_{10}1.31 、log_{10}2.19 、log_{10}5.63 、log_{10}6.08 の値を対数表から読み取る。
log_{10}(131 \times 219 \times 563 \times 608) \fallingdotseq 0.1173 + 0.3404 + 0.7505 + 0.7839 + 8 = 1.9921 + 8 = 0.9921 + 9 ここで、0.9921 に近い値を対数表から読み取ると、0.9921 \fallingdotseq log_{10}9.82 となる。
log_{10}(131 \times 219 \times 563 \times 608) \fallingdotseq log_{10}9.82 + 9 = log_{10}9.82 + log_{10}10^{9} = log_{10}(9.82 \times 10^{9})
こうして得られた log_{10}(131 \times 219 \times 563 \times 608) \fallingdotseq log_{10}(9.82 \times 10^{9}) の両辺を見比べると
131 \times 219 \times 563 \times 608 \fallingdotseq 9.82 \times 10^{9} = 9820000000 (実際の値は、 9820359456)

例 2^{29} を計算する場合

対数法則③を使います。
2^{29}は、10を底とする対数で表すと log_{10}2^{29} = 29 \times log_{10}2 となる。
log_{10}2 の値は、常用対数表から 0.3010 である。よって、log_{10}2^{29} = 29 \times 0.3010 = 8.7290 = 0.7290 + 8
常用対数の値が、0.7290 に近い真数の値を表から読み取ると、0.7290 \fallingdotseq log_{10}5.36 となる。
log_{10}2^{29} \fallingdotseq log_{10}5.36 + 8 \times log_{10}10 = log_{10}5.36 + log_{10}10^{8} = log_{10}(5.36 + 10^{8})
よって、2^{29}は、約 5.36 \times 10^{8} = 536000000 と計算できる。(実際の値は、536870912)

基本情報技術者試験の問題

ソートの計算量

計算量は一般に、O記法（オー記法、もしくはビッグオー記法と呼ぶ）を使って表します。O記法では計算量を、O(n) や O(n^{2}) のような形で記述します。
性能が良い : O(1) \lt O(\log n) \lt O(n) \lt O(n\log n) \lt O(n^{2}) \lt O(n^{3}) : 性能が悪い

  • O(1)
    データ数nにかかわらず一定なので、nが2倍、3倍、10倍となっても計算量は等しいままです。
  • O(\log n)
    計算量を求める場合は底を2として考えます。
    nが2倍、3倍と増えても、計算量は1、2と増える程度なのでデータ量が多い場合に有用です。
  • O(n)
    比例なので、nが2倍、3倍と増えると、計算量も同様に2倍、3倍となります。
  • O(n \log n)
    先述のO(\log n) にnの計算量がかけられたものです。
    例えば、nが2倍になると計算量は約2.5倍、nが10倍になると計算量は約20倍、nが100倍になると計算量は約300倍になります。
  • O(n^{2})
    n個の要素を用いた総当りの組み合わせを求める場合がこのパターンに当てはまります。
    二次関数なので、nが2倍、3倍と増えると、計算量は4倍、9倍と増えます。

O記法にはいくつかのルールがある。

  • 中身は一番大きな規模だけ残す(最大次数の項以外を除く) 例 2+n+2n^{2}→2n^{2}
  • 係数は除く 例 2n^{2}→n^{2}
  • 係数は1にする

例 1 ~ Nまでの整数の総和を求めよ

  • 普通に計算する → N-1回足し算 → O(N-1) → O(N) 　
  • 公式 \displaystyle S=\frac{(n+1)n}{2} を使う → 足し算と掛け算と割り算 → O(3) → O(1)

その中で対数が使われるのが下記となります。

O記法 概要 使用例
O(logN) <対数時間>
処理をひとつ行うたびに処理対象を何割か減らせるようなアルゴリズム。データ量が増えても計算時間がほとんど増えない。多くの場合、これ以上の改善はする必要はない。		 ソート済み配列の二分探索
O(NlogN) <準線形、線形対数時間>
ちょっと重いO(N)程度。マージソートのように二分木でデータを分割し、かつそれらをリニアサーチするようなアルゴリズムの計算量。二分木のオーダー(logN)×リニアサーチのオーダー(N)をかけあわせてNlogNになる。		 クイックソート、マージソート、ヒープソート

二分探索の計算量

ソートされている前提で、中央の値より小さい(左側)か大きい(右側)が分かるので、次からは半分は探索する必要はなくなる。
これを2の指数で表した場合、指数の数が1つずつ減ることになるので対数にすると、O(\log n) になる。
f:id:Yaju3D:20180122012235p:plain

クイックソートの計算量

分割統治法に基づくアルゴリズムです。
f:id:Yaju3D:20180122013630p:plain
毎回半分ずつにデータが分かれたと仮定すると、平均計算量は n個 \times \log_{2} n 段 → O(n\log n) となる。
f:id:Yaju3D:20180122213121p:plain

n進数での最大値

14けたの16進数の最大値は、10進数で表すと何けたか。ここで、 log_{10}2 ＝ 0.301 とする。

log_{10}(16^{14})
= 14 \times log_{10}(16)
= 14 \times log_{10}(2^{4})
= 14 \times 4 \times log_{10}(2)
= 14 \times 4 \times 0.301 = 16.856

繰り上げて答えは17桁になります。

身近な問題

地震の大きさや音階や星の等数や化石の年代測定や複利計算にも対数が使用します。他にも電気通信や音の単位に用いられるデシベル、音の大きさを表すホンでも使用されます。

地震のエネルギーの大きさ

地震のエネルギーの大きさを表す言葉として「マグネチュード」と言います。

マグニチュードと震度の違いは？
「マグニチュード」は、地震そのものの大きさ（規模）を表すものさしです。一方「震度」は、ある大きさの地震が起きた時のわたしたちが生活している場所での揺れの強さのことを表します。

マグネチュードを M、地震が発生するエネルギーを E (単位モジュール)とすると、logE = 4.8 + 1.5M という関係式があります。
この式より、E = 10^{4.8 + 1.5M} = 10^{4.8}10^{1.5M} となります。
したがって、 M が1増えるとエネルギーは10^{1.5}倍(およそ32倍)、 M が2増えるとエネルギーは10^{3}倍(およそ1000倍) になります。また、M が0.2増えるとエネルギーは10^{0.3}倍(およそ2倍) になります。

マグネチュードが2違うと1000倍で1違うと32倍って差が分かりにくいですが、x \times x = 1000 とすると x^{2} = 1000 = x = \sqrt{1000} = 31.6227 \fallingdotseq 32 となるわけです。

音階の違い

ピタゴラスの定理などで有名な古代ギリシャの数学者ピタゴラスが音階「ドレミファソラシド」を作りました。ピタゴラスが発見した音階は弦の長さで決められたものでしたが、我々が通常聞いている音階は振動数の比で決められています。

平均律音階は、1オクターブ高い音は振動数が2倍であり、その1オクターブを12音階に均等に分けたものである。
f:id:Yaju3D:20180128121348p:plain

音階 ド ド# レ レ# ミ ファ ファ# ソ ソ# ラ ラ# シ ド
周波数比 2^{0}=1 2^{\frac{1}{12}} 2^{\frac{2}{12}} 2^{\frac{3}{12}} 2^{\frac{4}{12}} 2^{\frac{5}{12}} 2^{\frac{6}{12}} 2^{\frac{7}{12}} 2^{\frac{8}{12}} 2^{\frac{9}{12}} 2^{\frac{10}{12}} 2^{\frac{11}{12}} 2^{\frac{12}{12}}=2
周波数 fr fr^{1} fr^{2} fr^{3} fr^{4} fr^{5} fr^{6} fr^{7} fr^{8} fr^{9} fr^{10} fr^{11} fr^{12}
近似値 262 277 294 311 330 349 370 392 415 440 466 494 523

※周波数 f=262

r^{12} = 2 で、これを満たす r を求めると 1.06^{12} \fallingdotseq 2 となります。半音上がることに振動数を 1.06倍されている。
\displaystyle r = \sqrt[12]{2} = 2^{\frac{1}{12}} = 1.0594631 \fallingdotseq 1.06

つまり、音階も対数変換(底1.06)ととらえることができる。\log_{1.06}(周波数)＝音階

星の等数

歴史上で星の明るさをランクづけしたのは、古代ギリシアの天文学者ヒッパルコスです。彼は、肉眼で見えるもっとも明るい星を1等星、かろうじて見える星を6等星として、1等星～6等星までの6段階に分けました。月日が経ち、1856年にイギリスの天文学者のノーマン・ロバート・ポグソンは、ジョン・ハーシェルの研究結果を元にし1 等星の明るさは 6 等星の 100 倍であり、1 等級ごとの明るさの違いは 100^{\frac{1}{5}} 倍であると定義しました。

f:id:Yaju3D:20180128204643p:plain

\displaystyle r = \sqrt[5]{100} = 100^{\frac{1}{5}} = 2.51188643151 \fallingdotseq 2.51

つまり、星の等数は対数(底2.51)となります。

化石の年代測定

化石の年代測定ってどうやっているのか不思議でした。
これは「放射性炭素年代測定法」というのが1974年にアメリカのリビーという人が発見した方法で、炭素14の性質を利用しています。
炭素14は放射性炭素ともいわれ、電子を放出して窒素14に変わる。この崩壊によって炭素14の数が半分になるまでの期間を半減期といい5730年となっている。
生きている生物(動物、植物)はこの炭素14を体内に取り込むので、体内の炭素14の割合は大気中の割合と同じとなる。また生物が死ぬと炭素14の供給がなくなり崩壊だけが続くので発見した資料の炭素14の割合を調べることで、その資料が死んでからの年数が推定できる。

ある木管に炭素14の割合を調べたら、75%に減っていた。
この場合、炭素14が1年で r 倍に現象するとして、この木簡が x 年前のものだとすると、
r^{x}=0.75 また r^{5730}=0.5

① x\log r=\log 0.75 、② 5730\log r=\log 0.5
① と ② より
\displaystyle x=\frac{\log0.75}{\log r} = \frac{5730}{\log 0.5} \times \log 0.75
\displaystyle =\frac{5730(\log3 - 2\log2)}{- \log 2}
= 5730 \times 0.4150 = 2378 年前

複利計算

サラリーローンで、1万円を1日1割の複利で2ヶ月借りる場合、
＜複利計算＞元利合計 = 元金 \times (1 + 利率)^{期間の数} の式となります。
1 \times (1+0.1)^{60} \fallingdotseq 305 なんと、約305万円にもなります。

これを対数で計算してみます。電卓がない状態で1.1を60回かけるんなんてうんざりですが、その場合は対数表を使えばいい。
\log_{10}x = \log_{10}1.1^{60}
対数の法則を使う
\log_{10}x = 60\log_{10}1.1
対数表では\log_{10}1.1 = 0.0414 となります。
\log_{10}x = 60 \times 0.414 = 2.484
このうち、小数部の 0.484 を対数表で見ると、3.05 です。
\log_{10}x = 2.484 = 2 + 0.484
= \log_{10}10^{2} + \log_{10}3.05
= \log_{10}100 \times 3.05
= \log_{10}305
x = 305

自然対数

自然対数は、2.71828182845\cdots という無理数を底にした対数となります。一方、常用対数は底を10にした対数となります。
常用対数は電卓やコンピューターの登場により使用されることがなくなりましたが、自然対数はコンピューターで使用する上で標準と底となっています。これは次記事に書きます。 yaju3d.hatenablog.jp

機械学習での役割

  • かけ算を足し算または割り算を引き算にすることで計算を楽にさせる。
    \displaystyle \prod_{n=1}^{N} f(n) = \sum_{n=1}^{N}\log f(n) ※\prodは、かけ算を繰り返す記号
  • 機械学習では非常に大きな数を扱います。プログラムでこのような数を扱うと、オーバーフローを起こしてしまうときがあります。そのような数は、対数 をとって扱うことで、オーバーフローを防ぐことができます。例えば、100000000 = 10^{8}、0.000000001 = 10^{-8} となります。
  • 機械学習では確率を使いますが、同時確率の場合には 1 以下の数の掛け算の連続になります。多い場合にはアンダーフローを起こしてしまうときがあります。そのような数は、対数 をとって扱うことで、アンダーフローを防ぐことができます。
    例としてサイコロを2回投げた場合、1回目に1の目が出て2回目に2の目が出る確率は、下記の計算となる。
    \displaystyle \frac{1}{6} \times \frac{1}{6} = \frac{1}{36}
  • 対数は単調増加関数であるため、ある関数 f(x) があって f(x) を最小にする値を求めた場合、対数をとった \log f(x) でも最小の値は同じになります。また、最大値を求める場合でも同じです。
    f(x)=(x-1)^{2}+2は、x=1の時、最小値をとる
    \log f(x)=\log ((x-1)^{2}+2) も、x=1の時、最小値をとる

単調増加関数

単調増加な関数は、任意な値 x1 と x2 があったとき 、 x1 \lt x2 ならば f(x1) \lt f(x2) となるような関数 f(x) であるいうことです。
例えば、 x1=3、x2 = 5 は x1 \lt x2 が成り立ちます。この関係性は対数に直しても \log x1 \lt \log x2 が成り立ちます。

「サザエさんのじゃんけん データ分析」の2017年の結果

機械学習

はじめに

明けましておめでとうございます。

昨年も東芝に暗雲が立ち込め、東証2部に降格して日経平均から外れ、ついには18年3月末にサザエさんからのスポンサー契約を降板する方向で調整とあいなりました。
さて、サザエさんのじゃんけんを長年記録しデータをWebスクレイピングさせて頂いた「サザエさんジャンケン学」が2017年6月25日を以て終了されていました。データを取得した際に7月以降がなかったので、あれっと思って見たらそういうことでした。長い間お疲れ様でしたm(_ _)m
www.asahi-net.or.jp

さてさて、2017年のサザエさんのじゃんけん結果はどうなったでしょう。 ちなみに、2016年のサザエさんのじゃんけん結果は、27勝11敗12分(勝率0.711)でした。

人工知能による予測化を断念

本来はDeepLearningを使用して予想をする予定でしたが、まだ勉強中のままです。あと少しだと思うので今年中にはなんとかなると思います。
昨年、R言語からPythonに切り替えました。
qiita.com

Python使うに際にはいつも Jyupter notebook だったので、この為に久しぶりに Visual Studio Code で Python を実行しようとしてバージョンアップ(1.19.1)したらタスク実行の設定がいろいろ変わって動かなくて、四苦八苦しながらなんとか動かせるようにしました。下記の記事を修正してあります。 yaju3d.hatenablog.jp

次の手の予測アルゴリズム

  • チョキが多いので、グー ＞ チョキ ＞ パーの優先順位とする
  • 前回と違う手を出すので、上記の優先順位で勝手を選ぶ
  • 二手前と一手前が違う手なら、残りの手を出すので勝手を選ぶ
  • 三手の中に同手がある場合、 残りの手を出すので勝手を選ぶ
  • 二手前と一手前が同じ手なら、勝手を出すので負手を選ぶ

2017年の勝敗結果

年月 サザエさんの手 予想の手 勝敗結果
01月01日 休み
01月08日 チョキ グー 勝ち
01月15日 チョキ グー 勝ち
01月22日 パー パー 引き分け
01月29日 グー パー 勝ち
02月05日 グー グー 引き分け
02月12日 パー グー 負け
02月19日 チョキ グー 勝ち
02月26日 チョキ パー 負け
03月05日 グー パー 勝ち
03月12日 パー チョキ 勝ち
03月19日 休み
03月26日 チョキ グー 勝ち
04月02日 チョキ パー 負け
04月09日 グー パー 勝ち
04月16日 パー チョキ 勝ち
04月23日 パー グー 負け
04月30日 チョキ グー 勝ち
05月07日 グー パー 勝ち
05月14日 パー チョキ 勝ち
05月21日 チョキ グー 勝ち
05月28日 グー パー 勝ち
06月04日 チョキ チョキ 引き分け
06月11日 パー チョキ 勝ち
06月18日 グー パー 勝ち
06月25日 グー グー 引き分け
07月02日 チョキ グー 勝ち
07月09日 パー チョキ 勝ち
07月16日 グー パー 勝ち
07月23日 パー グー 負け
07月30日 チョキ グー 勝ち
08月06日 パー パー 引き分け
08月13日 グー パー 勝ち
08月20日 チョキ グー 勝ち
08月27日 グー チョキ 負け
09月03日 パー チョキ 勝ち
09月10日 パー グー 負け
09月17日 チョキ グー 勝ち
09月24日 グー パー 勝ち
10月01日 チョキ チョキ 引き分け
10月15日 パー チョキ 勝ち
10月22日 グー パー 勝ち
10月29日 休み
11月05日 チョキ グー 勝ち
11月12日 チョキ チョキ 引き分け
11月19日 グー チョキ 負け
11月26日 パー チョキ 勝ち
12月03日 パー グー 負け
12月10日 チョキ グー 勝ち
12月17日 グー パー 勝ち
12月24日 パー チョキ 勝ち

結果は、32勝9敗7分(勝率0.780)となりました。

ちなみに、サザエさんじゃんけん研究所 公式ウェブサイトの サザエさんの手の予想と勝負結果(2017年)が29勝8敗11分(勝率0.783)でした。

今回は勝数では上回ったのですが、勝率では僅差で負けました。 勝率の計算は、「勝ち / (勝ち + 負け)」で行っているのですが、負けの1つが響いたわけです。

【2017/01/09追記】
２０１７年冬版 サザエさんじゃんけん白書によるとクール(四半期)の初回(1月、4月、7月、10月の初回)はチョキが出やすいとのことです。もし、これを取り入れた場合、04月02日(負け)と10月01日(引き分け)なので、34勝8敗6分(勝率0.809)になるかな。

年 データ分析 研究所公式
2013 24勝13敗12分(勝率0.649) 25勝9敗17分(勝率0.735)
2014 30勝10敗11分(勝率0.750) 30勝9敗12分(勝率0.769)
2015 32勝9敗9分(勝率0.780) 33勝9敗8分(勝率0.785)
2016 27勝11敗12分(勝率0.711) 22勝13敗15分(勝率0.628)
2017 32勝9敗7分(勝率0.780) 29勝8敗11分(勝率0.783)

スライド

2013年に静岡Developers勉強会で機械学習を学び、2014年1月にネタとしてSlideShareに公開しました。

最後に

今年中にはTensorFlowを使って人工知能のディープラーニングで予測手を作りたいと思います。まだ組むだけの理解度(n-gramモデル)が足りないためですが、あと少しで理解度が進むはずです。

【2018/01/07 追記】
r_stdさんが、機械学習を用いて検証してくれました。2016年に統計検定準1級に合格している方で幾つかの手法を使っています。その中で、naive bayesで33勝9敗6分と好成績を残しています。タイトルが①なので次回も期待しています。
r-std.hatenablog.com

ディープラーニング(深層学習)を理解してみる(勾配降下法：計算確認)

機械学習 tensorflow 人工知能

はじめに

前回の続きです。 yaju3d.hatenablog.jp

前回の計算が本当に合っているのか、Pythonを使って実証してみたいと思います。

プログラム

重みの更新

①現在の重みで推測値を求める

import numpy as np
a = np.array([10, 20])
b = np.array([[1,3,5],[3,1,7]])
print(np.dot(a,b))

[ 70, 50, 190]

全体の誤差を求める。
二乗誤差は、mean_squared_error の名前です。

def mean_squared_error(y, t):
    return 0.5 * np.sum((y - t)**2)

y = np.array([70, 50, 190])
t = np.array([190, 330, 660])

print(mean_squared_error(y, t))

156850.0

②勾配\Delta Eを計算する
勾配は、gradient の名前です。転置は x.T とTを付けるだけなんですね。

def gradient(y, t, x):
    return np.dot((y - t),x.T)
    
t = np.array([190, 330, 660])
y = np.array([70, 50, 190])
x = np.array([[1,3,5],[3,1,7]])

print(gradient(y, t, x))

[-3310 -3930]

③重みを更新する
②で求めた勾配\Delta E を用いて、現在の重みw を更新します。学習係数η は 0.02 としています。

def weight(w, lr, gd):
    return w - lr * gd

t = np.array([190, 330, 660])
y = np.array([70, 50, 190])
x = np.array([[1,3,5],[3,1,7]])
w = np.array([10, 20])

print(weight(w, 0.02, gradient(y, t, x)))

[ 76.2 98.6]

2回目の①現在の重みで推測値を求める

a = weight(w, 0.02, gradient(y, t, x))
b = np.array([[1,3,5],[3,1,7]])
print(np.dot(a,b))

[ 372. 327.2 1071.2]

全体の誤差を求める。
二乗誤差は、mean_squared_error の名前です。

y = np.dot(a,b)
t = np.array([190, 330, 660])

print(mean_squared_error(y, t))

101108.64

というのを繰り返す。

import numpy as np

# 二乗誤差を求める
def mean_squared_error(y, t):
    return 0.5 * np.sum((y - t)**2)

# 勾配を求める
def gradient(y, t, x):
    return np.dot((y - t),x.T)

# 重みを求める
def weight(w, lr, gd):
    return w - lr * gd

# 教師データ
t = np.array([190, 330, 660])
# 重み
w = np.array([10, 20])
# 入力データ(個数)
x = np.array([[1,3,5],[3,1,7]])
y = None
for i in range(100):
    if y is not None:
        # 勾配を求める
        grad = gradient(y, t, x)
        # 重みを求める
        w = weight(w, 0.02, grad)
    # 推測値を求める
    y = np.dot(w, x)
    # 誤差を求める
    l = mean_squared_error(y, t)
    if (i + 1) % 10 == 0:
        print("step=%3d, a1=%6.2f, a2=%6.2f, loss=%.2f" % (i + 1, w[0], w[1], l))
print("Estimated: a1=%6.2f, a2=%6.2f" % (w[0], w[1]))

step= 10, a1= 78.84, a2= 48.86, loss=4531.39
step= 20, a1= 89.41, a2= 32.85, loss=564.82
step= 30, a1= 94.92, a2= 27.91, loss=264.21
step= 40, a1= 97.30, a2= 26.05, loss=217.63
step= 50, a1= 98.28, a2= 25.30, loss=209.89
step= 60, a1= 98.68, a2= 25.00, loss=208.59
step= 70, a1= 98.84, a2= 24.88, loss=208.38
step= 80, a1= 98.91, a2= 24.83, loss=208.34
step= 90, a1= 98.94, a2= 24.81, loss=208.33
step=100, a1= 98.95, a2= 24.80, loss=208.33
Estimated: a1= 98.95, a2= 24.80

実行結果
リンゴ 98.95 ≒ 99、ミカン 24.80 ≒ 25 → 99 x 2 + 25 x 4 ≒ 297

以前、TensorFlowで組んだ際の結果とほぼ同じになったのでこれでいいかな。
リンゴ 98.81 ≒ 99、ミカン 24.90 ≒ 25 → 99 x 2 + 25 x 4 ≒ 297
yaju3d.hatenablog.jp

上記サイトでTensorFlowで組んだのに近付けるために、train_x と train_y を縦ベクトル に変更してみました。

f:id:Yaju3D:20160419005112p:plain f:id:Yaju3D:20160419002003p:plain

import numpy as np

# 二乗誤差を求める
def mean_squared_error(y, t):
    return 0.5 * np.sum((y - t)**2)

# 勾配を求める
def gradient(y, t, x):
    return np.dot((y.T - t.T), x)

# 重みを求める
def weight(w, lr, gd):
    return w - lr * gd

train_x = np.array([[1., 3.], [3., 1.], [5., 7.]])
train_y = np.array([190, 330, 660]).reshape(3, 1)
y = None
for i in range(100):
    if y is not None:
        # 勾配を求める
        grad = gradient(y, train_y, train_x)
        # 重みを求める(横ベクトルを一次元に変換)
        w = weight(w, 0.02, grad).reshape(-1,)
    else:
        # 初期重み
        w = np.array([10, 20])

    # 推測値を求める
    y = np.dot(w, train_x.T)
    # 誤差を求める
    l = mean_squared_error(y, train_y)
    if (i + 1) % 10 == 0:
        print("step=%3d, a1=%6.2f, a2=%6.2f, loss=%.2f" % (i + 1, w[0], w[1], l))
print("Estimated: a1=%6.2f, a2=%6.2f" % (w[0], w[1]))

ディープラーニング(深層学習)を理解してみる(勾配降下法：計算方法)

機械学習 tensorflow 人工知能

はじめに

前回の続きです。 yaju3d.hatenablog.jp

幾つかの人工知能関連の本やWebサイトを見ても、数式やプログラムのソースリストは記載されていても、数学が苦手な自分が理解できるようになるまでの説明が無い、そんな中でも下記3つの本(Kindle)がまだ理解できそうな感じで参考になりそうである。

　　　　

基礎的な知識から、やっと実際の計算方法に入っていきます。

勾配降下法

関数 z=x^2 + y^2 について、その最小値を与える x と y の値を勾配降下法で求めてみます。
ちなみに、正解は(x, y)=(0, 0)です。
f:id:Yaju3D:20171028224142p:plain

最初に勾配を求めておきましょう。
前回記事の偏微分で説明したように、関数 z=x^{2} + y^{2} を偏微分すると指数を係数にした \displaystyle\frac{\partial z}{\partial x}=2x、 \displaystyle\frac{\partial z}{\partial y}=2y となります。

①勾配式　\displaystyle\left(\frac{\partial z}{\partial x}, \frac{\partial z}{\partial y} \right) = (2x, 2y)

それでは、ステップを追って計算を進めます。

１.初期設定

初期位置と学習係数 \eta を適当に与えます。
今回は初期位置を(3.00,2.00)、学習係数 \eta = 0.1 とします。

No 位 置 勾 配 変位ベ クトル 関数値
i x_i y_i \partial z / \partial x \partial z / \partial y \Delta x \Delta y z
0 3.00 2.00

２．変位ベクトルを算出

現在位置 (x_i, y_i) に対して、勾配式から算出し、勾配降下法の基本式から変位ベクトル \Delta x = (\Delta x_i, \Delta y_i)を求めます。
前回記事の勾配降下法に適用で、2変数関数 z=f(x,y)の勾配降下法の基本式は次のように表しました。
基本式 \displaystyle({\Delta x}, {\Delta y}) = -\eta\nabla f(x,y)

これに①勾配式を当てはめたのが次の式となります。「\cdot」は掛け算の意味です。
②変位ベクトル (\Delta x_i, \Delta y_i) = - \eta(2x, 2y) = (- \eta\cdot2x, - \eta\cdot2y)

No 位 置 勾 配 変位ベ クトル 関数値
i x_i y_i \partial z / \partial x \partial z / \partial y \Delta x \Delta y z
0 3.00 2.00 6.00 4.00 -0.60 -0.40 13.00
各計算結果の求め方

勾配
\partial z / \partial x = 2 \cdot x_0 = 2 \times 3.00 = 6.00
\partial z / \partial y = 2 \cdot y_0 = 2 \times 2.00 = 4.00
変位ベクトル
\Delta x = - \eta \cdot \partial z / \partial x = -0.1 \times 6.00 = -0.60
\Delta y = - \eta \cdot \partial z / \partial y = -0.1 \times 4.00 = -0.40
関数値
z = x_0^2 + y_0^2 = 3.00^2 + 2.00^2 = 9.00 + 4.00 = 13.00

３．位置を更新

勾配降下法に従って、現在位置 (x_i, y_i) から移動先 (x_{i+1}, y_{i+1}) の点を次の式から求めます。
②移動先 (x_{i+1}, y_{i+1}) = (x_i, y_i) + (\Delta x_i, \Delta y_i)

No 位 置 勾 配 変位ベ クトル 関数値
i x_i y_i \partial z / \partial x \partial z / \partial y \Delta x \Delta y z
0 3.00 2.00 6.00 4.00 -0.60 -0.40 13.00
1 2.40 1.60

移動先
x_1 = x_0 + \Delta x_0 = 3.00 + (-0.60) = 2.40
y_1 = y_0 + \Delta y_0 = 2.00 + (-0.40) = 1.60

４．２と３の繰り返し

２と３の繰り返し(6～27は省略)、30回繰り返したときの座標 (x_{30}, y_{30}) の値です。
正解の (x, y) = (0, 0) と一致します。

No 位 置 勾 配 変位ベ クトル 関数値
i x_i y_i \partial z / \partial x \partial z / \partial y \Delta x \Delta y z
0 3.00 2.00 6.00 4.00 -0.60 -0.40 13.00
1 2.40 1.60 4.80 3.20 -0.48 -0.32 8.32
2 1.92 1.28 3.84 2.56 -0.38 -0.26 5.32
3 1.54 1.02 3.07 2.05 -0.31 -0.20 3.41
4 1.23 0.82 2.46 1.64 -0.25 -0.16 2.18
5 0.96 0.66 1.97 1.31 -0.20 -0.13 1.40
28 0.01 0.00 0.01 0.01 0.00 0.00 0.00
29 0.00 0.00 0.01 0.01 0.00 0.00 0.00
30 0.00 0.00 0.01 0.00 0.00 0.00 0.00

※Excel計算で小数誤差により微妙に値が違います。

バイアスについて

バイアス(bias)とは、一般に真値からの偏り、つまり系統的な誤差を指す。

切片(せっぺん)

ニューラルネットワークのパラメーターは重みとバイアスがセットとなります。
バイアスがイメージしやすいものとして、回帰直線があります。
f:id:Yaju3D:20171105230836p:plain
回帰直線は次のような1次式で表現されます。
回帰方程式 y = ax + b
f:id:Yaju3D:20171106001201p:plain
a を回帰係数(傾き)、bを切片の呼びます。切片が無いと必ず原点を通すことになってしまいます。
傾きが求まったところで切片で上下位置の調整をします。この切片がバイアスのことになります。

閾値(しきいち)

f:id:Yaju3D:20170528144340p:plainf:id:Yaju3D:20170528154127p:plain
左図が本物の神経細胞(ニューロン) で、右図が形式ニューロンです。
簡単に説明すると、入力が2つあり各入力に対して重みが掛け算され、その値が閾値を超えれば出力は「1」、そうでなければ出力は「0」となります。 たとえば、入力が(1,0)、重みが(0.5, 0.7)だとすると、1×0.5 ＋ 0×0.7 ＝ 0.5 を計算して閾値と比較します。
閾値未満だと出力無し(発火なし)、閾値を超えると出力有り(発火あり)となります。

出力信号無し (y=0):w_1x_1+w_2x_2 \lt \theta
出力信号有り (y=1):w_1x_1+w_2x_2 \geqq \theta

ここで \theta はニューロン固有の閾値です。 \theta を左に移行させた発火の式は次のように表現することができます。
発火の式 y=a(w_1x_1+w_2x_2 - \theta)
※a は活性化関数となります。

活性化関数にシグモイド曲線を使った場合の閾値 \theta は生物的にはニューロンの個性を表現する値です。
\theta が大きければ興奮しにくく(すなわち鈍感)、小さければ興奮しやすい(すなわち敏感)という感受性を表します。
f:id:Yaju3D:20171106013953p:plain

発火の式 y=a(w_1x_1+w_2x_2 - \theta) にて、\theta だけマイナス記号が付いているのは数学的に美しくありません。美しさが欠けることは数学が嫌うところです。また、マイナスは計算ミスを誘発しやすいという欠点を持ちます。そこで、-\theta を b と置き換えたのが次の式となります。
y=a(w_1x_1+w_2x_2 + b)
こうすれば式として美しく、計算ミスも起こりにくくなります。
この b をバイアス(bias)と呼びます。

バイアスの式表現

数式では b とするよりは w_0 として、次の式にします。
y=a(w_0 + w_1x_1+w_2x_2)

一般的に書き直すと、重みをwではなく \theta として、\theta_0 をバイアスとした場合、入力値 x とすると\theta_0 と x で次元が違うと扱いにくいので、最初の要素に 1 をセットする。
f:id:Yaju3D:20180103202510p:plain

それを x_0=1 と定義して、 x の最初の要素に x_0 を置くほうがより数学上では綺麗となる。
f:id:Yaju3D:20180103210343p:plain

\theta を転置したものと x を掛けたものを計算すると次の式になります。
\theta^Tx = \theta_0x_0 + \theta_1x_1+ \theta_2x_2 + \cdot + \theta_nx_n
更にこれを書き直すと 、f_\theta(x) = \theta^Tx と簡易的な表現の式となる。

やる夫で学ぶ機械学習 - 多項式回帰と重回帰 - · けんごのお屋敷

リンゴとミカン

リンゴとミカンを例題に勾配降下法を計算していきます。上記と違うのは、訓練データがあることです。 f:id:Yaju3D:20160417230234p:plain
f:id:Yaju3D:20160419005112p:plain f:id:Yaju3D:20160419002003p:plain

以前の記事を参考にします。
yaju3d.hatenablog.jp

前提

f:id:Yaju3D:20171029165625p:plain
1層の全結合ニューラルネットワークを用いて、勾配降下法による重みの更新例を示します。
入力層ユニット数は2、出力層はユニット数が1のシンプルなネットワークです。損失関数は二乗誤差を用います。

No 入力 データ 教師データ
i x_1 x_2 t
1 1 3 190
2 3 1 330
3 5 7 660

上表はトレーニングデータセットです。サンプル数は3で、1番目のサンプルを見ると、x_1=1、x_2=3が与えられたときの教師データ t=190 (正解)になっています。

勾配\Delta Eの計算式

出力層の出力値 y は、入力層の出力値 x を用いて、次のような一次多項式で表すことができます。
a_1x_1+a_2x_2+b=y
パラメーター　a_1,a_2,b

これにトレーニングデータセットを当てはめると、次のような連結方程式が出来上がります。
a_1+3a_2+b=y_1
3a_1+a_2+b=y_2
5a_1+7a_2+b=y_3

この連結方程式は行列とベクトルを用いて次のように表すことが出来ます。

\begin{pmatrix}a_1 & a_2 & b \end{pmatrix} 
\begin{bmatrix}
1 & 3 & 5 \\ 
3 & 1 & 7 \\ 
1 & 1 & 1 \\ 
\end{bmatrix} = \begin{pmatrix}y_1 & y_2 & y_3 \end{pmatrix}  
    w　　　　　X　　　　   Y
   重み　　入力データ　出力データ(推測値)

今回はバイアスbは 0 とするので除外しました。そうしないリンゴとミカンの金額(重み)を求めたいのにバイアスを含む3つ値が求まってしまいます。

\begin{pmatrix}a_1 & a_2 \end{pmatrix} 
\begin{bmatrix}
1 & 3 & 5 \\ 
3 & 1 & 7 \\ 
\end{bmatrix} = \begin{pmatrix}y_1 & y_2 & y_3 \end{pmatrix}  
    w　　　　　X　　　　   Y
   重み　　入力データ　出力データ(推測値)

誤差を含めた式と行列を表現すると下記のようなります。
y=a_1x_1 + a_2x_2 + 誤差 e

\begin{bmatrix}
190 \\ 
330 \\ 
660
\end{bmatrix}=
\begin{bmatrix}
1 & 3 \\ 
3 & 1 \\ 
5 & 7
\end{bmatrix}
\begin{bmatrix}
a_1\\ 
a_2
\end{bmatrix} +
\begin{bmatrix}
e_1\\ 
e_2\\ 
e_3
\end{bmatrix}

入力データを X、重み(パラメーター)を w、出力データ(推測値) Y としています。w と X を用いると、出力データ Y は次の式で表すことができます。
wX=Y
また、教師データは、 t を用いて次のように表します。
t=(190 \quad 330 \quad 660)

損失関数には二乗誤差を使用します。誤差は次のような式で表すことができます。
E=\displaystyle \frac{1}{2} (Y - t)^2
\quad=\displaystyle \frac{1}{2} (wX - t)^2
勾配\Delta Eは、誤差 E を重み w で微分すると、勾配\Delta Eを次のような式になります。
\Delta E = \displaystyle \frac{\partial E}{\partial w} = (Y - t)X^T
X^T：X の転置行列、(Y - t)：誤差信号 \delta

今回、(Y - t)を誤差信号と呼び、記号\delta (デルタ)で表します。
勾配\Delta Eは、誤差信号 \delta と、入力データX から求めることができます。入力データはすでにわかっている値なので、誤差信号の値を求めれば、勾配\Delta Eを求めることができ、そして勾配\Delta Eがわかれば、式 w \leftarrow w-\eta\Delta E で重みを w を更新することができます。\eta：学習係数

初期値の設定

初めに、重みw と学習係数\etaに適当な初期値を設定します。ここでは学習係数\etaは 0.02 とし、重みwの初期値を10円と20円から、バイアスは 0 にして次のように設定します。
w=(a_1 \quad a_2) = (10 \quad 20)

重みの更新

①～③の手順で重みの更新を行います。

①現在の重みで推測値を求める

入力データXと現在の重みwから、推測値Yを求めます。
推測値

wX=\begin{pmatrix}10 & 20\end{pmatrix} 
\begin{bmatrix}
1 & 3 & 5 \\ 
3 & 1 & 7 \\ 
\end{bmatrix}
\quad\quad=(1 \times 10 + 3 \times 20 \quad 3 \times 10 + 1 \times 20 \quad 5 \times 10 + 7 \times 20)  
\quad\quad=(70 \quad 50 \quad 190)  
\quad\quad=Y

ここで、全体の誤差をいったん計算してみます。損失関数には二乗誤差を使用するので、全体の誤差Eは次のように計算します。

現在の値
教師データ t=(190 \quad 330 \quad 660)
推測値　　 Y=( 70 \quad 50 \quad 190)
全体の誤差
E=\displaystyle \sum_{n=1}^{3} \frac{1}{2}(y_n - t_n)^2
\quad=\displaystyle \frac{1}{2}(y_1 - t_1)^2 + \frac{1}{2}(y_2 - t_2)^2 + \frac{1}{2}(y_3 - t_3)^2
\quad=\displaystyle \frac{1}{2}(70 - 190)^2 + \frac{1}{2}(50 - 330)^2 + \frac{1}{2}(190 - 660)^2
\quad=\displaystyle \frac{1}{2}(-120)^2 + \frac{1}{2}(-280)^2 + \frac{1}{2}(-470)^2
\quad=\displaystyle 7200 + 39200 + 110450
\quad=\displaystyle 156850

②勾配\Delta Eを計算する

現在の重みwに対する勾配\Delta Eを計算します。
入力データを転置行列X^Tにするのは、行列の掛け算の定義で「行列Aの列数(横の個数)と、行列Bの行数(縦の個数)が等しくないと掛け算できない」ためです。
行列A(Y-t) 1行3列 と 行列B(X^T) 3行2列 にして掛け算させます。行列の積

現在の値
教師データ t=(190 \quad 330 \quad 660)
推測値　　 Y=( 70 \quad 50 \quad 190)

入力データ　X^T = \begin{bmatrix}
1 & 3 \\ 
3 & 1 \\ 
5 & 7 \\ 
\end{bmatrix}  

勾配
\Delta E=(Y-t)X^T
\quad\quad=(70-190 \quad 50-330 \quad 190-660)\begin{bmatrix}
1 & 3 \\ 
3 & 1 \\ 
5 & 7 \\ 
\end{bmatrix}
\quad\quad=(-120 \quad -280 \quad -470)\begin{bmatrix}
1 & 3 \\ 
3 & 1 \\ 
5 & 7 \\ 
\end{bmatrix}
\quad\quad=(-120 \times 1 + -280 \times 3 + -470 \times 5 \quad -120 \times 3 + -280 \times 1 + -470 \times 7)
\quad\quad=(-3310 \quad -3930)
③重みを更新する

②で求めた勾配\Delta Eを用いて、現在の重みwを更新します。学習係数\etaは 0.02 としています。

現在の値
現在の重み w=(10 \quad 20)
勾配　　 \Delta E=(-3310 \quad -3930)
更新後の重み
w \leftarrow w-\eta\Delta E
\quad=(10 \quad 20) - 0.02 \times (-3310 \quad -3930)
\quad=(10 \quad 20) - (-66.2 \quad -78.6)
\quad=(10 - (-66.2) \quad 20 - (-78.6))
\quad=(76.2 \quad 98.6)

これで重みwが更新されました。ここで、更新後の重みwを使って推測値を求め、全体の誤差を再び計算してみます。

現在の値
教師データ t=(190 \quad 330 \quad 660)

推測値

wX=\begin{pmatrix}76.2 & 98.6 \end{pmatrix} 
\begin{bmatrix}
1 & 3 & 5 \\ 
3 & 1 & 7 \\ 
\end{bmatrix}
\quad\quad=(1 \times 76.2 + 3 \times 98.6  \quad 3 \times 76.2 + 1 \times 98.6 \quad 5 \times 76.2 + 7 \times 98.6  
\quad\quad=(372 \quad 327.2 \quad 1071.2)  
\quad\quad=Y

全体の誤差
E=\displaystyle \sum_{n=1}^{3} \frac{1}{2}(y_n - t_n)^2
\quad=\displaystyle \frac{1}{2}(y_1 - t_1)^2 + \frac{1}{2}(y_2 - t_2)^2 + \frac{1}{2}(y_3 - t_3)^2
\quad=\displaystyle \frac{1}{2}(372 - 190)^2 + \frac{1}{2}(327.2 - 330)^2 + \frac{1}{2}(1071.2 - 660)^2
\quad=\displaystyle \frac{1}{2}(182)^2 + \frac{1}{2}(-2.8)^2 + \frac{1}{2}(411.2)^2
\quad=\displaystyle 16562 + 3.92 + 84542.72
\quad=\displaystyle 101108.64

全体の誤差は 101108.64 となり、重みの誤差前の 156850 より小さくなっていることが分かります。

更新の繰り返し

以上の①～③の計算で、全サンプル(今回は2サンプル)の更新が1回終わりました。これで1エポックを終了したことになります。
更新後の重みwを用いて、①～③をもう一度繰り返せば、2エポックが終了となります。

最後に

計算の繰り返しはコンピューターの得意なところです。
次回はこの計算が本当に合っているのか、Pythonを使って実証してみたいと思います。

スポンサーリンク