CNTK 2.2 Python API 解説 (7) - 強化学習の基礎: DQN, Policy Gradient / CNTK v2.3 リリース

#CNTK 2.2 Python API 解説 (7) - 強化学習の基礎: DQN, Policy Gradient / CNTK v2.3 リリース


##0. はじめに

◆ CNTK ( [Microsoft Cognitive Toolkit](https://www.microsoft.com/en-us/cognitive-toolkit/) ) 2.2 の Python API 解説第７弾です。
今回は強化学習がテーマで、OpenAI gym の CartPole を題材にして DQN と Policy Gradient アプローチで CNTK で実装します。


◆ **CNTK v2.3 リリース**

本題に入る前に、(もうご存知の方も多いかと思いますが、)
待望の <font color="red">CNTK v2.3 が 11月22日にリリースされました</font> ので、簡単にリリースノートを見ておきましょう :

>[CNTK v2.3 Release Notes](https://docs.microsoft.com/en-us/cognitive-toolkit/ReleaseNotes/CNTK_2_3_Release_Notes)

>**このリリースのハイライト**

>- より良い ONNX サポート。
- 分散トレーニングにおけるより良いパフォーマンスのために NCCL2 へスイッチしています。
- 改良された C# API。
- OpenCV は CNTK をインストールするために必要ではありません、それは Tensorboard 画像特徴と画像リーダーのために必要なだけです。
- 各種パフォーマンス改良。
- ネットワーク最適化 API の追加。
- スパースのためのより高速な Adadelta。

チュートリアルも新たに追加されていました。是非トライしてみてください。


◆ さて、強化学習は行動主義心理学により触発された機械学習の領域です。
ソフトウェア・エージェントが報酬を最大化するために、与えられた環境においてどのように行動を取るべきかに関心を寄せます。

ある種の機械学習の設定では正解ラベルへのアクセスを持ちませんので、教師あり学習のテクニックに頼ることはできません。
けれども、相互作用可能な何かがあり、そしてそれが教えてくれる何某かのフィードバックを得られる場合には、
動作をどのように改良するかを学習するために強化学習が利用できます。

定番ですが題材として OpenAI gym シミュレーターからの CartPole 環境を使用します。これはカートにポールのバランスを取ることを教えます :

![cntk15_polecart.gif](https://qiita-image-store.s3.amazonaws.com/0/203151/f014eb25-8fb2-462b-e18a-f1fb9d15cd1a.gif)

以下の GIF アニメは、本記事でトレーニングした DQN モデルのデモ・サンプルの実行です :

![cntk15_demo.gif](https://qiita-image-store.s3.amazonaws.com/0/203151/61081b37-744c-7eb4-07be-4f86adb8fd35.gif)


本記事の内容は :

1. 動作環境と Jupyter Notebook について
2. 強化学習
3. パート 1 : Deep Q-Network (DQN)
4. パート 2 : Policy Gradient (PG)

本記事は以下の CNTK チュートリアルを参考にしています :

- [CNTK 203: Reinforcement Learning Basics](https://github.com/Microsoft/CNTK/blob/master/Tutorials/CNTK_203_Reinforcement_Learning_Basics.ipynb)



## 1. 動作環境と Jupyter Notebook について

**動作環境**

動作環境の構築が必要な場合には、[Cognitive Toolkit 2.2 を Azure Linux GPU 仮想マシンにインストール](https://qiita.com/masao-classcat/items/5a922da9e698ec626d25) を参考にしてください。Azure ポータルと Ubuntu Linux にある程度慣れていれば、30 分程度で以下のような環境が構築できるかと思います :

- [Azure NC 仮想マシン](https://blogs.technet.microsoft.com/jpitpro/2016/11/22/azure-n-series-general-availability-on-december-1/) with NVIDIA Tesla® K80 GPU
- Ubuntu 16.04 LTS
- NVIDIA CUDA 8.0 & cuDNN 6.0
- Anaconda 3 4.1.1
- CNTK 2.2 (for GPU)


**Jupyter Notebook**

また、本記事でも CNTK チュートリアルでも Jupyter Notebook を多用します。
Jupyter Notebook の利用方法については「CNTK 2.2 Python API 入門 (2)」の記事中の [Jupyter Notebook の活用](https://qiita.com/masao-classcat/items/3512f212805e06abda5f#1-jupyter-notebook-%E3%81%AE%E6%B4%BB%E7%94%A8) を参照してください。



##2. 強化学習

###2-1 基本概念

強化学習 (Reinforcement learning, RL) は行動主義心理学 (behaviorist psychology) により触発された機械学習の領域です。[ソフトウェア・エージェント](https://en.wikipedia.org/wiki/Software_agent) が、累積する報酬 (の概念) を最大化するために、環境においてどのように行動 (= [actions](https://en.wikipedia.org/wiki/Action_selection)) を取るべきかに関心を寄せます。
機械学習では、環境は典型的にはマルコフ決定過程 ([Markov decision process](https://en.wikipedia.org/wiki/Markov_decision_process), MDP) として定式化できます。
このコンテキストのための多くの強化学習アルゴリズムは動的計画法 (= [Dynamic programming](https://en.wikipedia.org/wiki/Dynamic_programming)) のテクニックを利用するからです。


ある種の機械学習の設定では、正解ラベルへの即座のアクセスを持ちませんので、教師あり学習のテクニックに頼ることはできません。けれども、相互作用可能な何かがあり、そしてそれによって教えてくれる何某かのフィードバックを得られる場合には、以前の動作が良くても悪くても、動作をどのように改良するかを学習するために強化学習が利用できます。


教師あり学習のケースとは異なり、強化学習ではラベル付けられた正しい入出力ペアは決して提示されません。
そして最善ではない行動も明示的には決して正されません。
これは (前に決して学習されていない条件や行動の) 探索と (以前の経験から既に学習された条件または行動の) 搾取 (= exploitation) の間のバランスを見出すことを伴います。

>少し分かりにくい表現ですが、ある程度良いと分かっている行動を通して報酬を「搾取」するだけでなく、
より良い報酬を求めて他の行動を「探索」することも必要で、そのバランスを取る必要があることを意味しています。
もちろん、初期状態では「探索」が中心になります。


本記事で扱う、多腕バンディット問題 (Multi-arm bandit problems) は強化学習アルゴリズムのカテゴリの一つです。
そこでは探索 vs 搾取のトレードオフが徹底的に学習されます。




###2-2 問題設定 : CartPole (カートポール)

OpenAI [gym](https://github.com/openai/gym) シミュレーターからの [CartPole](https://gym.openai.com/envs/CartPole-v0) 環境を使用します。これはカートにポールのバランスを取ることを教えます。
CartPole では、ポールがカートに非駆動関節 (= un-actuated joint) によって装着されていて、カートが摩擦がないトラック (進路) に沿って移動します。
システムはカートに +1 か -1 の力を適用することによって制御します。
ポールが直立している総ての時間ステップに対して +1 の報酬が与えられます。
エピソードは、ポールが垂直から 15 度以上傾くか、またはカートが中央から 2.4 ユニット以上移動するときに終了します。

※ **エピソード** は学習の開始から終了までを指しますが、1 ゲームと考えれば分かりやすいでしょう。


以下はカートポールの図面です。詳細は [参考](http://www.simongrant.org/pubs/thesis/3/2.html) を参照してください :

```python:In_[1]==>
from IPython.display import Image
```

```python:In_[2]==>
Image(url="https://cntk.ai/jup/polecart.gif", width=300, height=300)
```

![cntk15_polecart.gif](https://qiita-image-store.s3.amazonaws.com/0/203151/f014eb25-8fb2-462b-e18a-f1fb9d15cd1a.gif)



**目標**

目標は、開始状態として (カートに対して垂直に) 直立しているポールを伴ってカートが動くときに、ポールが倒れることを回避することです。より具体的には、ポールが垂直から 15 度以内で、そしてカートが中央から 2.4 ユニット内であれば、報酬を集めることができます。
本記事では、最後の 50 バッチ (= 50 エピソード) で 200 かそれ以上の平均報酬に導かれるようなポリシー (i.e. 行動のセット) を学習するまでトレーニングします。


強化学習用語で換言すれば、目標は、<font color="red">ある環境 (この場合はカート上でバランスを取るポール) との相互作用を通して
報酬 $r$ (フィードバック) を最大化する、ポリシー $a$ を見つけること</font> です。


経験のシリーズ : $$s \xrightarrow{a} r, s'$$ が与えられたとき、
与えられた状態 $s$ において時間とともに累積していく報酬 $r$ を最大化するために、どのように行動 $a$ を選択するかを学習することができます :

$$
Q(s,a) = r_0 + \gamma r_1 + \gamma^2 r_2 + \ldots = r_0 + \gamma \max_a Q^*(s',a)
$$

ここで $\gamma \in [0,1)$ はディスカウント・ファクターで、(離れた時間ステップの) 報酬をどのくらい評価するかを制御します。これは [*Bellmann*-equation](https://en.wikipedia.org/wiki/Bellman_equation) と呼ばれます。



◆ 状態空間をどのようにモデル化するか、そしてどの行動が最高の未来の報酬を生成するかを見出すために、受け取った報酬をどのように使用するかを示します。
本記事では２つの異なるポピュラーなアプローチを提示して CNTK で実装しますが、DQN を中心的に扱います :

**DQN (Deep Q-Network)** : DQN は、Atari ゲームを生ピクセルだけからどのようにプレーするかをトレーニングするために成功的に使用され、2015 年に有名になりました。$Q(s,a)$ 値を学習するためにニューラルネットワークをトレーニングします (そのために $Q$-ネットワークと呼ばれます)。これらの $Q$ 関数値から最善の行動を選択します。


**Policy Gradient (ポリシー勾配)** : この方法はネットワークのポリシー (行動セット) を直接的に推定します。
行動の順序付けられたセットの学習で、これは行動のサブセットを確率的に選択することにより報酬を最大化することへと導きます。本記事では、ポリシーを学習するために勾配降下アプローチを使用して行動を学習します。



###2-3 準備

**インポート**

```python:In_[3]==>
from __future__ import print_function
from __future__ import division
import matplotlib.pyplot as plt
from matplotlib import style
import numpy as np
import pandas as pd
import seaborn as sns

style.use('ggplot')
%matplotlib inline
```


OpenAI gym パッケージがインストールされていない場合には、以下のコード・ブロックでインストールできます :

```python:In_[4]==>
try:
    import gym
except:
    !pip install gym
    import gym
```


例によって Fast モードと Slow モードが用意されています。デフォルトは Fast モードです :

```python:In_[5]==>
isFast = True
```

- GPU 環境であれば Slow モードで実行しても問題ありません (つまり、`isFast = False` で良いです)。

- 以後、トレーニングの結果やモデルは Slow モードで生成しています。




###2-4 CartPole: データと環境

カートにポールのバランスを取ることを教えるために OpenAI の gym シミュレータからの CartPole 環境を使用しますが、仕様の要点をまとめておきます。

総ての時間ステップで、エージェントは :

- 観測 $(x, \dot{x}, \theta, \dot{\theta})$ を得ます。それぞれ、カート位置、カート速度、(垂直からの) ポール角度、ポール角速度に相当します。

- `LEFT` または `RIGHT` の行動を遂行します、そして

- 以下を受け取ります :

  - もう一つの時間ステップを生存したことへの +1 の報酬、そして
  - 新しい状態 $(x', \dot{x}', \theta', \dot{\theta}')$


もし以下の条件に該当すれば、エピソードは終了します :

- ポールが垂直から 15 度以上、そして/または
- カートが中央から 2.4 ユニット以上移動する。


もし以下の条件に該当すれば、タスクは完了したとして考えられます :

- エージェントが最後の 50 エピソードで 200 の平均報酬を獲得。 


※ Fast モードではこれらのターゲットは緩和されます。




##3. パート 1 : DQN (Deep Q-Network)

###3-1 要点と準備

遷移 $(s,a,r,s')$ の後、価値関数 $Q(s,a)$ をターゲット $r+\gamma \max_{a'}Q(s',a')$ に近づけようとします、
ここで $\gamma$ は未来への報酬のためのディスカウント・ファクターで 0 と 1 の範囲の値を取ります。


**DQN**

- 観測 (状態, 行動) を `スコア` にマップする Q-関数を学習します。

- メモリ・リプレーを使用します。

- (本記事では扱いませんが、) 学習を安定させるために第２のネットワークを使用することもできます。




###3-2 DQN モデル

CNTK で DQN モデルを実装していきますが、ここでは簡単な浅いネットワーク - 2 層完全結合ネットワークを選択します :

$$
l_1 = relu( x W_1 + b_1)
$$ $$
Q(s,a) = l_1 W_2 + b_2 \\
$$

```python:入力不要
        l1 = C.layers.Dense(H, activation=C.relu)  # H = 64
        l2 = C.layers.Dense(ACTION_COUNT) # ACTION_COUNT = 2
```

※ 浅いモデルをより深いネットワークで置き換えることによりアプローチを拡張できます。

追加のインポートです :

```python:In_[6]==>
import numpy
import math 
import os 
import random

import cntk as C
```


最初に CartPole の環境を確認しておきます :

- STATE_COUNT = 4 ( $(x, \dot{x}, \theta, \dot{\theta})$ に相当 )

- ACTION_COUNT = 2 ( `LEFT` または `RIGHT` に相当 )

```python:In_[7]==>
env = gym.make('CartPole-v0')

STATE_COUNT  = env.observation_space.shape[0]
ACTION_COUNT = env.action_space.n

STATE_COUNT, ACTION_COUNT
```
```bash
(4, 2)
```


次に `Brain` クラスを定義します。このクラスはモデル定義に加えて、
損失・オプティマイザーの作成からの trainer オブジェクトの生成を含む、CNTK の基本的なプリミティブの単なるラッパーです :

```python:In_[8]==>
# ターゲット報酬
REWARD_TARGET = 30 if isFast else 200

# Averaged over these these many episodes
BATCH_SIZE_BASELINE = 20 if isFast else 50

H = 64 # hidden layer size

class Brain:
    def __init__(self):
        self.params = {}
        self.model, self.trainer, self.loss = self._create()
        # self.model.load_weights("cartpole-basic.h5")

    def _create(self):
        observation = C.sequence.input_variable(STATE_COUNT, np.float32, name="s")
        q_target = C.sequence.input_variable(ACTION_COUNT, np.float32, name="q")

        # モデル
        l1 = C.layers.Dense(H, activation=C.relu)
        l2 = C.layers.Dense(ACTION_COUNT)
        unbound_model = C.layers.Sequential([l1, l2])
        model = unbound_model(observation)

        self.params = dict(W1=l1.W, b1=l1.b, W2=l2.W, b2=l2.b)

        # 損失関数: loss='mse'
        loss = C.reduce_mean(C.square(model - q_target), axis=0)
        meas = C.reduce_mean(C.square(model - q_target), axis=0)

        # オプティマイザー
        lr = 0.00025
        lr_schedule = C.learning_rate_schedule(lr, C.UnitType.minibatch)
        learner = C.sgd(model.parameters, lr_schedule, gradient_clipping_threshold_per_sample=10)
        trainer = C.Trainer(model, (loss, meas), learner)

        # CNTK: trainer と loss を返します。
        return model, trainer, loss

    def train(self, x, y, epoch=1, verbose=0):
        #self.model.fit(x, y, batch_size=64, nb_epoch=epoch, verbose=verbose)
        arguments = dict(zip(self.loss.arguments, [x,y]))
        updated, results =self.trainer.train_minibatch(arguments, outputs=[self.loss.output])

    def predict(self, s):
        return self.model.eval([s])
```


`Memory` クラスは異なる状態 ($s$)、行動 ($a$) そして報酬 ($r$) をストアします :

```python:In_[9]==>
class Memory:   # stored as ( s, a, r, s_ )
    samples = []

    def __init__(self, capacity):
        self.capacity = capacity

    def add(self, sample):
        self.samples.append(sample)

        if len(self.samples) > self.capacity:
            self.samples.pop(0)

    def sample(self, n):
        n = min(n, len(self.samples))
        return random.sample(self.samples, n)
```


`Agent` クラスは `Brain` オブジェクトと `Memory` オブジェクトを使用して、
報酬を最大化する行動の最適なセットを選択するために過去の行動をリプレーします :

```python:In_[10]==>
MEMORY_CAPACITY = 100000
BATCH_SIZE = 64

GAMMA = 0.99 # ディスカウント・ファクター

MAX_EPSILON = 1
MIN_EPSILON = 0.01 # stay a bit curious even when getting old
LAMBDA = 0.0001    # speed of decay

class Agent:
    steps = 0
    epsilon = MAX_EPSILON

    def __init__(self):
        self.brain = Brain()
        self.memory = Memory(MEMORY_CAPACITY)

    def act(self, s):
        if random.random() < self.epsilon:
            return random.randint(0, ACTION_COUNT-1)
        else:
            return numpy.argmax(self.brain.predict(s))

    def observe(self, sample):  # in (s, a, r, s_) format
        self.memory.add(sample)

        # slowly decrease Epsilon based on our eperience
        self.steps += 1
        self.epsilon = MIN_EPSILON + (MAX_EPSILON - MIN_EPSILON) * math.exp(-LAMBDA * self.steps)

    def replay(self):
        batch = self.memory.sample(BATCH_SIZE)
        batchLen = len(batch)

        no_state = numpy.zeros(STATE_COUNT)


        # CNTK: float32 に明示的に設定します。
        states = numpy.array([ o[0] for o in batch ], dtype=np.float32)
        states_ = numpy.array([(no_state if o[3] is None else o[3]) for o in batch ], dtype=np.float32)

        p = agent.brain.predict(states)
        p_ = agent.brain.predict(states_)

        # CNTK: float32 に明示的に設定します。
        x = numpy.zeros((batchLen, STATE_COUNT)).astype(np.float32)
        y = numpy.zeros((batchLen, ACTION_COUNT)).astype(np.float32)

        for i in range(batchLen):
            s, a, r, s_ = batch[i]

            # CNTK: [0] because of sequence dimension
            t = p[0][i]
            if s_ is None:
                t[a] = r
            else:
                t[a] = r + GAMMA * numpy.amax(p_[0][i])

            x[i] = s
            y[i] = t

        self.brain.train(x, y)
```



###3-3 トレーニング

行動の初期状態は見当違いの荒々しい (= wild) 探検となり、そして反復を繰り返していくことによって
より長い実行を生成してより報酬を収集するような行動の範囲をシステムが学習することを期待します。
本記事では、下で [epsilon-greedy](https://en.wikipedia.org/wiki/Reinforcement_learning) (a.k.a. $\epsilon$-greedy) アプローチを実装します。


以下のコード・ブロックは重みのヒートマップを可視化します :

```python:In_[11]==>
def plot_weights(weights, figsize=(7,5)):
    '''Heat map of weights to see which neurons play which role'''
    sns.set(style="white")
    f, ax = plt.subplots(len(weights), figsize=figsize)
    cmap = sns.diverging_palette(220, 10, as_cmap=True)
    
    for i, data in enumerate(weights):
        axi = ax if len(weights)==1 else ax[i]
        if isinstance(data, tuple):
            w, title = data
            axi.set_title(title)
        else:
            w = data
                
        sns.heatmap(w.asarray(), cmap=cmap, square=True, center=True, #annot=True,
                    linewidths=.5, cbar_kws={"shrink": .25}, ax=axi)
```


**探検 (Exploration) - 搾取 (exploitation) のトレードオフ**

初期 $\epsilon$ は、完全に探検を行なう 1 に設定されます。
しかし、ステップが増えるにつれて探検を減らして報酬を収集する (= 搾取する) ために学習された空間を活用し始めます :

```python:In_[12]==>
def epsilon(steps):
    return MIN_EPSILON + (MAX_EPSILON - MIN_EPSILON) * np.exp(-LAMBDA * steps)
plt.plot(range(10000), [epsilon(x) for x in range(10000)], 'r')
plt.xlabel('step');plt.ylabel('$\epsilon$')
```

![cntk15_epsilon.png](https://qiita-image-store.s3.amazonaws.com/0/203151/fb8f01c4-e37d-4e0c-92bf-109541d5355f.png)


◆ DQN を使用してエージェントをトレーニングする準備ができましたので、実行してみましょう。モデルは `dqn.mod` としてセーブされます :

```python:In_[13]==>
TOTAL_EPISODES = 2000 if isFast else 3000

def run(agent):
    s = env.reset()
    R = 0

    while True:
        # Uncomment the line below to visualize the cartpole
        # env.render()

        # CNTK: explicitly setting to float32
        a = agent.act(s.astype(np.float32))

        s_, r, done, info = env.step(a)

        if done: # terminal state
            s_ = None

        agent.observe((s, a, r, s_))
        agent.replay()

        s = s_
        R += r

        if done:
            return R

agent = Agent()

episode_number = 0
reward_sum = 0
while episode_number < TOTAL_EPISODES:
    reward_sum += run(agent)
    episode_number += 1
    if episode_number % BATCH_SIZE_BASELINE == 0:
        print('Episode: %d, Average reward for episode %f.' % (episode_number,
                                                               reward_sum / BATCH_SIZE_BASELINE))
        if episode_number%200==0:
            plot_weights([(agent.brain.params['W1'], 'Episode %i $W_1$'%episode_number)], figsize=(14,5))
        if reward_sum / BATCH_SIZE_BASELINE > REWARD_TARGET:
            print('Task solved in %d episodes' % episode_number)
            plot_weights([(agent.brain.params['W1'], 'Episode %i $W_1$'%episode_number)], figsize=(14,5))
            break
        reward_sum = 0
agent.brain.model.save('dqn.mod')
```

```bash:Slow_モード
Episode: 50, Average reward for episode 21.080000.
Episode: 100, Average reward for episode 19.240000.
Episode: 150, Average reward for episode 18.680000.
Episode: 200, Average reward for episode 18.080000.
Episode: 250, Average reward for episode 17.760000.
Episode: 300, Average reward for episode 15.320000.
Episode: 350, Average reward for episode 14.800000.
Episode: 400, Average reward for episode 14.040000.
Episode: 450, Average reward for episode 13.160000.
Episode: 500, Average reward for episode 13.260000.
Episode: 550, Average reward for episode 12.200000.
Episode: 600, Average reward for episode 13.260000.
Episode: 650, Average reward for episode 12.180000.
Episode: 700, Average reward for episode 12.260000.
Episode: 750, Average reward for episode 12.040000.
Episode: 800, Average reward for episode 12.300000.
Episode: 850, Average reward for episode 14.740000.
Episode: 900, Average reward for episode 33.980000.
Episode: 950, Average reward for episode 13.340000.
Episode: 1000, Average reward for episode 10.660000.
Episode: 1050, Average reward for episode 10.380000.
Episode: 1100, Average reward for episode 11.000000.
Episode: 1150, Average reward for episode 11.180000.
Episode: 1200, Average reward for episode 11.680000.
Episode: 1250, Average reward for episode 11.140000.
Episode: 1300, Average reward for episode 10.500000.
Episode: 1350, Average reward for episode 10.580000.
Episode: 1400, Average reward for episode 10.420000.
Episode: 1450, Average reward for episode 10.680000.
Episode: 1500, Average reward for episode 39.440000.
Episode: 1550, Average reward for episode 47.720000.
Episode: 1600, Average reward for episode 28.560000.
Episode: 1650, Average reward for episode 26.040000.
Episode: 1700, Average reward for episode 24.760000.
Episode: 1750, Average reward for episode 26.300000.
Episode: 1800, Average reward for episode 26.300000.
Episode: 1850, Average reward for episode 27.620000.
Episode: 1900, Average reward for episode 28.820000.
Episode: 1950, Average reward for episode 33.060000.
Episode: 2000, Average reward for episode 35.340000.
Episode: 2050, Average reward for episode 42.400000.
Episode: 2100, Average reward for episode 53.060000.
Episode: 2150, Average reward for episode 62.140000.
Episode: 2200, Average reward for episode 57.060000.
Episode: 2250, Average reward for episode 74.300000.
Episode: 2300, Average reward for episode 89.100000.
Episode: 2350, Average reward for episode 100.820000.
Episode: 2400, Average reward for episode 114.280000.
Episode: 2450, Average reward for episode 147.280000.
Episode: 2500, Average reward for episode 188.360000.
Episode: 2550, Average reward for episode 198.220000.
Episode: 2600, Average reward for episode 199.980000.
Episode: 2650, Average reward for episode 200.000000.
Episode: 2700, Average reward for episode 200.000000.
Episode: 2750, Average reward for episode 200.000000.
Episode: 2800, Average reward for episode 200.000000.
Episode: 2850, Average reward for episode 200.000000.
Episode: 2900, Average reward for episode 200.000000.
Episode: 2950, Average reward for episode 200.000000.
Episode: 3000, Average reward for episode 200.000000.
CPU times: user 10min 38s, sys: 14 s, total: 10min 52s
Wall time: 10min 56s
```

報酬が 200 に達していることが分かります。
また、以下のプロットは重みのヒートマップの抜粋で、トレーニングの最後に併せて出力されます :

![cntk15_episode200.png](https://qiita-image-store.s3.amazonaws.com/0/203151/9bd9b1f4-f982-0a7a-caf7-5442c2d1b217.png)

![cntk15_episode600.png](https://qiita-image-store.s3.amazonaws.com/0/203151/24930c68-914a-3af2-b049-7114864240c4.png)

![cntk15_episode1000.png](https://qiita-image-store.s3.amazonaws.com/0/203151/76cf81e6-22fb-d5f1-72d3-65c90e2d7152.png)

![cntk15_episode1600.png](https://qiita-image-store.s3.amazonaws.com/0/203151/8243bf86-e78a-14ee-8db0-fb12cef4dc34.png)

![cntk15_episode2000.png](https://qiita-image-store.s3.amazonaws.com/0/203151/63b4cfb6-2a15-f109-e62a-eadd648f6811.png)

![cntk15_episode2600.png](https://qiita-image-store.s3.amazonaws.com/0/203151/84cece75-ed01-6b6c-d450-1627d37782b5.png)

![cntk15_episode3000.png](https://qiita-image-store.s3.amazonaws.com/0/203151/9be371dc-5e86-4944-7939-2801e925307a.png)




###3-4 DQN モデルを実行する

セーブされたモデル `dqn.mod` をロードすれば CartPole のデモを実行することができます。
但し、グラフィカルなデモ・システムですので、(X-) Window システムが必要になりますのでご注意ください :

```python:In_[14]==>
env = gym.make('CartPole-v0')

num_episodes = 10  # number of episodes to run

modelPath = 'dqn.mod'
#modelPath = 'dqn.mod'
root = C.load_model(modelPath)

for i_episode in range(num_episodes):
    print(i_episode)
    observation = env.reset()  # reset environment for new episode
    done = False
    while not done: 
        if not 'TEST_DEVICE' in os.environ:
            env.render()
        action = np.argmax(root.eval([observation.astype(np.float32)]))
        observation, reward, done, info  = env.step(action)    
```
```
0
1
2
3
4
5
6
7
8
9
```


以下のデモ・サンプルは、GPU 装備のクラウド仮想マシン上でトレーニングしたモデル `dqn.mod` を
ローカルで 10 エピソード実行した動画をキャプチャーしています :

![cntk15_demo.gif](https://qiita-image-store.s3.amazonaws.com/0/203151/61081b37-744c-7eb4-07be-4f86adb8fd35.gif)




##4. パート 2 : Policy Gradient (PG)

Policy Gradient については要点だけを簡単に示していきます。


###4-1 要点と準備

**目標 :**
$$
\text{maximize } E [R | \pi_\theta]
$$


**アプローチ** :

1. 経験を収集します (i.e $(s,a)$ 空間を通して多くの軌跡をサンプリングします)。

2. 良い経験がより起こりうるようにポリシーを更新します。


**DQN との違い** :

- 単一の $(s,a,r,s')$ 遷移は考えません、しかし勾配更新のためにエピソード全体を使用します。

- DQN では価値関数をモデル化 (出力は生のスコア) するところで、PG ではパラメータはポリシーを直接的にモデル化します (出力は行動確率)。


**報酬** :
総ての時間ステップのために +1 報酬を得ます、そこでは依然としてゲーム中です。
>問題点 : どの行為がゲームの継続に繋がるか、そしてどれが実際に悪い選択であったかを通常は知り得ません。
単純な発見的手法 (heuristic) : エピソードの始まりのほうの行為は概して良く、終わりに向かうにつれてそれらは悪くなりがちです (それらは結局、ゲームを失敗することに繋がりますから)。

報酬の 1D 浮動小数点配列を取り、ディスカウントされた報酬を計算します :

```python:In_[15]==>
def discount_rewards(r, gamma=0.999):
    """Take 1D float array of rewards and compute discounted reward """
    discounted_r = np.zeros_like(r)
    running_add = 0
    for t in reversed(range(0, r.size)):
        running_add = running_add * gamma + r[t]
        discounted_r[t] = running_add
    return discounted_r
```

プロットしてみます :

```python:In_[16]==>
discounted_epr = discount_rewards(np.ones(10))
f, ax = plt.subplots(1, figsize=(5,2))
sns.barplot(list(range(10)), discounted_epr, color="steelblue")
```

![cntk15_pg1.png](https://qiita-image-store.s3.amazonaws.com/0/203151/630ec709-be21-a1d8-a537-9fdbac470958.png)


報酬を正規化します、これによって終わりに向けてゼロより下にまで落ちます。
gamma は報酬が落ちるのをどの程度遅くするかを制御します :

```python:In_[17]==>
discounted_epr_cent = discounted_epr - np.mean(discounted_epr)
discounted_epr_norm = discounted_epr_cent/np.std(discounted_epr_cent)
f, ax = plt.subplots(1, figsize=(5,2))
sns.barplot(list(range(10)), discounted_epr_norm, color="steelblue")
```

![cntk15_pg2.png](https://qiita-image-store.s3.amazonaws.com/0/203151/15e64b73-a649-c9ca-51ea-391f70c41d99.png)

```python:In_[18]==>
discounted_epr = discount_rewards(np.ones(10), gamma=0.5)
discounted_epr_cent = discounted_epr - np.mean(discounted_epr)
discounted_epr_norm = discounted_epr_cent/np.std(discounted_epr_cent)
f, ax = plt.subplots(2, figsize=(5,3))
sns.barplot(list(range(10)), discounted_epr, color="steelblue", ax=ax[0])
sns.barplot(list(range(10)), discounted_epr_norm, color="steelblue", ax=ax[1])
```

![cntk15_pg3.png](https://qiita-image-store.s3.amazonaws.com/0/203151/4cee77a7-c0af-749d-6608-89dd1d8aa322.png)




###4-2 Policy Gradient モデル

policy gradient アプローチでは、Dense 層の出力は sigmod 関数を通して 0-1 の範囲にマップされます :

$$
l_1 = relu( x W_1 + b_1)
$$$$
l_2 = l_1 W_2 + b_2
$$$$
\pi(a|s) = sigmoid(l_2)
$$

```python:In_[19]==>
TOTAL_EPISODES = 2000 if isFast else 10000

D = 4  # 入力次元
H = 10 # 隠れ層ニューロンの数

observations = C.sequence.input_variable(STATE_COUNT, np.float32, name="obs")

W1 = C.parameter(shape=(STATE_COUNT, H), init=C.glorot_uniform(), name="W1")
b1 = C.parameter(shape=H, name="b1")
layer1 = C.relu(C.times(observations, W1) + b1)

W2 = C.parameter(shape=(H, ACTION_COUNT), init=C.glorot_uniform(), name="W2")
b2 = C.parameter(shape=ACTION_COUNT, name="b2")
score = C.times(layer1, W2) + b2
# Until here it was similar to DQN

probability = C.sigmoid(score, name="prob")
```



###4-3 PG モデルのトレーニングと実行

**ポリシー検索 (= Policy Search)** : 

最適なポリシー検索は勾配フリーなアプローチか、$\theta$ でパラメータ化されるポリシー空間 ($\pi_\theta$) に渡る勾配を計算することにより遂行されます。

本記事では、パラメータ化された空間 $\theta$ に渡るエラーの古典的な forward (`loss.forward`) と back (`loss.backward`) 伝播を使用します。この場合、$\theta = \{W_1, b_1, W_2, b_2\}$ がモデル・パラメータです :

```python:In_[20]==>
input_y = C.sequence.input_variable(1, np.float32, name="input_y")
advantages = C.sequence.input_variable(1, np.float32, name="advt")

loss = -C.reduce_mean(C.log(C.square(input_y - probability) + 1e-4) * advantages, axis=0, name='loss')

lr = 0.001
lr_schedule = C.learning_rate_schedule(lr, C.UnitType.sample)
sgd = C.sgd([W1, W2], lr_schedule)

gradBuffer = dict((var.name, np.zeros(shape=var.shape)) for var in loss.parameters if var.name in ['W1', 'W2', 'b1', 'b2'])

xs, hs, label, drs = [], [], [], []
running_reward = None
reward_sum = 0
episode_number = 1

observation = env.reset()

while episode_number <= TOTAL_EPISODES:
    x = np.reshape(observation, [1, STATE_COUNT]).astype(np.float32)

    # Run the policy network and get an action to take.
    prob = probability.eval(arguments={observations: x})[0][0][0]
    action = 1 if np.random.uniform() < prob else 0

    xs.append(x)  # 観測
    # grad that encourages the action that was taken to be taken

    y = 1 if action == 0 else 0  # a "fake label"
    label.append(y)

    # step the environment and get new measurements
    observation, reward, done, info = env.step(action)
    reward_sum += float(reward)

    # Record reward (has to be done after we call step() to get reward for previous action)
    drs.append(float(reward))

    if done:
        # Stack together all inputs, hidden states, action gradients, and rewards for this episode
        epx = np.vstack(xs)
        epl = np.vstack(label).astype(np.float32)
        epr = np.vstack(drs).astype(np.float32)
        xs, label, drs = [], [], []  # reset array memory

        # Compute the discounted reward backwards through time.
        discounted_epr = discount_rewards(epr)
        # Size the rewards to be unit normal (helps control the gradient estimator variance)
        discounted_epr -= np.mean(discounted_epr)
        discounted_epr /= np.std(discounted_epr)

        # Forward パス
        arguments = {observations: epx, input_y: epl, advantages: discounted_epr}
        state, outputs_map = loss.forward(arguments, outputs=loss.outputs,
                                          keep_for_backward=loss.outputs)

        # Backward パス
        root_gradients = {v: np.ones_like(o) for v, o in outputs_map.items()}
        vargrads_map = loss.backward(state, root_gradients, variables=set([W1, W2]))

        for var, grad in vargrads_map.items():
            gradBuffer[var.name] += grad

        # Wait for some batches to finish to reduce noise
        if episode_number % BATCH_SIZE_BASELINE == 0:
            grads = {W1: gradBuffer['W1'].astype(np.float32),
                     W2: gradBuffer['W2'].astype(np.float32)}
            updated = sgd.update(grads, BATCH_SIZE_BASELINE)

            # reset the gradBuffer
            gradBuffer = dict((var.name, np.zeros(shape=var.shape))
                              for var in loss.parameters if var.name in ['W1', 'W2', 'b1', 'b2'])

            print('Episode: %d. Average reward for episode %f.' % (episode_number, reward_sum / BATCH_SIZE_BASELINE))

            if reward_sum / BATCH_SIZE_BASELINE > REWARD_TARGET:
                print('Task solved in: %d ' % episode_number)
                break

            reward_sum = 0

        observation = env.reset()  # reset env
        episode_number += 1
probability.save('pg.mod')
```

```bash
Episode: 50. Average reward for episode 21.600000.
Episode: 100. Average reward for episode 21.880000.
Episode: 150. Average reward for episode 18.360000.
Episode: 200. Average reward for episode 18.840000.
Episode: 250. Average reward for episode 18.620000.
Episode: 300. Average reward for episode 21.040000.
Episode: 350. Average reward for episode 17.920000.
Episode: 400. Average reward for episode 18.200000.
Episode: 450. Average reward for episode 19.120000.
Episode: 500. Average reward for episode 15.220000.
Episode: 550. Average reward for episode 19.760000.
Episode: 600. Average reward for episode 20.600000.
Episode: 650. Average reward for episode 21.900000.
Episode: 700. Average reward for episode 23.380000.
Episode: 750. Average reward for episode 20.480000.
Episode: 800. Average reward for episode 26.620000.
Episode: 850. Average reward for episode 24.260000.
Episode: 900. Average reward for episode 27.740000.
Episode: 950. Average reward for episode 25.100000.
Episode: 1000. Average reward for episode 29.080000.
Episode: 1050. Average reward for episode 27.860000.
Episode: 1100. Average reward for episode 33.480000.
Episode: 1150. Average reward for episode 36.520000.
Episode: 1200. Average reward for episode 31.600000.
Episode: 1250. Average reward for episode 40.880000.
Episode: 1300. Average reward for episode 34.680000.
Episode: 1350. Average reward for episode 38.900000.
Episode: 1400. Average reward for episode 39.820000.
Episode: 1450. Average reward for episode 40.260000.
Episode: 1500. Average reward for episode 37.200000.
Episode: 1550. Average reward for episode 41.060000.
Episode: 1600. Average reward for episode 43.300000.
Episode: 1650. Average reward for episode 43.000000.
Episode: 1700. Average reward for episode 40.360000.
Episode: 1750. Average reward for episode 41.380000.
Episode: 1800. Average reward for episode 47.300000.
Episode: 1850. Average reward for episode 44.820000.
Episode: 1900. Average reward for episode 50.360000.
Episode: 1950. Average reward for episode 50.600000.
Episode: 2000. Average reward for episode 49.660000.
Episode: 2050. Average reward for episode 55.280000.
Episode: 2100. Average reward for episode 51.900000.
Episode: 2150. Average reward for episode 57.820000.
Episode: 2200. Average reward for episode 57.040000.
Episode: 2250. Average reward for episode 65.420000.
Episode: 2300. Average reward for episode 59.720000.
Episode: 2350. Average reward for episode 61.760000.
Episode: 2400. Average reward for episode 81.320000.
Episode: 2450. Average reward for episode 84.480000.
Episode: 2500. Average reward for episode 87.220000.
Episode: 2550. Average reward for episode 104.680000.
Episode: 2600. Average reward for episode 104.540000.
Episode: 2650. Average reward for episode 103.840000.
Episode: 2700. Average reward for episode 116.360000.
Episode: 2750. Average reward for episode 107.080000.
Episode: 2800. Average reward for episode 110.260000.
Episode: 2850. Average reward for episode 118.660000.
Episode: 2900. Average reward for episode 101.440000.
Episode: 2950. Average reward for episode 109.580000.
Episode: 3000. Average reward for episode 110.580000.
Episode: 3050. Average reward for episode 122.700000.
Episode: 3100. Average reward for episode 141.980000.
Episode: 3150. Average reward for episode 134.940000.
Episode: 3200. Average reward for episode 143.680000.
Episode: 3250. Average reward for episode 113.560000.
Episode: 3300. Average reward for episode 129.640000.
Episode: 3350. Average reward for episode 132.640000.
Episode: 3400. Average reward for episode 132.680000.
Episode: 3450. Average reward for episode 112.440000.
Episode: 3500. Average reward for episode 127.060000.
Episode: 3550. Average reward for episode 149.100000.
Episode: 3600. Average reward for episode 136.200000.
Episode: 3650. Average reward for episode 135.740000.
Episode: 3700. Average reward for episode 139.880000.
Episode: 3750. Average reward for episode 144.800000.
Episode: 3800. Average reward for episode 160.540000.
Episode: 3850. Average reward for episode 148.280000.
Episode: 3900. Average reward for episode 150.080000.
Episode: 3950. Average reward for episode 167.840000.
Episode: 4000. Average reward for episode 145.940000.
Episode: 4050. Average reward for episode 171.320000.
Episode: 4100. Average reward for episode 155.880000.
Episode: 4150. Average reward for episode 162.380000.
Episode: 4200. Average reward for episode 154.480000.
Episode: 4250. Average reward for episode 131.020000.
Episode: 4300. Average reward for episode 141.140000.
Episode: 4350. Average reward for episode 170.740000.
Episode: 4400. Average reward for episode 172.400000.
Episode: 4450. Average reward for episode 147.060000.
Episode: 4500. Average reward for episode 176.340000.
Episode: 4550. Average reward for episode 147.500000.
Episode: 4600. Average reward for episode 150.980000.
Episode: 4650. Average reward for episode 149.740000.
Episode: 4700. Average reward for episode 168.460000.
Episode: 4750. Average reward for episode 161.460000.
Episode: 4800. Average reward for episode 170.540000.
Episode: 4850. Average reward for episode 167.380000.
Episode: 4900. Average reward for episode 174.600000.
Episode: 4950. Average reward for episode 171.900000.
Episode: 5000. Average reward for episode 130.220000.
Episode: 5050. Average reward for episode 139.820000.
Episode: 5100. Average reward for episode 106.500000.
Episode: 5150. Average reward for episode 101.440000.
Episode: 5200. Average reward for episode 107.980000.
Episode: 5250. Average reward for episode 149.660000.
Episode: 5300. Average reward for episode 143.660000.
Episode: 5350. Average reward for episode 175.520000.
Episode: 5400. Average reward for episode 184.560000.
Episode: 5450. Average reward for episode 176.700000.
Episode: 5500. Average reward for episode 179.180000.
Episode: 5550. Average reward for episode 154.060000.
Episode: 5600. Average reward for episode 169.300000.
Episode: 5650. Average reward for episode 159.600000.
Episode: 5700. Average reward for episode 176.100000.
Episode: 5750. Average reward for episode 159.140000.
Episode: 5800. Average reward for episode 172.440000.
Episode: 5850. Average reward for episode 173.020000.
Episode: 5900. Average reward for episode 186.980000.
Episode: 5950. Average reward for episode 188.340000.
Episode: 6000. Average reward for episode 187.960000.
Episode: 6050. Average reward for episode 181.680000.
Episode: 6100. Average reward for episode 190.560000.
Episode: 6150. Average reward for episode 185.620000.
Episode: 6200. Average reward for episode 168.300000.
Episode: 6250. Average reward for episode 156.600000.
Episode: 6300. Average reward for episode 124.980000.
Episode: 6350. Average reward for episode 138.660000.
Episode: 6400. Average reward for episode 120.520000.
Episode: 6450. Average reward for episode 121.580000.
Episode: 6500. Average reward for episode 118.020000.
Episode: 6550. Average reward for episode 95.720000.
Episode: 6600. Average reward for episode 103.520000.
Episode: 6650. Average reward for episode 123.300000.
Episode: 6700. Average reward for episode 128.640000.
Episode: 6750. Average reward for episode 123.580000.
Episode: 6800. Average reward for episode 136.840000.
Episode: 6850. Average reward for episode 147.520000.
Episode: 6900. Average reward for episode 173.120000.
Episode: 6950. Average reward for episode 192.300000.
Episode: 7000. Average reward for episode 187.580000.
Episode: 7050. Average reward for episode 186.200000.
Episode: 7100. Average reward for episode 191.940000.
Episode: 7150. Average reward for episode 195.480000.
Episode: 7200. Average reward for episode 186.140000.
Episode: 7250. Average reward for episode 152.060000.
Episode: 7300. Average reward for episode 133.280000.
Episode: 7350. Average reward for episode 155.800000.
Episode: 7400. Average reward for episode 155.260000.
Episode: 7450. Average reward for episode 173.620000.
Episode: 7500. Average reward for episode 190.020000.
Episode: 7550. Average reward for episode 166.600000.
Episode: 7600. Average reward for episode 139.360000.
Episode: 7650. Average reward for episode 137.860000.
Episode: 7700. Average reward for episode 132.740000.
Episode: 7750. Average reward for episode 143.900000.
Episode: 7800. Average reward for episode 123.340000.
Episode: 7850. Average reward for episode 146.860000.
Episode: 7900. Average reward for episode 184.620000.
Episode: 7950. Average reward for episode 171.320000.
Episode: 8000. Average reward for episode 140.620000.
Episode: 8050. Average reward for episode 175.820000.
Episode: 8100. Average reward for episode 193.440000.
Episode: 8150. Average reward for episode 190.380000.
Episode: 8200. Average reward for episode 190.440000.
Episode: 8250. Average reward for episode 196.720000.
Episode: 8300. Average reward for episode 196.180000.
Episode: 8350. Average reward for episode 190.640000.
Episode: 8400. Average reward for episode 186.120000.
Episode: 8450. Average reward for episode 198.800000.
Episode: 8500. Average reward for episode 197.900000.
Episode: 8550. Average reward for episode 198.720000.
Episode: 8600. Average reward for episode 195.440000.
Episode: 8650. Average reward for episode 197.820000.
Episode: 8700. Average reward for episode 171.900000.
Episode: 8750. Average reward for episode 150.220000.
Episode: 8800. Average reward for episode 159.620000.
Episode: 8850. Average reward for episode 138.580000.
Episode: 8900. Average reward for episode 137.540000.
Episode: 8950. Average reward for episode 122.820000.
Episode: 9000. Average reward for episode 136.600000.
Episode: 9050. Average reward for episode 149.400000.
Episode: 9100. Average reward for episode 154.920000.
Episode: 9150. Average reward for episode 129.560000.
Episode: 9200. Average reward for episode 163.900000.
Episode: 9250. Average reward for episode 131.200000.
Episode: 9300. Average reward for episode 142.480000.
Episode: 9350. Average reward for episode 131.220000.
Episode: 9400. Average reward for episode 132.500000.
Episode: 9450. Average reward for episode 160.940000.
Episode: 9500. Average reward for episode 182.020000.
Episode: 9550. Average reward for episode 198.180000.
Episode: 9600. Average reward for episode 198.060000.
Episode: 9650. Average reward for episode 200.000000.
Episode: 9700. Average reward for episode 200.000000.
Episode: 9750. Average reward for episode 199.840000.
Episode: 9800. Average reward for episode 198.140000.
Episode: 9850. Average reward for episode 194.600000.
Episode: 9900. Average reward for episode 197.020000.
Episode: 9950. Average reward for episode 196.220000.
Episode: 10000. Average reward for episode 195.840000.
```

**ソリューション**

```python:In_[21]==>
observation = C.sequence.input_variable(STATE_COUNT, np.float32, name="s")

W1 = C.parameter(shape=(STATE_COUNT, H), init=C.glorot_uniform(), name="W1")
b1 = C.parameter(shape=H, name="b1")
layer1 = C.relu(C.times(observation, W1) + b1)

W2 = C.parameter(shape=(H, ACTION_COUNT), init=C.glorot_uniform(), name="W2")
b2 = C.parameter(shape=ACTION_COUNT, name="b2")
model = C.times(layer1, W2) + b2
W1.shape, b1.shape, W2.shape, b2.shape, model.shape
```
```bash
((4, 10), (10,), (10, 2), (2,), (2,))
```

```python:In_[22]==>
# Correspoding layers implementation - Preferred solution
def create_model(input):
    with C.layers.default_options(init=C.glorot_uniform()):
        z = C.layers.Sequential([C.layers.Dense(H, name="layer1"),
                                 C.layers.Dense(ACTION_COUNT, name="layer2")])
        return z(input)
    
model = create_model(observation)

model.layer1.W.shape, model.layer1.b.shape, model.layer2.W.shape, model.layer2.b.shape, model.shape
```
```bash
((4, 10), (10,), (10, 2), (2,), (2,))
```

**以上**