反応器計算におけるQ学習
はじめに
簡単な例として、管型反応器において反応温度を調節して所望の反応率を得るという問題について強化学習（Q学習）を適用してみる。
管型反応器の例題
対象とする問題は下記の通り。
- 成分Aを反応して生成物Bを得る管型反応器がある
- 反応器は連続運転しており、1日に1回だけ反応率を確認でき、その際に反応温度を変更することができる
- 反応温度は1回につき0.5℃変更できる（上げることも、下げることも、変えないことも可能）
- 反応温度を変更することで反応率を変化させることができる
ある初期温度で運転を開始したのちに所望の反応率に調整することを学習できるか。
ここで反応率と反応温度の関係は図中の式で表されるものとし、本例題ではk=0.01とする。 反応率と温度は下記グラフの関係となる。
管型反応器のモデル
反応器の状態と行動を次の表の通り整理する。ここでは反応率の目標を0.9599〜0.96としている。 状態は目標より低い、同じ、高い、の3つとし、行動は温度を上げる、何もしない、温度を下げる、の３つとする。 反応率が目標と同じで、かつ何もしないときに報酬が得られるものとする。
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
初期温度 340
メモ: 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に変わった
メモ: PyTorchの自動微分について
はじめに
PyTorchの自動微分についてのメモを記す。
自動微分
アルゴリズムによって定義された関数からその関数の偏導関数値を計算するアルゴリズムを導出するための技術。 数値微分や数式微分とも異なる方法で、これら方法の問題点（誤差や桁落ちによる精度低下、プログラムから微分表現を導けない）を 解消できる手段とのこと。
PyTorchでの自動微分の例
- tensorを作るときに、
requires_grad=Trueとする。 - テンソル演算する。 (ex. y = x **2)
y.backward()で勾配を計算する。x.gradで勾配値が得られる
import torch x = torch.tensor(1.0, requires_grad=True) # require_grad=Trueが必要。defaultはFalse y = x ** 2 y.backward() # 勾配を求める print(x.grad) # xにおける勾配の値がx.gradに入っている x = torch.tensor(1.0, requires_grad=True) y = 2 ** x y.backward() print(x.grad) x1 = torch.tensor(1.0, requires_grad=True) x2 = torch.tensor(1.0, requires_grad=True) y = x1**2 + torch.sqrt(x1*x2) # 偏微分のケース y.backward() print(x1.grad, x2.grad)
実行結果
tensor(2.) tensor(1.3863) tensor(2.5000) tensor(0.5000)
参考
スキナー箱のQ学習
はじめに
ネズミがボタンを押して餌を得るという問題について強化学習（Q学習）を適用してみる。 本記事の内容に関して下記書籍を参考にした。
ネズミの学習問題
対象とする問題は下記の通り。
- 箱の中にネズミが1匹いる
- 箱の中にネズミ餌のディスペンサーが1つある
- ディスペンサーには「電源ボタン」と「餌ボタン」がある
- 「電源ボタン」をネズミが押すたびに電源のON/OFFが切り替わる
- ディスペンサーの電源がONのときにネズミが「餌ボタン」を押すと餌が出てくる
ネズミは餌を得る手順を学習できるか。
ディスペンサーのモデル
まず、電源のON/OFF状態を保持し、ON状態で「餌ボタン」が押されると餌（報酬）を返すクラスを作成する。
import numpy class Dispenser(object): def __init__(self, init_state): """ 初期のON/OFF状態を設定する init_state: 0->電源OFF、1->電源ON """ self._state = init_state self._button = self.powerbutton(init_state) def powerbutton(self, init_state): """ ON/OFFを切り替えるgenerator 0: 電源OFF、1: 電源ON を返す。内部的には-1, 1で切り替えている """ state = -1 if init_state==0 else init_state while True: state = -1 * state yield numpy.maximum(0, state) def push_powerbutton(self): self._state = self._button.__next__() def step(self, action): """ action: 0->電源ボタンを押す 1->餌ボタンを押す """ if action == 0: self.push_powerbutton() reward = 0 else: reward = 1 if self._state==1 else 0 return self._state, reward
挙動としては下記。電源ボタンを押すとON/OFFが切り替わり、ONのとき餌ボタンを押すと餌(報酬 1)が得られる。 具体的には(状態, 報酬)のタプルが返る。これだけ。
>>> env = Dispenser(0) # 電源OFFを初期状態としてディスペンサーオブジェクトを作成 >>> env.step(0) # 0:電源ボタンを押す (1, 0) # 状態:電源ON, 報酬:0 が返る >>> env.step(1) # 1: 餌ボタンを押す (1, 1) # 状態:電源ON, 報酬:1が返る
Q値
各状態および行動に対するQ値の表を次の通り決める。
Q(0, 0): 状態 電源OFF、行動 電源ボタンを押す
Q(0, 1): 状態 電源OFF、行動 餌ボタンを押す
Q(1, 0): 状態 電源ON、 行動 電源ボタンを押す
Q(1, 1): 状態 電源ON、 行動 餌ボタンを押す
各Q値は下記の通り配列で保持することとする。初期値は全て0。
qtable = numpy.zeros((2, 2))
Q値に従って次の行動を決める
Q値および現在の状態から次に取る行動を決める。基本的にはQ値が最大の行動を取るが、
時折ランダムな行動を取る - グリーディー法を取り入れることとする。
def get_action(qtable, state, episode): """ 現在の状態とQ値から次の行動を決める episode回数が増えるごとにランダム行動を取らなくする """ ## ε- greedy法 行動を重ねる(episodeが増える)毎にランダム行動を取らなくする epsilon = 0.5 * (1 / (episode + 1)) if epsilon <= numpy.random.uniform(0, 1): ## まずQ値が最大の要素を確認(必ずしも1つでないので) ## 最大値がTrue、それ以外はFalseとなっているリストを得る qmax_bool_list = qtable[state]==qtable[state].max() ## Qが最大値の要素のインデックスを得る(複数あるときはその中からランダムに選ぶことになる) choices = numpy.where(qmax_bool_list)[0] else: ## ランダム行動。全要素のインデックスを得る choices = numpy.arange(len(qtable[state])) return numpy.random.choice(choices) # 候補の中から1つランダムに選ぶ
Q値を更新する
状態、行動および行動を取った結果の報酬から下記式に従いQ値を更新する。
は状態、
は行動、
は報酬、
は次の状態での最大のQ値。
は学習率、
は割引率と呼ばれる。
def update_qtable(qtable, state, action, next_state, reward): 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あたり5回行動することとして、10episode学習させる。 各episodeの初期のディスペンサー電源ON/OFFはランダムに決める。
MAX_STEP_NUMBER = 5 # 1episodeあたり5回行動することとする MAX_EPISODE_NUMBER = 10 # episodeは10までとする for episode in range(MAX_EPISODE_NUMBER): total_reward = 0 state = numpy.random.choice([0, 1]) # ディスペンサー電源初期状態はランダム env = Dispenser(state) for step in range(MAX_STEP_NUMBER): action = get_action(qtable, state, episode) # q値、状態から行動を決める next_state, reward = env.step(action) # 行動を取った結果を得る print(state, action, reward) total_reward += reward qtable = update_qtable(qtable, state, action, next_state, reward) # Q値の更新 state = next_state print('episode : %d total reward %d' %(episode+1, total_reward)) print(qtable)
学習結果
初期が電源OFFならtotal rewardは最大4、ONなら最大5となるが、episode 9, 10ではいずれもrewardは最大となっている。 初期ON/OFFがランダムでも正しく行動が取られている。
>>> 0 1 0 0 0 0 1 1 1 1 1 1 1 1 1 episode : 1 total reward 3 [[ 0. 0. ] [ 0. 1.42625]] 0 0 0 1 0 0 0 0 0 1 1 1 1 0 0 episode : 2 total reward 1 [[ 0.96271875 0. ] [ 0.57763125 1.8549375 ]] . . . 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 episode : 9 total reward 4 [[ 6.00931629 1.23142328] [ 1.25417229 7.85361236]] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 episode : 10 total reward 5 [[ 6.00931629 1.23142328] # 電源OFF時のQ値、電源ボタンを押すQ値が大きい [ 1.25417229 8.33916616]] # 電源ON時のQ値、餌ボタンを押すQ値が大きい
全コード
import numpy class Dispenser(object): def __init__(self, init_state): """ 初期のON/OFF状態を設定する init_state: 0->電源OFF、1->電源ON """ self._state = init_state self._button = self.powerbutton(init_state) def powerbutton(self, init_state): """ ON/OFFを切り替えるgenerator 0: 電源OFF、1: 電源ON を返す。内部的には-1, 1で切り替えている """ state = -1 if init_state==0 else init_state while True: state = -1 * state yield numpy.maximum(0, state) def push_powerbutton(self): self._state = self._button.__next__() def step(self, action): """ action: 0->電源ボタンを押す 1->餌ボタンを押す """ if action == 0: self.push_powerbutton() reward = 0 else: reward = 1 if self._state==1 else 0 return self._state, reward def get_action(qtable, state, episode): ## ε- greedy法 行動を重ねる(episodeが増える)毎にランダム行動を取らなくする epsilon = 0.5 * (1 / (episode + 1)) if epsilon <= numpy.random.uniform(0, 1): ## Q値に従った行動。Q値が最大の行動を取る ## まずQ値が最大の要素を確認(必ずしも1つでないので) qmax_bool_list = qtable[state]==qtable[state].max() # 複数あるときはその中からランダムに選ぶ choices = numpy.where(qmax_bool_list)[0] else: ## ランダム行動 choices = numpy.arange(len(qtable[state])) return numpy.random.choice(choices) def update_qtable(qtable, state, action, next_state, reward): 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 qtable = numpy.zeros((2, 2)) MAX_STEP_NUMBER = 5 # 1episodeあたり5回行動することとする MAX_EPISODE_NUMBER = 10 # episodeは10までとする for episode in range(MAX_EPISODE_NUMBER): total_reward = 0 state = numpy.random.choice([0, 1]) # ディスペンサー電源初期状態はランダム env = Dispenser(state) for step in range(MAX_STEP_NUMBER): action = get_action(qtable, state, episode) # q値、状態から行動を決める next_state, reward = env.step(action) # 行動を取った結果を得る print(state, action, reward) total_reward += reward qtable = update_qtable(qtable, state, action, next_state, reward) # Q値の更新 state = next_state print('episode : %d total reward %d' %(episode+1, total_reward)) print(qtable)
ニューラルネットワークによる手書き数字の認識
はじめに
scikit-learnライブラリに含まれる手書き数字データを用いて、ニューラルネットワークによる 手書き数字認識をやってみる。
データの準備
手書き数字データはサイズ 8 x 8で、グレースケールの階調が17段階となっている。具体的な手書き数字画像は scikit-learnのexampleで確認できる。 入力の個数は 8 x 8 = 64個で、出力は0〜9の10個となっている。データをPyTorchで扱うために、入力はFloatTensor、出力はLongTensorに 変換している。
## 必要なモジュールのインポート import numpy import torch import torch.nn as nn import torch.nn.functional as F from matplotlib import pyplot as plt from sklearn.datasets import load_digits # 手書き数字データ読み込み from sklearn.model_selection import train_test_split # 手書き数字データを学習用とテスト用に分類する ## データの準備 digits = load_digits() x_train, x_test, y_train, y_test = train_test_split(digits.data, digits.target, test_size=0.2) # 2割がテストデータ x_train = torch.FloatTensor(x_train) y_train = torch.LongTensor(y_train) x_test = torch.FloatTensor(x_test) y_test = torch.LongTensor(y_test)
ニューラルネットワークモデルの定義
入力64個、出力10個。中間層は2層とし、ノード数は100個とした。 中間層の活性化関数はReLUとし、出力層はl3の計算値をそのまま出力とした。
class Model(nn.Module): def __init__(self): super(Model, self).__init__() self.l1 = nn.Linear(64, 100) self.l2 = nn.Linear(100, 100) self.l3 = nn.Linear(100, 10) def forward(self, x): x = F.relu(self.l1(x)) x = F.relu(self.l2(x)) x = self.l3(x) return x
損失関数、最適化関数
損失関数にnn.CrossEntropyLossを使用した。learning_rateは試行錯誤して
良さそうなところに決めた。
model = Model()
criterion = nn.CrossEntropyLoss()
learning_rate = 1e-2
optimizer = torch.optim.Adam(model.parameters(), lr=learning_rate)
学習
epoch数は試行錯誤して良さそうなところに設定した。
epoch = 200 loss_log = [] # 学習状況のプロット用 for t in range(epoch): optimizer.zero_grad() y_model = model(x_train) loss = criterion(y_model, y_train) loss.backward() optimizer.step() loss_log.append(loss.item())
学習結果の確認
## 学習状況のプロット plt.plot(loss_log) plt.xlabel('epoch') plt.ylabel('loss') plt.yscale('log') plt.show() ## 正答率の確認 outputs = model(x_test) _, predicted = torch.max(outputs.data, 1) correct = (y_test == predicted).sum().item() # Tensorの比較で正答数を確認 print("正答率: {} %".format(correct/len(predicted)*100.0))
結果
epoch数5000までやってみたが、学習状況はなかなか落ち着いてこなかった。いっぽうで、正答率はepoch数が100でも5000でも ほとんど変わらなかった。
学習状況
正答率
>>> 正答率: 98.61111111111111 % # epoch数 100 >>> 正答率: 98.05555555555556 % # epoch数 5000
所感
人でも実際に見間違うことはあるだろうから、そんなに悪くない正答率なのではないかと思う。 畳み込みニューラルネットワークを利用したら精度が上がる？画像が 8 x 8 と小さいのであまり効果はなさそうに思われる。
参考
ニューラルネットワークによるばらつきのあるデータへの近似
はじめに
前回は関数値そのものを学習させたが、今回は関数値をもとに作成したばらつきのあるデータを学習させて 関数をニューラルネットワークで表現できるかを見てみる。
以下の記事で同様の内容がtensorflowを用いて紹介されているが、 ここではPyTorchを使う。
ニューラルネットワークは任意の関数を表現できるのか？ - Qiita
データの準備
ここでの関数はnumpy.sin(x)としている。xの値はランダムに選んだ。yの値は関数値に、正規分布に従う
乱数を加えてデータをばらつかせた。
## 必要なモジュールのインポート import numpy import torch import torch.nn as nn import torch.nn.functional as F from matplotlib import pyplot as plt ## データ作成 data_num = 200 x = numpy.random.rand(data_num)*6 # 0 - 6 の範囲の一様乱数 y = numpy.sin(x) + numpy.random.randn(data_num)*0.3 # 正規分布に従う乱数を足す x_train = torch.FloatTensor(x.reshape(data_num, 1)) y_train = torch.FloatTensor(y.reshape(data_num, 1))
ニューラルネットワークモデルの定義
入力1つ、出力1つ。中間層は1層とし、ノード数は2個とした。 中間層の活性化関数はReLUとし、出力層はl2の計算値をそのまま出力させた。。
class Model(nn.Module): def __init__(self): super(Model, self).__init__() self.l1 = nn.Linear(1, 2) self.l2 = nn.Linear(2, 1) def forward(self, x): x = F.relu(self.l1(x)) x = self.l2(x) return x
ニューラルネットワークモデルの登録
損失関数に平均二乗誤差nn.MSELossを使用した。learning_rateは試行錯誤して
良さそうなところに決めた。
model = Model()
criterion = nn.MSELoss()
learning_rate = 1e-1
optimizer = torch.optim.Adam(model.parameters(), lr=learning_rate)
学習
epoch数は試行錯誤して設定したが、計算の初期値次第で学習状況は多少変わる模様。
epoch = 100 for t in range(epoch): optimizer.zero_grad() y_model = model(x_train) loss = criterion(y_model, y_train) loss.backward() optimizer.step() print(t, loss.item()) # 学習状況の表示
学習結果の確認
plt.plot(x, y, '.', label='data') x = numpy.linspace(0, 6, data_num) x_model = torch.FloatTensor(x.reshape(data_num, 1)) y_model = model(x_model).detach().numpy().reshape(1, data_num)[0] plt.plot(x, y_model, label='model', linewidth=3) plt.legend() plt.xlabel('x') plt.ylabel('y') plt.show()
結果
計算させるたびに多少結果が変わるが、下図のようにだいぶ無理矢理な感じとなった。epoch数を増やしてもほとんど状況に変化なかった。
なお、活性化関数をsigmoidに変えると無理矢理感は低下した。
さらに中間層を2層（ノード数 20 ,40)、epoch数=1000で計算するとフィット感は改善した。
epoch数を5000回まで増やすと合わせすぎた感じとなり、いわゆる過学習状態になっていると思われる。。
所感
epoch数、中間層ノード数、層の深さなど、適切に選ぶにはそこそこの知識・経験が要りそう。
参考
メモ: PyTorch TensorDataset、DataLoader について
はじめに
PyTorchのtorch.utils.data.TensorDataset、torch.utils.data.DataLoaderの使い方についてのメモを記す。
torch.utils.data.TensorDataset
import numpy import torch import torch.utils.data x = numpy.array([[1], [2], [3], [4], [5], [6]]) y = numpy.array([[10], [20], [30], [40], [50], [60]]) x = torch.tensor(x) y = torch.tensor(y) dataset = torch.utils.data.TensorDataset(x, y)
中身をそのまま表示させると、TensorDatasetというオブジェクトであることが示される。
>>> dataset <torch.utils.data.dataset.TensorDataset object at 0x7f4e91fdca20>
[ ]で要素を取り出せる。Tensorの組がタプルになっている。
>>> dataset[0] (tensor([1.]), tensor([10.])) >>> for x, y in dataset: ... print(x, y) ... tensor([1.]) tensor([10.]) tensor([2.]) tensor([20.]) tensor([3.]) tensor([30.]) tensor([4.]) tensor([40.]) tensor([5.]) tensor([50.]) tensor([6.]) tensor([60.])
torch.utils.data.DataLoader
torch.utils.data.DataLoaderにTensorDatasetを渡すと得られ、ミニバッチを返すiterableなオブジェクトとなる。
batch_sizeでミニバッチのデータの数を指定できる(defaultは1)。shuffleをTrueとするとTensorDatasetの中身の
組の順序がシャッフルされる(defaultはFalse)。
batch_size = 3 # ミニバッチのデータの数 data_loader = torch.utils.data.DataLoader(dataset, batch_size=batch_size, shuffle=True)
forで中身を取り出すと次のようにミニバッチサイズ数ごとにデータが取り出される(この例では3つ)。順序が
シャッフルされていることが確認できる。
>>> for x, y in data_loader: ... print(x, y) ... tensor([[5.], [4.], [1.]]) tensor([[50.], [40.], [10.]]) tensor([[3.], [6.], [2.]]) tensor([[30.], [60.], [20.]])
ミニバッチ数を4、shuffle=Falseとすると次の通りとなる。データは全部で6組なので、2回めに取り出す際は残りの2組だけとなっている。
データの順序はシャッフルされずにもとのままとなっている。
>>> data_loader = torch.utils.data.DataLoader(dataset, ... batch_size=4, shuffle=False) >>> for x, y in data_loader: ... print(x, y) ... tensor([[1.], [2.], [3.], [4.]]) tensor([[10.], [20.], [30.], [40.]]) tensor([[5.], [6.]]) tensor([[50.], [60.]])
DataLoaderを用いた例
ニューラルネットワークによる関数近似 - 化学系エンジニアがAIを学ぶ
