化学系エンジニアがAIを学ぶ

具体的にはPyTorch、深層強化学習を学び、主に化学工学の問題に取り組みます

メモ: PyTorch 重み、バイアスの初期化について

はじめに

PyTorchのニューラルネットワークの重み・バイアスの初期化についてのメモを記す。

重み

重みの内容は次のようにして確認できる。

>>> import torch.nn as nn
>>> l = nn.Linear(1, 3)
>>> l.weight
Parameter containing:
tensor([[ 0.6204],
        [-0.5651],
        [-0.6809]], requires_grad=True)

重みの初期化は次のようにnn.initモジュールの関数を使用することで実施できる。

>>> nn.init.xavier_uniform_(l.weight) # Xavierの初期値
Parameter containing:
tensor([[ 0.5318],
        [-0.8438],
        [-1.1093]], requires_grad=True)

>>> nn.init.kaiming_uniform_(l.weight) # Heの初期値
Parameter containing:
tensor([[ 2.3349],
        [-0.9218],
        [-1.3357]], requires_grad=True)

単純な正規分布・一様分布や定数を設定できる関数もnn.init内に準備されている。

詳しくはこちらを参照: https://pytorch.org/docs/master/nn.init.html

バイアス

バイアスの内容は重みと同様に確認できる。

>>> l.bias
Parameter containing:
tensor([-0.9920,  0.6406, -0.6558], requires_grad=True)

バイアスの初期化も重みと同様のやり方。

>>> nn.init.uniform_(l.bias, 0, 1) # 一様分布
Parameter containing:
tensor([0.8032, 0.7810, 0.1619], requires_grad=True)

>>> nn.init.constant_(l.bias, 0) # 定数
Parameter containing:
tensor([0., 0., 0.], requires_grad=True)

デフォルトの初期値

nn.Linearのソースを見ると重みについてはHeの初期値が入っている模様。

pytorch/linear.py at master · pytorch/pytorch · GitHub

def reset_parameters(self):
    init.kaiming_uniform_(self.weight, a=math.sqrt(5))
    if self.bias is not None:
        fan_in, _ = init._calculate_fan_in_and_fan_out(self.weight)
        bound = 1 / math.sqrt(fan_in)
        init.uniform_(self.bias, -bound, bound)

参考

まるまるにっき | pytorchの重みの初期化の方法

python - How to initialize weights in PyTorch? - Stack Overflow

What is the default initialization of a conv2d layer and linear layer? - PyTorch Forums

ニュース: Preferred Networks、JXTGエネルギーと石油精製プラントの最適化・自動化に関する共同研究を開始

Preferred NetworksがJXTGエネルギーと石油精製プラントの最適化・自動化に関する共同研究を開始したそうです。10億円もの資金が調達されるそうで、一気に技術開発が進むことが期待されます。近い業界に関わる者としては、得られた知見などを可能な範囲でオープンにしてもらえたらと思いますが、難しいでしょうね。。。

it.impressbm.co.jp

iotnews.jp

深層強化学習(Deep Q Network, DQN)の簡単な例 〜Experience Replay追加〜

はじめに

前回の記事でOpenAI Gymを使わず非常に簡単な問題を対象にDQNを適用してみたが、"Experience Replay"を入れていなかった。今回は前回の問題にExperience Replayを追加してみる。なおこれを実施するにあたり、下記サイトを参考にした。

第15回 CartPole課題で深層強化学習DQNを実装|Tech Book Zone Manatee

Reinforcement Learning (DQN) Tutorial — PyTorch Tutorials 1.1.0 documentation

Experience Replay

行動とその結果を経験として記録しておき、その過去の経験をサンプリングして学習データとする方法。通常の逐次で学習する方法に比べ、時間的に離れたデータを使用することができ、偏りを抑えた安定した学習ができるとのこと。複数データをまとめてニューラルネットワークをミニバッチ学習させることもできる。

具体的には「現在の状態」、「行動」、「行動後の状態」、「報酬」のデータの組をある一定の数だけメモリに残しておいて、その中から学習に使うデータをいくつかランダムにサンプリングし、これらを1つのバッチとしてニューラルネットワークに学習させる。

対象とする問題

前回記事と同一の報酬払出装置を考える。プログラムは前回のものにいくつか追加・修正してExperience Replayを適用する。

深層強化学習(Deep Q Network, DQN)の簡単な例 - 化学系エンジニアがAIを学ぶ

Experience Replayの実装

下記のpytorch.orgのチュートリアルにあるように、namedtupleを用いて経験を記録するクラスを追加する。

Reinforcement Learning (DQN) Tutorial — PyTorch Tutorials 1.1.0 documentation

ここでは500回分の行動の結果を記録することにしている(数は適当)。500を超えると古い方から順番に上書きされる形となっている。サンプリング時は指定したバッチサイズ数だけのデータが取り出される。

from collections import namedtuple

# 「現在の状態」、「行動」、「行動後の状態」、「報酬」を記録するnamedtuple
Transition = namedtuple('Transition', ('state', 'action', 'next_state', 'reward'))

import random

class ReplayMemory(object):
    def __init__(self, capacity):
        self.capacity = capacity # メモリの許容サイズ
        self.memory = []         # これまでの経験を保持
        self.index = 0           # 保存した経験のindex

    def push(self, state, action, next_state, reward):
        if len(self.memory) < self.capacity: # 許容サイズ以下なら追加
            self.memory.append(None)
        self.memory[self.index] = Transition(state, action, next_state, reward)
        self.index = (self.index + 1) % self.capacity # indexを1つ進める。許容サイズ超えたら古い順から上書き

    def sample(self, batch_size):
        return random.sample(self.memory, batch_size) # バッチサイズ分のランダムサンプリング

    def __len__(self):
        return len(self.memory)

CAPACITY = 500

memory = ReplayMemory(CAPACITY)

DQNの学習

上記の行動を記録したクラスを利用できる形に、DQNの学習部分を修正する。

BATCH_SIZE = 50 # サンプリングするデータ数

def update_dqn(replay_memory): # 行動を記録したクラスを引数とする
    ## メモリがバッチサイズより小さいときはまだ学習しない
    if len(replay_memory) < BATCH_SIZE:
        return

    ## バッチ取得
    transitions = replay_memory.sample(BATCH_SIZE)
    ## (状態、行動、次の状態、報酬)✕バッチサイズ を (状態xバッチサイズ、行動✕バッチサイズ、、、)に変える
    batch = Transition(*zip(*transitions))
    
    ## 各値をtensorに変換
    state = torch.FloatTensor(batch.state).unsqueeze(1)
    action = torch.LongTensor(batch.action).unsqueeze(1)
    next_state = torch.FloatTensor(batch.next_state).unsqueeze(1)
    reward = torch.FloatTensor(batch.reward).unsqueeze(1)

    ## Q値の算出
    q_now = dqn(state).gather(-1, action) # 今の状態のQ値
    max_q_next = dqn(next_state).max(-1)[0].unsqueeze(1).detach() # 状態移行後の最大のQ値
    gamma = 0.9
    q_target = reward + gamma * max_q_next # 目標のQ値

    # DQNパラメータ更新
    optimizer.zero_grad()
    loss = criterion(q_now, q_target) # 今のQ値と目標のQ値で誤差を取る
    loss.backward()
    optimizer.step()
    
    return loss.item() # lossの確認のために返す

学習結果

他のコードは前回記事と同じである。下図に学習の推移を示す。左は誤差の推移、右は1エピソードにおけるトータル報酬の推移である。実際に計算してみるとわかるが逐次で学習させるより、Experience Replayを適用するほうが学習が安定しやすい。

f:id:schemer1341:20190609202605p:plain

コード

chemical-engineer-learns-ai/dispenser_experience_replay.py at master · nakamura-13/chemical-engineer-learns-ai · GitHub

参考

第15回 CartPole課題で深層強化学習DQNを実装|Tech Book Zone Manatee

Reinforcement Learning (DQN) Tutorial — PyTorch Tutorials 1.1.0 documentation

深層強化学習(Deep Q Network, DQN)の簡単な例

はじめに

DQNを学ぼうとして色々と調べたが、どこもかしこもOpenAI Gymを使っていて、まずそれの扱いから考えないといけないのでつらい。ここではOpenAI Gymを使わず、非常に簡単な問題を対象にDQNを適用してみることとする。Python、PyTorchを用いる。

また、DQNをやろうとすると出てくる"Experience Replay"。これも入れようとすると考えるのがつらいので、入れない。

下記記事にDQN全般やExperience Replayについてちょっと解説がある。

Deep Q Network (DQN) - DeepLearningを勉強する人

対象とする問題

下図のような報酬払出装置を考える。装置には「電源」ボタンと「払出」ボタンが設置されている。「電源」ボタンを押すと装置電源のON/OFFが切り替わる。電源ONのときに「払出」ボタンを押すと報酬が払い出される。電源OFFのときに「払出」ボタンを押しても報酬は払い出されない。

f:id:schemer1341:20190502015126p:plain:w250

1回の行動で、いずれかのボタンを1度押すことができるものとする。5回の行動を行う場合に、報酬を最大化する手順を学習することができるか?

なお、この問題の内容に関して下記書籍を参考にした。

この払出装置(Dispenser)をプログラムにすると次の通りとなる。

class Dispenser(object):
    def __init__(self, init_state):
        """
        初期のON/OFF状態を設定する
        init_state: 0->電源OFF、1->電源ON
        """
        self.state  = init_state

    def powerbutton(self):
        """
        電源ボタンを押し、ON/OFFを切り替える
        """
        if self.state == 0:
            self.state = 1
        else:
            self.state = 0

    def step(self, action):
        """
        払出機を操作する
        action: 0->電源ボタンを押す 1->払出ボタンを押す
        状態と報酬が返る
        """
        if action == 0: # 電源ボタンを押した場合
            self.powerbutton() # 電源ON/OFF切り替え
            reward = 0 # 報酬はない
        else:           # 払出ボタンを押した場合
            if self.state == 1:
                reward = 1 # 電源ONのときのみ報酬あり
            else:
                reward = 0
        return self.state, reward

Deep Q Network (DQN)

報酬払出機の状態は2つ、取れる行動は2つであるので、これを表すQテーブルは次のようになる。

f:id:schemer1341:20190502090108p:plain:w400

Q学習ではこのQテーブルのQ値(Q(0,0), Q(0,1), Q(1,0), Q(1,1))を逐次更新していくが、DQNではこれらのQ値を得る関数を逐次更新していく。その関数を表すのにニューラルネットワークを用いるが、次のようなイメージとなる。

f:id:schemer1341:20190502141215p:plain:w550

状態を入力として、Q値を出力として得る形である。ここでは上図の通りのニューラルネットワークを用いることとする。(あまりディープではないが・・・)

import torch
import torch.nn as nn
import torch.nn.functional as F

class DQN(nn.Module):
    def __init__(self):
        super(DQN, self).__init__()
        self.l1 = nn.Linear(1, 3)
        self.l2 = nn.Linear(3, 3)
        self.l3 = nn.Linear(3, 2)

    def forward(self, x):
        x = F.relu(self.l1(x))
        x = F.relu(self.l2(x))
        x = self.l3(x)
        return x

dqn = DQN()
optimizer = torch.optim.SGD(dqn.parameters(), lr=0.01)
criterion = nn.MSELoss()

DQNの学習方法

では、どのような情報を使ってこのニューラルネットワークを学習させていくのか、ということになるが、これにはQ学習のQ値の更新式を用いる。

f:id:schemer1341:20190502150605p:plain:w500

この式は今のQ値( Q(s_t, a_t))を、今の推定の目標値( r_{t+1}+\gamma {\rm max}Q)に近づくよう更新するものである。すなわち、

f:id:schemer1341:20190502155057p:plain:w290

となることを目指している。よって Q(s_t, a_t) r_{t+1}+\gamma {\rm max}Qの誤差を小さくするようにニューラルネットワークを学習させるという形にすればよいことになる。いずれの値も現在の状態、行動、次の状態、報酬が分かれば算出可能であり、DQNの学習をプログラムにすると次のようになる。

def update_dqn(state, action, next_state, reward):
    ## 各変数をtensorに変換
    state = torch.FloatTensor([state])
    action = torch.LongTensor([action]) # indexとして使うのでLong型
    next_state = torch.FloatTensor([next_state])

    ## Q値の算出
    q_now = dqn(state).gather(-1, action) # 今の状態のQ値
    max_q_next = dqn(next_state).max(-1)[0].detach() # 状態移行後の最大のQ値
    gamma = 0.9
    q_target = reward + gamma * max_q_next # 目標のQ値

    # DQNパラメータ更新
    optimizer.zero_grad()
    loss = criterion(q_now, q_target) # 今のQ値と目標のQ値で誤差を取る
    loss.backward()
    optimizer.step()
    
    return loss.item() # lossの確認のために返す

DQNからの行動の決定

DQNは入力された状態に対し、取りうるすべての行動のQ値のtensorを返す(1次元tensor)。そこから最大のQ値のindexを取り出して(0または1)、それをそのまま行動としている。

import numpy as np

EPS_START = 0.9
EPS_END = 0.0
EPS_DECAY = 200

def decide_action(state, episode):
    state = torch.FloatTensor([state]) # 状態を1次元tensorに変換
    ## ε-グリーディー法
    eps = EPS_END + (EPS_START - EPS_END) * np.exp(-episode / EPS_DECAY)
    if eps <= np.random.uniform(0, 1):
        with torch.no_grad():
            action = dqn(state).max(-1)[1] # Q値が最大のindexが得られる
            action = action.item() # 0次元tensorを通常の数値に変換
    else:
        num_actions = len(dqn(state)) # action数取得
        action = np.random.choice(np.arange(num_actions)) # ランダム行動
    return action

学習する

現在の状態からDQNより行動を決めて、環境から次の状態・報酬を取得し、これらの結果を用いてDQNを更新する、を繰り返す。エピソードの数や上記のoptimizerの学習率、ε-グリーディー法のepsの変化のさせ方はうまく学習が進むように調節して決めている。(適度な値を決めるのにかなり苦労した。)

import matplotlib.pyplot as plt

NUM_EPISODES = 1200
NUM_STEPS = 5

log = [] # 結果のプロット用

for episode in range(NUM_EPISODES):
    env = Dispenser(0)
    total_reward = 0 # 1エピソードでの報酬の合計を保持する

    for s in range(NUM_STEPS):
        ## 現在の状態を確認
        state = env.state
        ## 行動を決める
        action = decide_action(state, episode)
        ## 決めた行動に従いステップを進める。次の状態、報酬を得る
        next_state, reward = env.step(action)
        ## DQNを更新
        loss = update_dqn(state, action, next_state, reward)
        total_reward += reward

    log.append([total_reward, loss])

## 結果表示
r, l = np.array(log).T
fig = plt.figure(figsize=(11,4))
ax1 = fig.add_subplot(121)
ax2 = fig.add_subplot(122)
ax1.set_xlabel("Episode")
ax1.set_ylabel("Loss")
ax1.set_yscale("log")
ax2.set_xlabel("Episode")
ax2.set_ylabel("Total reward")
ax1.plot(l)
ax2.plot(r)
plt.show()

学習結果

下図に学習の推移を示す。左は誤差の推移、右は1エピソードにおけるトータル報酬の推移である。払出機は初め電源OFFのため、5回の行動で得られる最大のトータル報酬は4であるが、右下のグラフから最終的にトータル報酬4が得られる行動を適切に学習できていることを確認できる。

f:id:schemer1341:20190503161744p:plain

学習後のDQNから各状態におけるQ値を確認できるが、今回の学習での値は下記である。

>>> dqn(torch.FloatTensor([0]))
tensor([9.0000, 8.1200], grad_fn=<AddBackward0>)
>>> dqn(torch.FloatTensor([1]))
tensor([8.1060, 10.0000], grad_fn=<AddBackward0>)

計算の詳細は示さないが、今回と同一の問題にQ学習を適用するとこれらとほぼ同じQ値が得られ、DQNによってQ値を適切に表現する関数を構築できていることを確認できる。

さいごに

DQNを適用するのが全く意味がないような簡単な問題にDQNを適用した。"Experience Replay"も考えなかった(状態✕行動の組み合わせが4パターンしかない問題なのであまり効果はないのかも)。簡単すぎる問題であるがDQNを適用する上での基礎的なポイントは掴めたのではないかと思う。

コード

ここで示した全コードおよびQ学習を適用した場合のコードはこちら:

chemical-engineer-learns-ai/20190504_simple_example_of_DQN at master · nakamura-13/chemical-engineer-learns-ai · GitHub

参考

Reinforcement Learning (DQN) Tutorial — PyTorch Tutorials 1.1.0 documentation

・牧野浩二、西崎博光, 『Pythonによる深層強化学習入門』, オーム社 (2018)

反応器計算におけるQ学習2

はじめに

以前の記事で管型反応器において所望の反応率を得るという問題について強化学習(Q学習)を適用したが、 今回は反応率に応じて目的成分の選択率が変化するという条件で、収率を最大化する問題を対象とする。

管型反応器の例題

対象とする問題は下記の通り。

  • 成分Aを反応して目的の生成物Bを得る管型反応器がある
  • 反応器は連続運転しており、1日に1回だけAの反応率とBの収率を確認でき、その際に反応温度を変更することができる
  • 反応温度は1回につき0.5℃変更できる(上げることも、下げることも、変えないことも可能)
  • 反応温度を変更することで反応率を変化させることができる
  • 生成物Bと一緒に副生物Cも生成し、反応率に応じて生成物Bの選択率は変化する

ある初期温度で運転を開始したのちに、収率を最大化することを学習できるか?

f:id:schemer1341:20190430211358p:plain:w400

ここで反応温度・反応率・選択率・収率の関係は図中の式で表されるものとし、 本例題ではk=0.01、a=100とする。 反応温度と反応率、反応率と選択率の関係は下記グラフとなる。

f:id:schemer1341:20190430213228p:plain

なお、収率が最大となる条件はこの時点で計算可能であり、反応率0.955、収率0.945である。 これを強化学習で推算できるかがポイントとなる。

管型反応器のモデル

反応器の状態と行動を次の表の通り整理する。 状態は反応率を離散化して表現する。ここでは反応率0.9未満を状態0とし、0.9から1.0を100分割して番号を振った。 (収率が最大となる反応率が少なくとも0.9以上とわかっている前提)

f:id:schemer1341:20190430215220p:plain:w500

報酬は-(1-収率)とし、1に満たない分をマイナスの報酬(ペナルティ)として与える形とした。

import numpy as np
import matplotlib.pyplot as plt

class Tubular_selec(object):
    def __init__(self, temp):
        self._k = 0.01  # kを0.01とする
        self._a = 100.0 # aを100とする
        init_conv  = self.calc_conv(self._k, temp)
        init_yield = self.calc_yield(self._a, init_conv)
        ## ログとして[温度, 反応率, 収率]を残す
        self._log = [[temp, init_conv, init_yield]]

    def calc_conv(self, k, temp):
        ## 反応率計算
        return 1 - np.exp(-k * temp)

    def calc_selc(self, a, conv):
        ## 選択率計算
        return 1 - conv**a

    def calc_yield(self, a, conv):
        ## 収率計算
        return conv * self.calc_selc(a, conv)

    def state(self):
        ## 状態を返す 状態を転化率0.9-1.0の間で100分割する
        ## 0.9未満を状態0とするので、状態は全部で101個
        _, conv, _ = self._log[-1] # 一番最新のlogからデータ取得
        n = 100
        if conv < 0.9: # 反応率 0.1未満は状態0とする
            st = 0
        else:
            st = int(np.ceil((conv - 0.9)*10*n)) # 0.9~1.0をn分割
        return st

    def step(self, action):
        temp, _, _ = self._log[-1] # 最新の温度を取得
        ## 行動0で温度UP、行動2で温度Down、行動1は何もしない
        if action == 0:
            temp += 0.5
        elif action == 2:
            temp -= 0.5
        conv = self.calc_conv(self._k, temp)
        yild = self.calc_yield(self._a, conv)
        self._log.append([temp, conv, yild])
        # 1.0に満たない分を負の報酬(ペナルティ)として与える
        reward = - abs(1.0 - yild)
        return reward

Q値に従った次の行動の決定

def decide_action(qtable, state, episode):
    eps = 1 - episode/240 # ε-グリーディー法
    if eps <= np.random.uniform(0, 1): # グリーディー行動
        t = np.where(qtable[state]==qtable[state].max())[0] # Q値が最大のインデックスを返す
    else:                              # ランダム行動
        t = np.arange(qtable.shape[-1]) # 取りうる行動すべて
    return np.random.choice(t) # 行動の候補からランダムに選ぶ

Q値の更新

def update_qtable(qtable, state, action, reward, next_state):
    gamma = 0.9
    alpha = 0.5
    next_qmax = max(qtable[next_state])
    qtable[state, action] = (1 - alpha) * qtable[state, action] + alpha * (
        reward + gamma * next_qmax)
    return qtable

学習する

1episodeあたり100日運転することとして、300episode学習させる。

OPERATION_PERIOD = 100 # 100日分
EPISODE_NUM = 300
qtable = np.zeros((101, 3)) # サイズは状態数✕行動数

## 学習状況確認用グラフ設定
fig = plt.figure(figsize=(11,4))
ax1 = fig.add_subplot(121)
ax2 = fig.add_subplot(122)

## 学習
for episode in range(EPISODE_NUM):
    #t = np.random.uniform(280, 350) # 初期温度をランダムにする場合
    env = Tubular_selec(280)
    
    for d in range(OPERATION_PERIOD):
        ## まず現在の状態を確認
        state = env.state()
        ## 次に行動を決める
        action = decide_action(qtable, state, episode)
        ## 決めた行動に従い運転を1日進める.また報酬を環境から得る
        reward = env.step(action)
        ## 新しい状態を確認
        next_state = env.state()
        ## Qテーブル更新
        qtable = update_qtable(qtable, state, action, reward, next_state)
    ## 10episodeにつき1回のグラフ描画
    if (episode+1)%10 == 0:
        t, c, y = np.array(env._log).T
        ax1.cla()
        ax2.cla()
        ax1.set_xlabel("Days")
        ax1.set_ylabel("Temperature")
        ax1.set_ylim(270,320)
        ax2.set_xlabel("Days")
        ax2.set_ylabel("Conversion, Yield")
        ax2.set_ylim(0.93,0.96)
        ax1.plot(t)
        ax2.plot(c, label="Conversion")
        ax2.plot(y, label="Yield")
        ax2.legend()
        plt.pause(0.1)

学習結果

収率最大となるのは反応率0.955、収率0.945のときであるが、 反応温度を上げて、まっすぐその状態に近づいていることが確認できる。

f:id:schemer1341:20190430221901p:plain

反応器計算におけるQ学習

はじめに

簡単な例として、管型反応器において反応温度を調節して所望の反応率を得るという問題について強化学習(Q学習)を適用してみる。

管型反応器の例題

対象とする問題は下記の通り。

  • 成分Aを反応して生成物Bを得る管型反応器がある
  • 反応器は連続運転しており、1日に1回だけ反応率を確認でき、その際に反応温度を変更することができる
  • 反応温度は1回につき0.5℃変更できる(上げることも、下げることも、変えないことも可能)
  • 反応温度を変更することで反応率を変化させることができる

ある初期温度で運転を開始したのちに所望の反応率に調整することを学習できるか。

f:id:schemer1341:20190427154711p:plain:w400

ここで反応率と反応温度の関係は図中の式で表されるものとし、本例題ではk=0.01とする。 反応率と温度は下記グラフの関係となる。

f:id:schemer1341:20190427144449p:plain:w400

管型反応器のモデル

反応器の状態と行動を次の表の通り整理する。ここでは反応率の目標を0.9599〜0.96としている。 状態は目標より低い、同じ、高い、の3つとし、行動は温度を上げる、何もしない、温度を下げる、の3つとする。 反応率が目標と同じで、かつ何もしないときに報酬が得られるものとする。

f:id:schemer1341:20190427153801j:plain:w400
import numpy as np
import matplotlib.pyplot as plt

class Tubular(object):
    def __init__(self, temp):
        self._k = 0.01 # kを0.01とする
        ## ログとして[反応温度, 反応率]を残す
        self._log = [[temp, self.calc_conv(self._k, temp)]]

    def calc_conv(self, k, temp):
        ## 反応率を計算
        return 1 - np.exp(-k * temp)

    def state(self):
        ## 状態を返す 反応率が目標内なら1、それより下は0、上は2
        _, conv = self._log[-1] # 一番最新のlogからデータ取得
        if conv < 0.9599:
            st = 0
        elif conv > 0.9600:
            st = 2
        else:
            st = 1
        return st

    def step(self, action):
        ## 計算を1日分進めるステップ計算
        ## まず報酬計算  状態1で行動1(何もしない)のとき報酬あり
        if self.state()==1 and action==1:
            reward = 1
        else:
            reward = 0
        ## 1日分進める
        temp, _  = self._log[-1] # 一番最新のlogからデータ取得
        if action == 0:
            temp += 0.5
        elif action == 2:
            temp -= 0.5
        conv = self.calc_conv(self._k, temp)
        self._log.append([temp, conv])
        return reward

Q値に従った次の行動の決定

def decide_action(qtable, state, episode):
    eps = 0.8 - 0.8*episode/800 # ε-グリーディー法
    if eps <= np.random.uniform(0, 1):
        ## Q値が最大のインデックスを返す
        t = np.where(qtable[state]==qtable[state].max())[0]
    else:
        ## 取りうる行動すべてを返す
        t = np.arange(qtable.shape[-1])
    return np.random.choice(t) # 行動の候補からランダムに選ぶ

Q値の更新

def update_qtable(qtable, state, action, reward, next_state):
    gamma = 0.9
    alpha = 0.5
    next_qmax = max(qtable[next_state])
    qtable[state, action] = (1 - alpha) * qtable[state, action] + alpha * (
        reward + gamma * next_qmax)
    return qtable

学習する

1episodeあたり100日運転することとして、1000episode学習させる。 各episodeの初期の反応温度はランダムに決める。

operation_period = 100 # 100日分
episode_num = 1000
qtable = np.zeros((3, 3)) # サイズは状態数✕行動数

for episode in range(episode_num):
    t = np.random.uniform(300, 350)
    env = Tubular(t)

    for d in range(operation_period):
        ## まず現在の状態を確認
        state = env.state()
        ## 次に行動を決める
        action = decide_action(qtable, state, episode)
        ## 決めた行動に従い運転を1日進める.また報酬を環境から得る
        reward = env.step(action)
        ## 新しい状態を確認
        next_state = env.state()
        ## Qテーブル更新
        qtable = update_qtable(qtable, state, action, reward, next_state)

学習結果

Qテーブルの一例。反応率が目標未満なら温度を上げ、目標より高いなら温度を下げるような値になっている。

#          温度上げる       何もしない      温度下げる
array([[1.97888446e-02, 3.12142611e-12, 1.27102555e-18], # 反応率目標未満
       [2.34887193e+00, 1.00000000e+01, 5.69672696e+00], # 目標内
       [7.11430684e-27, 7.22125479e-27, 1.91366261e-02]]). # 目標より高い

次の結果確認用コードを用いて各種初期温度を計算すると、いずれも最終的に所望の反応率に到達していることを確認できる。

def suisan(init_temp):
    env = Tubular(init_temp)
    for d in range(operation_period):
        state = env.state()
        action = decide_action(qtable, state, 10000) # episode数は適当に大きい値を入れる。探索をしないため。
        _ = env.step(action)
    t, c = np.array(env._log).T
    fig = plt.figure(figsize=(11,4))
    ax1 = fig.add_subplot(121)
    ax2 = fig.add_subplot(122)
    ax1.set_xlabel("Days")
    ax2.set_xlabel("Days")
    ax1.set_ylabel("Temperature")
    ax2.set_ylabel("Conversion")
    ax1.plot(t)
    ax2.plot(c)
    plt.show()

初期温度 310

f:id:schemer1341:20190427161633p:plain

初期温度 340

f:id:schemer1341:20190427161650p:plain

メモ: PyTorch tensor requires_gradのTrue/False確認、切り替え

requires_gradのTrue/False確認

属性 requires_gradを参照して確認できる。

import torch

a = torch.tensor([0, 1], dtype=torch.float32)
b = torch.tensor([2, 3], dtype=torch.float32, requires_grad=True)

a.requires_grad
b.requires_grad

実行結果

False
True

requires_gradのTrue/False切り替え

メソッド requires_grad()で切り替えられる

a.requires_grad_(True)
a.requires_grad
True # <- FalseからTrueに変わった