剖析强化学习 - 第一部分

作者：Massimiliano Patacchiola

前言 [本文是对强化学习的介绍，适合已经有一些机器学习背景，并且懂一些数学和Python的读者。当我研究一种新算法时，我总是希望了解底层机制，从这个意义上讲，使用一种编程语言从头开始实现算法对理解算法是有帮助的。我在这篇文章中采用了这种方法，虽然需要花更长时间阅读但值得这样。我不是以英语为母语的人，所以如果你发现一些难以理解的错误句子，请在下面的评论中说明，我会解决它]

当我开始学习强化学习时，我没有找到任何从基础上解释什么是强化学习的好的在线资源。大多数（非常好）博客都关注现代方法（DeepReinforcement Learning），都会介绍Bellman方程但没有令人满意的解释。我把注意力转向这本书中学习，Russel和Norvig的《人工智能：一种现代的方法》。

本文基于第二版的第17章，可以被认为是本章的延伸阅读。我将使用与作者相同的数学符号，通过这种方式，您可以使用这本书来覆盖本文一些缺失的部分，反之亦然。你可以在我的github版本库中找到这篇文章中使用的完整代码，并附带pdf版本的帖子。在下一节我将介绍马尔可夫链，如果你已经知道这个概念，可以跳过下一节...

Andrey Markov

Andrey Markov是一位研究随机过程的俄罗斯数学家，他对遵循一系列相关事件的系统特别感兴趣。1906年，马尔可夫对他称之为链（Chain）的离散过程的结果产生了兴趣。一个马尔可夫链具有一组状态 S={s0,s1,...,sm}以及可以从一个状态连续移动到另一个状态的过程，每一次移动都是一个步骤（Step），并基于转移模型 T。您应该努力记住粗体字，因为我们将在文章的其余部分广泛使用它们。总结马尔可夫链的定义如下：

1.一组可能的状态：S={s0,s1,...,sm}

2.初始状态：s0

3.转移模型：T(s,s′)

马尔科夫链中有一个我没有提及的重要特性，马尔可夫链是基于马尔可夫属性，当前状态的未来变化条件上独立于过去的状态，就是过程中当前的状态只依赖于它在t-1时的状态。一个例子可以简化马尔可夫链的理解，假设我们有一个只有两个状态s0和s1，其中s0是初始状态。该过程在s0保持90％的时间，剩下的10％的时间它移动到s1，当过程处于状态s1时它将保持50％的时间。根据这些数据，我们可以创建一个如下的转移矩阵T ：

转移矩阵总是一个方形矩阵，并且由于我们处理的是概率分布，因此所有条目都在0和1之间，单个行总和为1。我们可以用图形表示马尔可夫链。在下面的表示中，链的每个状态都是一个节点，转移概率用边表示，最高的概率具有最厚边：

到目前为止，我们没有提及时间，但我们必须这样做，因为马尔可夫链是随时间发展的动态过程。假设我们猜测这个过程在3个步骤之后和50个步骤之后的状态是什么，我们该怎么做？我们注意到这是个具有有限数量的状态并且时间均匀的链，这意味着转移矩阵不随时间变化。考虑到这些假设，我们可以将k步转移概率计算为转移矩阵的k次幂，让我们在Numpy中计算：

 

import numpy as np

#Declaring the Transition Matrix T
T = np.array([[0.90, 0.10],
              [0.50, 0.50]])

#Obtaining T after 3 steps
T_3 = np.linalg.matrix_power(T, 3)
#Obtaining T after 50 steps
T_50 = np.linalg.matrix_power(T, 50)
#Obtaining T after 100 steps
T_100 = np.linalg.matrix_power(T, 100)

#Printing the matrices
print("T: " + str(T))
print("T_3: " + str(T_3))
print("T_50: " + str(T_50))
print("T_100: " + str(T_100))

T: [[ 0.9  0.1]
    [ 0.5  0.5]]

T_3: [[ 0.844  0.156]
      [ 0.78   0.22 ]]

T_50: [[ 0.83333333  0.16666667]
       [ 0.83333333  0.16666667]]

T_100: [[ 0.83333333  0.16666667]
        [ 0.83333333  0.16666667]]

现在我们定义代表k = 0时系统状态的初始分布，系统由两个状态组成，我们可以将初始分布建模为具有两个元素的向量，向量的第一个元素表示停留在状态s0的概率，第二个元素是停留在状态s1的概率。我们假设从s0开始，向量v 代表下列形式的初始分布：

v=（1，0）

我们可以计算出在k次迭代之后的特定状态的概率：初始分布和转移矩阵相乘：v⋅Tk。让我们在Numpy中实现这一点：

 

import numpy as np

#Declaring the initial distribution
v = np.array([[1.0, 0.0]])
#Declaring the Transition Matrix T
T = np.array([[0.90, 0.10],
              [0.50, 0.50]])

#Obtaining T after 3 steps
T_3 = np.linalg.matrix_power(T, 3)
#Obtaining T after 50 steps
T_50 = np.linalg.matrix_power(T, 50)
#Obtaining T after 100 steps
T_100 = np.linalg.matrix_power(T, 100)

#Printing the initial distribution
print("v: " + str(v))
print("v_1: " + str(np.dot(v,T)))
print("v_3: " + str(np.dot(v,T_3)))
print("v_50: " + str(np.dot(v,T_50)))
print("v_100: " + str(np.dot(v,T_100)))

v: [[ 1.  0.]]

v_1: [[ 0.9  0.1]]

v_3: [[ 0.844  0.156]]

v_50: [[ 0.83333333  0.16666667]]

v_100: [[ 0.83333333  0.16666667]]

发生了什么？该过程从s0开始经过一次迭代后，我们有90%确定它仍然处于这个状态。这很容易理解，我们的转移模型说明，该过程有90％的概率可以留在s0状态，没有什么新意。展望k = 3时的状态分布，我们注意到有一些不同，未来有更多的可能性进入不同的状态。如果我们想找到经过三次迭代后状态处于s0的概率，我们应该把所有可能导致s0的分支概率加起来。一张图片胜过千言万语：

在k = 3时处于s0的概率由（0.729 + 0.045 +0.045 + 0.025）给出，它等于0.844，我们得到了相同的结果。现在假设在开始时我们对于过程的开始状态有一些不确定性，让我们定义另一个起始向量，如下所示：

v =（0.5，0.5）

就是说，有50％的概率，我们可以从s0开始。再次运行Python脚本，我们在1次，3次，50次和100次迭代后打印结果：

 

v: [[ 0.5, 0.5]]

v_1: [[ 0.7  0.3]]

v_3: [[ 0.812  0.188]]

v_50: [[ 0.83333333  0.16666667]]

v_100: [[ 0.83333333  0.16666667]]

这一次在k = 3时进入s0的概率较低（0.812），但在长时间执行后，我们有相同的结果（0.8333333）。长时间运行发生了什么？经过50和100次迭代的结果是相同的，不管开始是什么概率始终v_50=v_100。该链收敛于平衡状态，意味着随着时间的推移，它忘记了开始的分配。但我们必须小心，收敛并不总能得到保证，马尔可夫链的动态性可能非常复杂，特别是可能存在瞬态和循环状态。对于我们的讨论范围，这些就已经足够了。我建议你看看setosa.io博客，因为他们有一个用于马尔可夫链可视化的交互式页面。

马尔科夫决策过程（MDP）

在强化学习中，经常使用与马尔可夫链有关的概念，下面我们说说马尔可夫决策过程（MDP）。MDP是马尔可夫链的重新释义，包括Agent和决策过程。MDP由这些组件定义：

1.  可能的状态集合: S={s0,s1,...,sm}

2.  初始状态: s0

3.  可能的动作（action）集合: A={a0,a1,...,an}

4.  转移模型: T(s,a,s′)

5.  奖励（Reward）函数: R(s)

正如您所看到的，我们正在引入一些关于马尔可夫链的新元素，特别是依赖于当前状态、下一个状态和Agent动作的转移模型。转移模型返回在状态s完成动作a到达状态s′的概率。但是给定s和a的该模型有条件地独立于以前的所有状态和行为（马尔可夫属性）。此外还有奖励函数 R(s)，每当Agent从一个状态移动到另一个状态时，它都会返回一个实数值（注意：将奖励函数定义为仅依赖于s可能会引起混淆，Russel和Norvig在书中使用这种表示法来简化描述，但并没有改变问题的本质）。由于有奖励功能，我们可以说一些状态更期望某些状态，因为当Agent在这些状态移动时，它会得到更高的奖励。相反，有些状态恰恰相反，因为当Agent移动到这些状态时，会收到负的奖励。

·        问题：Agent必须最大化获得的奖励，避免哪些返回负值的状态，而是选择返回正值的状态。

·        解决方案：找到一个策略 π(s) 它返回最高奖励的动作。

Agent可以尝试不同的策略，但只有其中一个可以被认为是最优策略，用π*表示，它会产生最高的预期效用(utility)。现在开始介绍我将在所有帖子中使用的例子，这个例子的灵感来自Russel和Norving在他们书中第17.1章中提出的简单环境，假设我们有一台必须到达充电站的清洁机器人，我们的简单世界是一个4×3矩阵，其中起点s0在（1,1）位置，（4,3）位置是充电站，（4,2）位置是危险楼梯，以及（2,2）位置是障碍物。机器人必须找到到达充电站的最佳方式（奖励+1）并避免跌落楼梯（奖励-1）。每当机器人做出决定时，就有可能受到随机因素（例如地面滑溜，一只恶猫刺痛机器人）的干扰，这会使机器人在20％的时间内偏离原始路径。机器人向前走时将在10％的情况下向左，10％的情况向右的状态，如果机器人撞到墙壁或障碍物，它会弹回到原来的位置。这个世界的主要特点如下：

·        离散的时间和空间

·        完全可观察

·        无限的视野

·        已知的转移模型

环境是完全可观察的，这意味着机器人总是知道它处于哪个状态。应该进一步解释术语无限时域，无限时域意味着，没有一个固定的时间限制。如果Agent有相同的两种终止状态的策略，它将永远持续下去，这个假设并不意味着每个剧集(episode )都需要Agent走过一系列无限状态，当两个终止状态中的一个到达时，该剧集结束。下面说明了这个世界和转移模型的表示，小心Russell和Norvig使用的状态索引，这可能会让人困惑，他们从左下角开始按列和行来命名世界的每个状态。

机器人的目标是找到到达充电站的最佳方式，这意味着什么才是最好的方式呢？根据机器人在每个中间状态下获得的奖励类型，我们可以有不同的最优策略π*。假设我们正在设计机器人的固件，根据电池级别，我们在每个时间步骤给予不同的奖励，两个终止状态的奖励保持不变（充电器= + 1，楼梯= -1），（2,2）处的障碍不是有效的状态，因此没有与其相关的奖励。鉴于这些假设，我们可以有四种不同的情况：

1.R(s)≤−1.6284R(s)≤−1.6284 极低的电池

2.−0.4278≤R(s)≤−0.085−0.4278≤R(s)≤−0.085 较低的电池

3.−0.0221≤R(s)≤0−0.0221≤R(s)≤0 稍低的电池

4.R(s)>0R(s)>0 完全充满

对于这些条件中的每一个，我们都可以尝试猜测Agent将选择哪种策略。在极低的电池情况下，Agent会受到如此高的惩罚，以至于只想尽快让伤痛停止，生活是如此痛苦，以至于坠下楼梯是一个不错的选择。在非常低的电池情况下，Agent将选择到充电站的最短路径，它不关心是否跌倒。在电池电量稍低的情况下，机器人完全不承担风险，并躲避楼梯避免撞墙的代价。最后，在充满电的情况下，Agent避免了两个终止状态并保持在稳定状态，在每个时间步骤都获得积极回报。

到目前为止，我们知道哪些策略可以在具有确定奖励的特定环境中出现，但我仍然没有谈到：Agent如何选择最佳策略？

Bellman方程

上一节结束时提出了一个问题：Agent如何选择最佳策略？为了回答这个问题，我将介绍Bellman方程。首先，必须找到一种方法来比较两种策略，我们可以使用每个状态给出的奖励来获得状态序列效用的度量。我们定义状态历史 h的效用如：

Uh=R(s0)+γR(s1)+γ2R(s2)+...+γnR(sn)

前面的公式定义了状态序列的折扣奖励，其中的 γ∈[0,1]被称为折扣因子。折扣因子描述了Agent关于当前奖励对于未来奖励的偏好，折扣因子为1.0会导致先前公式退化为累加奖励。折扣奖励不是我们估算效用的唯一方式，但它是减少问题的方法。例如，在无限状态序列的情况下，折扣奖励会产生确定的效用（使用无限序列的奖励汇总），而且我们还可以使用每个时间步骤获得的平均奖励来比较无限序列。如何比较单一状态的效用？效用U(s)可以定义为：

让我们记住效用是根据策略π定义的，为了简洁起见我上面没有提到。一旦我们有了效用，我们如何为下一个状态选择最佳的动作？利用最大期望效用原则，即理性的Agent应该选择一种最大化其预期效用的动作。我们离Bellman公式更近了一步，一个状态s的效用和与它相邻的状态 s′的效用相关，即：

我们刚刚推导了Bellman方程！使用Bellman方程，Agent可以估计最佳的动作并找出最佳策略。我们试着剖析这个等式,首先，术语R(s)是我们必须在等式中加上的东西，我们知道效用必须考虑在状态s给予该状态的奖励；其次，请注意，该公式使用转移模型T乘以下一个状态s的效用。如果你认为这是有道理的，那么发生概率较低的状态（如我们简化世界中左侧和右侧移动的概率为10％）在总和中将具有最低的权重。

我们将在简化的4x3环境中使用我们的清洁机器人，以经验方式测试Bellman方程。在这个例子中，每个非终止状态的奖励是 R(s)=−0.04。我们可以想当然每个状态的效用值，因为目前你不需要知道我们如何得到这些值，可以想象它们是神奇的。以同样神奇的方式，我们为世界获得了最佳策略（重新检查我们从Bellman方程中获得的结果是否合理）。这张图片非常重要，牢记在心。

在我们的例子中，我们假设机器人从状态（1,1）开始，使用Bellman方程，我们必须找到 UP、LEFT、DOWN和RIGHT中效用最高的动作。我们没有最优策略，但我们有每个状态的转移模型和效用值。你必须记住我们环境的两个主要规则，首先，如果机器人碰到墙，它会回到之前的状态；其次，根据转移模型，选定的动作仅以80％的概率执行。我不想说明简陋的数字，而是想向您展示可能结果的图形说明：

对于每个可能的结果，我都说明了效用和转移模型给出的概率，这对应于Bellman方程的第一部分。下一步是计算效用和转移概率之间的乘积，然后总结每个动作的值。

我们发现，对于状态（1,1），UP的行为具有最高的价值。这符合我们神奇般获得的最佳策略。Bellman方程的这一部分返回最大化后续状态的期望效用的动作，这是最优策略应该做的事情：

现在我们有了所有的元素，并且我们可以将值带入到Bellman方程得到状态（1,1）的效用：

U(s11)=−0.04+1.0×0.7456=0.7056

这就是Bellman方程的工作原理！下面说明在我们的模拟世界中使用的等式的Python实现，我们将使用前面部分的相同术语。我们的世界有4x3 = 12个可能的状态，起始矢量包含12个值，转移矩阵是一个巨大的12x12x4矩阵（12个起始状态，12个下一个状态，4个动作），其中大部分值为零（我们只能从一个状态移动到其邻近状态）。我使用脚本生成转移矩阵，并将其保存为Numpy矩阵（可以在此处下载）。在脚本中我定义了这个函数return_state_utility()，它是Bellman方程的一个实现。使用这个函数，我们将打印状态（1,1）的效用并检查它是否与我们之前找到的一样：

 

import numpy as np

def return_state_utility(v, T, u, reward, gamma):
    """Return the state utility.

    @param v the state vector
    @param T transition matrix
    @param u utility vector
    @param reward for that state
    @param gamma discount factor
    @return the utility of the state
    """
    action_array = np.zeros(4)
    for action in range(0, 4):
        action_array[action] = np.sum(np.multiply(u, np.dot(v, T[:,:,action])))
    return reward + gamma * np.max(action_array)

def main():
    #Starting state vector
    #The agent starts from (1, 1)
    v = np.array([[0.0, 0.0, 0.0, 0.0, 
                   0.0, 0.0, 0.0, 0.0, 
                   1.0, 0.0, 0.0, 0.0]])

    #Transition matrix loaded from file
    #(It is too big to write here)
    T = np.load("T.npy")

    #Utility vector
    u = np.array([[0.812, 0.868, 0.918,   1.0,
                   0.762,   0.0, 0.660,  -1.0,
                   0.705, 0.655, 0.611, 0.388]])

    #Defining the reward for state (1,1)
    reward = -0.04
    #Assuming that the discount factor is equal to 1.0
    gamma = 1.0

    #Use the Bellman equation to find the utility of state (1,1)
    utility_11 = return_state_utility(v, T, u, reward, gamma)
    print("Utility of state (1,1): " + str(utility_11))

if __name__ == "__main__":
    main()

Utility of state (1,1): 0.7056

很好，我们获得了完全相同的价值！到现在为止，我们认为效用价值神奇般地出现。我们不想使用魔法，而是想找到一种算法来获得这些值。有一个问题，对于n个可能的状态有n个Bellman方程，每个方程包含n个未知数。使用任何线性代数包将有可能解决这些方程，问题是因为max的操作导致它们不是线性的。该怎么办？我们可以使用值迭代算法...

值迭代算法

Bellman方程是用于求解MDP的值迭代算法的核心。我们的目标是找到每个状态的效用（也称为价值）。正如前面所说的，不能使用线性代数库，我们需要一种迭代方法，从任意初始效用值（通常为零）开始，然后使用Bellman方程计算一个状态的效用，并将其分配给状态，这个迭代过程被称为Bellman更新。无限地应用Bellman更新可以确保达到平衡，一旦我们达到了平衡状态，我们就可以找到我们正在寻找的效用值，并且我们可以用它们来估计每个状态的最佳动作。如何理解算法何时达到平衡？我们需要一个停止标准，考虑到两次连续迭代之间的效用，我们可以在没有状态变化的情况下停止算法。

这个结果是收缩属性的结果，我会跳过它，因为它在书中的第17.2章已经很好地解释了。好的，是时候用Python实现算法了。我将重新使用return_state_utility()函数来更新效用向量u。

 

import numpy as np

def return_state_utility(v, T, u, reward, gamma):
    """Return the state utility.

    @param v the value vector
    @param T transition matrix
    @param u utility vector
    @param reward for that state
    @param gamma discount factor
    @return the utility of the state
    """
    action_array = np.zeros(4)
    for action in range(0, 4):
        action_array[action] = np.sum(np.multiply(u, np.dot(v, T[:,:,action])))
    return reward + gamma * np.max(action_array)

def main():
    #Change as you want
    tot_states = 12
    gamma = 0.999 #Discount factor
    iteration = 0 #Iteration counter
    epsilon = 0.01 #Stopping criteria small value

    #List containing the data for each iteation
    graph_list = list()

    #Transition matrix loaded from file (It is too big to write here)
    T = np.load("T.npy")

    #Reward vector
    r = np.array([-0.04, -0.04, -0.04,  +1.0,
                  -0.04,   0.0, -0.04,  -1.0,
                  -0.04, -0.04, -0.04, -0.04])    

    #Utility vectors
    u = np.array([0.0, 0.0, 0.0,  0.0,
                   0.0, 0.0, 0.0,  0.0,
                   0.0, 0.0, 0.0,  0.0])
    u1 = np.array([0.0, 0.0, 0.0,  0.0,
                    0.0, 0.0, 0.0,  0.0,
                    0.0, 0.0, 0.0,  0.0])

    while True:
        delta = 0
        u = u1.copy()
        iteration += 1
        graph_list.append(u)
        for s in range(tot_states):
            reward = r[s]
            v = np.zeros((1,tot_states))
            v[0,s] = 1.0
            u1[s] = return_state_utility(v, T, u, reward, gamma)
            delta = max(delta, np.abs(u1[s] - u[s])) #Stopping criteria       
        if delta < epsilon * (1 - gamma) / gamma:
                print("=================== FINAL RESULT ==================")
                print("Iterations: " + str(iteration))
                print("Delta: " + str(delta))
                print("Gamma: " + str(gamma))
                print("Epsilon: " + str(epsilon))
                print("===================================================")
                print(u[0:4])
                print(u[4:8])
                print(u[8:12])
                print("===================================================")
                break

if __name__ == "__main__":
    main()

在收敛期间看看每个效用的稳定是很有趣的，我使用matplotlib绘制每个状态25次迭代的效用值。

使用相同的代码，针对折扣因子gamma使用不同值执行不同的模拟，当折扣因子接近1.0时，我们对效用的预测变得更加精确，在 gamma = 1.0 的极限情况下，算法将永远不会结束，因为我们永远不会达到停止规则。

 

=================== FINAL RESULT ==================
Iterations: 9
Delta: 0.000304045
Gamma: 0.5
Epsilon: 0.001
===================================================
[[ 0.00854086  0.12551955  0.38243452  1.        ]]
[[-0.04081336  0.          0.06628399 -1.        ]]
[[-0.06241921 -0.05337728 -0.01991461 -0.07463402]]
===================================================

=================== FINAL RESULT ==================
Iterations: 16
Delta: 0.000104779638547
Gamma: 0.9
Epsilon: 0.001
===================================================
[[ 0.50939438  0.64958568  0.79536209  1.        ]]
[[ 0.39844322  0.          0.48644002 -1.        ]]
[[ 0.29628832  0.253867    0.34475423  0.12987275]]
===================================================

=================== FINAL RESULT ==================
Iterations: 29
Delta: 9.97973302774e-07
Gamma: 0.999
Epsilon: 0.001
===================================================
[[ 0.80796344  0.86539911  0.91653199  1.        ]]
[[ 0.75696623  0.          0.65836281 -1.        ]]
[[ 0.69968285  0.64882069  0.6047189   0.38150244]]
===================================================

还有另一种算法可以让我们找到效用向量，同时还有一个最优策略，那就是策略迭代算法。

策略迭代算法

利用值迭代算法，我们可以估计每个状态的效用。我们仍然想找到估算最佳策略的另一种方法，在本节中，我将向您展示如何使用策略迭代算法来寻找最大化期望奖励的最优策略，没有其它策略会产生比最优策略π*更多的奖励。策略迭代保证收敛，当前策略及其效用函数是最优策略和最优效用函数。首先，我们定义一个策略π它为每个状态分配一个动作，我们可以为此策略分配随机动作，这并不重要。使用return_state_utility()函数（Bellman方程），我们可以计算策略的预期效用。有一个好消息,我们并不需要Bellman方程的完整版本，即：

由于我们有一个策略和策略中每个状态相关联的动作，我们可以去掉max运算符并使用Bellman等式的简化版本：

当评估策略时，我们可以改进它，该策略的改进是算法的第二个和最后一个步骤。我们的环境具有有限数量的状态及有限数量的策略，每次迭代都会产生更好的策略。我实现了一个函数return_policy_evaluation()，它包含Bellman方程的简化版本。此外，我们需要函数return_expected_action()返回基于当前值u和T具有最高效用的动作。为了检查发生了什么，我还创建了一个print函数，它将策略向量p中包含的每个动作映射到一个符号并将其打印在终端上。

 

import numpy as np

def return_policy_evaluation(p, u, r, T, gamma):
  """Return the policy utility.

  @param p policy vector
  @param u utility vector
  @param r reward vector
  @param T transition matrix
  @param gamma discount factor
  @return the utility vector u
  """
  for s in range(12):
    if not np.isnan(p[s]):
      v = np.zeros((1,12))
      v[0,s] = 1.0
      action = int(p[s])
      u[s] = r[s] + gamma * np.sum(np.multiply(u, np.dot(v, T[:,:,action])))
  return u

def return_expected_action(u, T, v):
    """Return the expected action.

    It returns an action based on the
    expected utility of doing a in state s, 
    according to T and u. This action is
    the one that maximize the expected
    utility.
    @param u utility vector
    @param T transition matrix
    @param v starting vector
    @return expected action (int)
    """
    actions_array = np.zeros(4)
    for action in range(4):
       #Expected utility of doing a in state s, according to T and u.
       actions_array[action] = np.sum(np.multiply(u, np.dot(v, T[:,:,action])))
    return np.argmax(actions_array)

def print_policy(p, shape):
    """Printing utility.

    Print the policy actions using symbols:
    ^, v, <, > up, down, left, right
    * terminal states
    # obstacles
    """
    counter = 0
    policy_string = ""
    for row in range(shape[0]):
        for col in range(shape[1]):
            if(p[counter] == -1): policy_string += " *  "            
            elif(p[counter] == 0): policy_string += " ^  "
            elif(p[counter] == 1): policy_string += " <  "
            elif(p[counter] == 2): policy_string += " v  "           
            elif(p[counter] == 3): policy_string += " >  "
            elif(np.isnan(p[counter])): policy_string += " #  "
            counter += 1
        policy_string += '\n'
    print(policy_string)

现在我将在main函数的循环中使用这些函数，这是一个策略迭代算法的实现。我声明了一个新的向量p，其中包含每个状态的动作，该算法的停止条件是两次连续迭代的效用向量之间的差，当步骤中对效用的改进不会发生变化（或发生非常小的变化）时，算法终止。

def main():

 

    gamma = 0.999
    epsilon = 0.0001
    iteration = 0
    T = np.load("T.npy")
    #Generate the first policy randomly
    # NaN=Nothing, -1=Terminal, 0=Up, 1=Left, 2=Down, 3=Right
    p = np.random.randint(0, 4, size=(12)).astype(np.float32)
    p[5] = np.NaN
    p[3] = p[7] = -1
    #Utility vectors
    u = np.array([0.0, 0.0, 0.0,  0.0,
                  0.0, 0.0, 0.0,  0.0,
                  0.0, 0.0, 0.0,  0.0])
    #Reward vector
    r = np.array([-0.04, -0.04, -0.04,  +1.0,
                  -0.04,   0.0, -0.04,  -1.0,
                  -0.04, -0.04, -0.04, -0.04])

    while True:
        iteration += 1
        #1- Policy evaluation
        u_0 = u.copy()
        u = return_policy_evaluation(p, u, r, T, gamma)
        #Stopping criteria
        delta = np.absolute(u - u_0).max()
        if delta < epsilon * (1 - gamma) / gamma: break
        for s in range(12):
            if not np.isnan(p[s]) and not p[s]==-1:
                v = np.zeros((1,12))
                v[0,s] = 1.0
                #2- Policy improvement
                a = return_expected_action(u, T, v)         
                if a != p[s]: p[s] = a
        print_policy(p, shape=(3,4))

    print("=================== FINAL RESULT ==================")
    print("Iterations: " + str(iteration))
    print("Delta: " + str(delta))
    print("Gamma: " + str(gamma))
    print("Epsilon: " + str(epsilon))
    print("===================================================")
    print(u[0:4])
    print(u[4:8])
    print(u[8:12])
    print("===================================================")
    print_policy(p, shape=(3,4))
    print("===================================================")

if __name__ == "__main__":
    main()

使用gamma=0.999和epsilon=0.0001运行脚本，在22次迭代后收敛，我们得到下面的结果:

 

=================== FINAL RESULT ==================
Iterations: 22
Delta: 9.03617490833e-08
Gamma: 0.999
Epsilon: 0.0001
===================================================
[ 0.80796344  0.86539911  0.91653199  1.        ]
[ 0.75696624  0.          0.65836281 -1.        ]
[ 0.69968295  0.64882105  0.60471972  0.38150427]
===================================================
 >   >   >   *  
 ^   #   ^   *  
 ^   <   <   <  
===================================================

算法返回的最终策略等于最优策略。此外，使用简化的Bellman方程，该算法设法找到效用向量的最优值。如果我们看看策略的演变，我们会注意到一些有趣的现象，开始时，策略是随机生成的，经过四次迭代后，算法找到次最优策略并坚持到第10次迭代，直到找到最优策略。从迭代10到迭代22，算法根本不改变策略。一个次优的策略在无模型强化学习中可能是一个问题，因为贪婪的Agent可能停留在次优策略，就目前对我们来说它不是一个问题。

策略迭代和值迭代，哪个更好？如果您有很多动作，或者您从合理的策略开始，请选择策略迭代；如果您有很少的动作，并且转移是非循环的，那么选择值迭代。如果你想从两个中获得最好的结果，那么请看看修改过的策略迭代算法。

使用线性代数的策略评估

我说过，从Bellman方程中消除max操作使得我们的世界变得更加简单，因为我们可以使用任何线性代数包来计算效用。在最后一节中，我想说明如何使用线性代数方法达到相同的结果。在Bellman方程中，有一个带n个变量和n个约束的线性系统，我们需要处理矩阵和向量。给定一个策略p、与状态s相关的动作、奖励向量r、转移矩阵T和折扣因子gamma，我们可以用一行代码来估计效用：

u[s]= np.linalg.solve(np.identity(12)- gamma*T[:,:,p[s]],r)[s]

我使用了Numpy方法np.linalg.solve，它将系数矩阵A和一个数组b作为参数输入，该方法返回给系统A x = b的解。对于参数矩阵A我传递单位矩阵I和gamma * T之间的差，对于参数数组b我传递奖励矢量r。为什么用I -gamma*T作为第一个参数传递 ？我们可以从简化的Bellman方程推导出这个值：

事实上，我们可以在Numpy中实现最后一个等式得到u：

u[s]=np.dot(np.linalg.inv(np.identity(12)-gamma*T[:,:,p[s]]),r)[s]

不过，我更喜欢使用np.linalg.solve，因为它更具可读性。

结论

在第一部分中，我总结了强化学习的基本思想。例如，我使用了具有预定义转移模型的有限环境。如果我们没有转移模型，会发生什么？在接下来的部分中，我将介绍无模型强化学习，它用一组新的有趣工具回答这个问题。你可以在我的github存储库中找到完整的代码。 

 

索引

1. [第一篇]马尔科夫决策过程，贝尔曼方程，值迭代和策略迭代算法。

2.  [第二篇]蒙特卡罗概念，蒙特卡洛方法，预测与控制，广义策略迭代，Q函数。

3.  [第三篇]时间差分概念，动物学习，TD(0),TD(λ)和资格痕迹，SARSA，Q-learning。

4.  [第四篇]Actor-Critic方法背后的神经生物学，计算Actor-Critic方法，Actor-only和Critic-only方法。

5.  [第五篇]进化算法介绍，强化学习中的遗传算法，遗传算法的策略选择。

6.  [第六篇]强化学习应用，多臂老虎机，山地车，倒立摆，无人机着陆，难题。

7.  [第七篇]函数逼近概念，线性逼近器，应用，高阶逼近器。

8. [第八篇] 非线性函数逼近，感知器，多层感知器，应用，政策梯度。

资源

The dissecting-reinforcement-learning repository.

The setosa blog containing a good-lookingsimulator for Markov chains.

Official github repository for the book “ArtificialIntelligence: a Modern Approach”.

参考

Bellman,R. (1957). A Markovian decision process (No. P-1066). RAND CORP SANTA MONICACA.

Russell,S. J., Norvig, P., Canny, J. F., Malik, J. M., & Edwards, D. D. (2003).Artificial intelligence: a modern approach (Vol. 2). Upper Saddle River:Prentice hall.