{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 第十六讲：马尔可夫决策过程\n",
    "\n",
    "# 第十三部分：强化学习及控制\n",
    "\n",
    "接下来我们将要学习的是**强化学习（RL: Reinforcement Learning）**以及**自适应控制（Adaptive Control）**。\n",
    "\n",
    "回顾前面学过的内容：\n",
    "\n",
    "* 在监督学习中，算法尝试模仿训练集中的数据，给新的输入$x$贴上一个合适的标签$y$。在这种情形下，$y$会明确的指出每个输入$x$所对应的那个“正确的”标签；\n",
    "* 在无监督学习中，算法尝试发掘数据中潜在的固有结构，使用模型将这种结构表示出来，之后再对数据进行标记，使其存在于某个结构中。\n",
    "\n",
    "但是在很多需要作出渐进决策及控制的问题中，我们很难提供这种明确的监督，去告诉算法什么是正确、什么是错误。比如我们建造了一个四腿机器人，然后尝试通过编程让它行走。从一开始我们就无法定义对于行走来说，什么是正确的动作。所以也就无法为算法提供一个明确的监督方案令其模仿了。\n",
    "\n",
    "在强化学习中，我们只会给算法提供一个**奖励函数（reword function）**，这个函数会告诉算法在什么情况下是做的好，在什么情况下是做的不好，比如在四腿机器人的例子中，当机器人向前行走时奖励函数会给算法正面的反馈，而当机器人无故后退或翻倒时函数则会给算法负面的反馈。而学习算法的任务就是自主发现“通过做出怎样的动作才能获得更多的奖励”。\n",
    "\n",
    "强化学习算法在直升机自动驾驶、机器人腿部移动、手机网络路由、市场策略选择、工厂控制、高效网页索引等领域都有非常成功的应用案例。而我们对强化学习算法的介绍将从**[马尔可夫决策过程](https://zh.wikipedia.org/wiki/%E9%A6%AC%E5%8F%AF%E5%A4%AB%E6%B1%BA%E7%AD%96%E9%81%8E%E7%A8%8B)（[MDP: Markov Decision Process](https://en.wikipedia.org/wiki/Markov_decision_process)）**开始，MDP将形式化RL算法经常遇到的问题。\n",
    "\n",
    "## 1. 马尔可夫决策过程\n",
    "\n",
    "一个马尔可夫决策过程就是一个元组$(S,A,\\{P_{sa},\\gamma,R\\})$，其中：\n",
    "\n",
    "* $S$是一个**状态（states）**的集合；（在直升机自动驾驶的例子中，就代表了直升机所有可能的位置和姿态。）\n",
    "* $A$是一个**动作（action）**的集合；（在直升机自动驾驶的例子中，就代表了通过操作杆可以对直升机做出的所有动作。）\n",
    "* $P_{sa}$是状态转换概率，对每个状态$s\\in S$和动作$a\\in A$，$P_{sa}$是一个在状态空间上的分布。在后面我们会详细讨论这一点，简单的说，就是如果在状态$s$下发生了动作$a$，“那么下一步所将变为哪一个状态”的概率分布就由$P_{sa}$确定；\n",
    "* $\\gamma\\in[0,1)$称为**折扣因子（discount factor）**，用来调整未来奖励与当作奖励之间的权重；\n",
    "* $R:S\\times A\\to\\mathbb R$就是奖励函数。（有时只把奖励函数当做关于状态$S$的函数$R:S\\to\\mathbb R$。）\n",
    "\n",
    "简要的描述一下MDP的动态过程：首先，我们从状态$s_0$开始，选择在MDP中做一个动作$a_0\\in A$。动作的结果就是产生一个随机的状态转换，转换到下一个状态$s_1$，而$s_1\\sim P_{s_0a_0}$（也就是事件“转换到下一个状态”服从由当前状态$s_0,a_0$确定的分布）。接下来我们继续在MDP中做出下一个动作$a_1$。同样，$a_1$的发生使得MDP转换到下一个状态$s_2\\sim P_{s_1a_0}$。之后我们继续做出动作$a_2$，等等……我们可以用下图描述上面的文字：\n",
    "\n",
    "$$\n",
    "\\require{AMScd}\n",
    "\\begin{CD}\n",
    "s_0@>{a_0}>>s_1@>{a_1}>>s_2@>{a_2}>>s_3@>{a_3}>>\\cdots\\\\\n",
    "\\end{CD}\n",
    "$$\n",
    "\n",
    "MDP在做出动作$a_0,a_1,\\cdots$后从状态$s_0,s_1$中得到的“总收益”为：\n",
    "\n",
    "$$\n",
    "R(s_0,a_0)+\\gamma R(s_1,a_1)+\\gamma^2R(s_2,a_2)+\\cdots\n",
    "$$\n",
    "\n",
    "也可以写成只与状态相关的函数：\n",
    "\n",
    "$$\n",
    "R(s_0)+\\gamma R(s_1)+\\gamma^2R(s_2)+\\cdots\n",
    "$$\n",
    "\n",
    "对于大多数情况，我们都使用较为简单的状态奖励函数$R(s)$，尽管更加一般化的状态-动作奖励函数$R(s,a)$也并没有增加什么计算难度。\n",
    "\n",
    "在强化学习中，我们的目标是按照时间顺序依次选择动作，使得整个MDP获得的“总收益”最大：\n",
    "\n",
    "$$\n",
    "\\mathrm E\\left[R(s_0)+\\gamma R(s_1)+\\gamma^2R(s_2)+\\cdots\\right]\n",
    "$$\n",
    "\n",
    "注意到第$t$步的奖励被折扣因子$\\gamma^t$降低了（而且随着步数的推移折扣会越来越大），于是，为了保证期望值足够大，我们倾向于尽可能早的获得正反馈（并让负反馈到来的越晚越好）。在经济学应用领域，$R(\\cdot)$代表赚到的钱，$\\gamma$也有自然的解释——利率（也就是今天的美元比明天的美元更值钱）。\n",
    "\n",
    "**策略（policy）**是一个函数$\\pi:S\\to A$，它是一个从状态到动作的映射。不论何时，只要处于状态$s$，就采用动作$a=\\pi(s)$，我们称这个过程为**执行（executing）**策略$\\pi$（在复杂模型中，策略不止通过当前状态，还可能通过过去的状态及动作，对下一步动作做出判断；而在本课程中，我们仅使用当前状态作为策略的输入）。为策略$\\pi$定义**价值函数（value function）**：\n",
    "\n",
    "$$\n",
    "V^\\pi(s)=\\mathrm E\\left[\\underbrace{R(s_0)}_\\textrm{immediate reward}+\\underbrace{\\gamma R(s_1)+\\gamma^2R(s_2)+\\cdots}_\\textrm{future reward}\\mid s_0=s,\\pi\\right]\n",
    "$$\n",
    "\n",
    "$V^\\pi(s)$是从状态$s$开始依照策略$\\pi$所得到的折扣奖励的预期总收益。（式中的这种标记法理论上说是不正确的，因为$\\pi$并不是随机变量，不过这种写法在文献中是标准记法。）这个式子由两部分组成：\n",
    "* 我们通常称前半部分为**即时奖励（immediate reward）**$R(s)$，即从状态$s$开始时直接得到的奖励；\n",
    "* 称后半部分为**未来奖励（future reward）**，即当前状态之后的所有步骤决策所产生的折扣奖励的总收益。\n",
    "\n",
    "容易看出，这个式子可以写成$\\displaystyle V^\\pi(s)=\\mathrm E\\left[R(s_0)+\\gamma\\underbrace{\\left(R(s_1)+\\gamma R(s_2)+\\cdots\\right)}_{V^\\pi\\left(s'\\right)}\\mid s_0=s,\\pi\\right]$（为了方便，我们暂时定义映射$s_0\\to s,\\ s_1\\to s'$），这是一个递归定义，可以按照递归式写作：\n",
    "\n",
    "$$\n",
    "V^\\pi(s)=R(s)+\\gamma\\sum_{s'\\in S}P_{s\\pi(s)}\\left(s'\\right)V^\\pi\\left(s'\\right)\n",
    "$$\n",
    "\n",
    "我们之所以不把后半部分直接写成$\\gamma V^\\pi\\left(s'\\right)$，是因为在位于前一个状态时，下一个状态是一个随机变量。所以上式的未来奖励部分写成了期望的定义式$\\displaystyle\\operatorname*{E}_{s'\\sim P_{s\\pi(s)}}\\left[V^\\pi\\left(s'\\right)\\right]$——即“位于状态$s$时，执行策略$\\pi(s)$，到达状态$s'$的概率：$P_{s\\pi(s)}\\left(s'\\right)$（可以看做是$P_{sa}\\left(s'\\right)$，而$a=\\pi(s)$）”与“状态$s'$下的价值函数：$V^\\pi\\left(s'\\right)$”之积在所有状态$s'\\in S$上求和——这是从状态$s'$开始依照策略$\\pi$所得到的折扣奖励的预期总收益。其中$s'$是一个与概率$P_{s\\pi(s)}$相关的分布，这是一个在MDP中从状态$s$开始，按照策略$\\pi(s)$执行动作后，落在下一个状态$s'$上时对应的分布。因此，第二项就是MDP在执行第一步之后，在未来能够获得的折扣奖励的预期总收益。\n",
    "\n",
    "上面这个式子也称作**贝尔曼方程（Bellman Equation）**，这是我们在解决MDP问题时用到的主要方程之一。我们可以说，对于给定的策略$\\pi$，其价值函数$V^\\pi$满足贝尔曼方程。\n",
    "\n",
    "贝尔曼方程可以高效的解出$V^\\pi$，尤其是对有限状态的MDP（$\\lvert S\\rvert\\lt\\infty$），我们可以为每一个状态$s$计算一次$V^\\pi(s)$。如此就可以得到一个由“关于$V^\\pi(s)$的有$\\lvert S\\rvert$个变量（每个状态$s\\in S$都需要一个预期总收益$V^\\pi(s)$）同时具有$\\lvert S\\rvert$个方程（每个状态$s\\in S$的价值函数$V^\\pi(s)$都有一个方程）的线性方程”组成线性方程组，进而通过这些方程高效的解出每一个$V^\\pi(s)$。\n",
    "\n",
    "我们再给出**最优价值函数（optimal value function）**的定义：\n",
    "\n",
    "$$\n",
    "V^*(s)=\\max_\\pi V^\\pi(s)\\tag{1}\n",
    "$$\n",
    "\n",
    "换句话说，这是在状态$s$时，对于所有可能的策略，能够使折扣奖励的总预期收益最大化的策略$\\pi$下，得到的总收益。当然，这个最优价值函数也有一个对应的贝尔曼方程：\n",
    "\n",
    "$$\n",
    "V^*(s)=R(s)+\\max_{a\\in A}\\gamma\\sum_{s'\\in S}P_{s\\pi(s)}\\left(s'\\right)V^*\\left(s'\\right)\\tag{2}\n",
    "$$\n",
    "\n",
    "第一项没有变，仍旧是策略的立即奖励；而第二项是对于所有可能的动作，在选择执行动作$a$后，使获得的未来折扣奖励的预期总收益最大化。\n",
    "\n",
    "这也引出了最佳策略的定义$\\pi^*:S\\to A$为：\n",
    "\n",
    "$$\n",
    "\\pi^*(s)=\\arg\\max_{a\\in A}\\sum_{s'\\in S}P_{sa}\\left(s'\\right)V^*\\left(s'\\right)\\tag{3}\n",
    "$$\n",
    "\n",
    "也就是说$\\pi^*(s)$给了我们一个使$(2)$式中使$\\max$部分取到最大值时$a$的取值（我们这里省略了$\\gamma$，因为对最大化的式子来说它只是一个常数）。（忘记$\\max f(x)$与$\\arg\\max f(x)$区别的话可以参考[What is the difference between $\\arg\\max$ and $\\max$?](http://math.stackexchange.com/questions/312012/what-is-the-difference-between-arg-max-and-max)，简单地说$\\max f(x)$返回函数$f(x)$的极值，而$\\arg\\max f(x)$返回函数取到极值点时参数的取值。）\n",
    "\n",
    "那么，对于所有状态$s$以及所有策略$\\pi$，有：\n",
    "\n",
    "$$\n",
    "V^*(s)=V^{\\pi^*}(s)\\geq V^\\pi(s)\n",
    "$$\n",
    "\n",
    "前面的等式说明$V^{\\pi^*}$（即策略$\\pi^*$的价值函数）等于“对于所有状态$s$取到最优的价值函数$V^*$”。后面的不等式说明$\\pi^*$的值比任何策略的值都要大。也就是说$(3)$式定义的$\\pi^*$是最优策略。\n",
    "\n",
    "注意到$\\pi^*$的一个有趣的属性——它是对于所有状态$s$的最优策略。这并不是说，如果从状态$s$开始，就有一个针对$s$的最优策略；如果从状态$s'$开始，就会采取别的针对$s'$的最优策略。而是说令$(1)$式取到最大值的的$\\pi^*$是对所有状态$s$而言的，这意味着不论MDP从什么状态开始，我们都可以用策略$\\pi^*$。\n",
    "\n",
    "## 2. 值迭代与策略迭代\n",
    "\n",
    "接下来我们将介绍两种针对有限状态的MDP的高效算法。我们假设MDP具有有限的状态以及有限的动作空间（$\\lvert S\\rvert\\leq\\infty,\\lvert A\\rvert\\leq\\infty$）。\n",
    "\n",
    "要介绍的第一个算法是**值迭代（value iteration）**：\n",
    "1. 对于每个状态$s$，以$V(s):=0$初始化；（相当于在内存中创建了大小为$\\lvert S\\rvert$的状态向量，用$\\vec0$初始化。）\n",
    "2. 重复直到收敛：`{`\n",
    "\n",
    "  * 对于每个状态，更新$\\displaystyle V(s):=R(s)+\\max_{a\\in A}\\gamma\\sum_{s'\\in S}P_{sa}\\left(s'\\right)V\\left(s'\\right)$。（注意，本小节介绍的算法都是在状态转换概率$P_{sa}$及奖励函数$R$已知时使用的。）\n",
    "  \n",
    "  `}`\n",
    "\n",
    "这个算法可以看做是不停的尝试使用贝尔曼方程更新价值函数。对于内部的循环，有两种更新方法：\n",
    "1. 可以先计算每个$s$的$V(s)$，然后用这些新算出来的值代替所有的旧值（用新计算的状态向量更新第1步中的状态向量），这也称作**同步（synchronous）**更新。这种方法可以看做是一种“Bellman backup operator”的实现，用对当前价值函数的估值映射到新的估值上；\n",
    "2. 可以遍历状态$s$（按某种固定顺序），每次更新一个值（每次$V(s)$的计算都从状态向量中取最新的分量，得出结果后直接更新状态向量中想要的分量），这也称作**异步（asynchronous）**更新。\n",
    "\n",
    "不论是使用同步还是异步的更新（通常异步会快一点），对值的迭代可以使得$V$最终收敛于$V^*$。得到$V^*$后就可以使用$(3)$式计算最优策略了。\n",
    "\n",
    "在MDP中还有一种用于计算最优策略的标准算法，称为**策略迭代（policy iteration）**：\n",
    "1. 随机初始化策略$\\pi$；\n",
    "2. 重复直到收敛：`{`\n",
    "  * (a) 令$V:=V^\\pi$；（通过贝尔曼方程组求解。）\n",
    "  * (b) 对每一个状态$s$，令$\\displaystyle\\pi(s):=\\arg\\max_{a\\in A}\\sum_{s'\\in S}P_{sa}\\left(s'\\right)V\\left(s'\\right)$。\n",
    "  \n",
    "  `}`\n",
    "\n",
    "内部的循环会不停的计算当前策略下的价值函数，然后使用当前的价值函数更新策略。（(b)步骤中求$\\pi$也称为**关于$V$的贪心策略（greedy with respect to $V$）**）。而步骤(a)可以通过求解贝尔曼方程组得到——也就是在策略$\\pi$给定时，有$\\lvert S\\rvert$个变量$V^\\pi(s)$，以及$\\lvert S\\rvert$个方程的线性方程组。\n",
    "\n",
    "在有限的迭代之后，$V$会收敛于$V^*$且$\\pi$会收敛于$\\pi^*$。\n",
    "\n",
    "值迭代与策略迭代都是求解MDP的标准算法，目前大家并没有对于哪个算法更好达成共识。对于小型MDP，策略迭代通常更快并会在很少的几步迭代之后收敛。然而在求解具有很多状态的大型MDP时，对$V^\\pi$的求解将涉及到解大型线性方程组，这通常会很困难。在这种情况下，值迭代可能会是更好的选择。也正是因为这个原因，我们通常在现实中更喜欢使用值迭代。\n",
    "\n",
    "## 3. MDP的模型估计\n",
    "\n",
    "到目前为止，我们讨论的MDP以及MDP算法都是建立在状态转换概率以及奖励函数已知的前提下。但是，在现实问题中，状态转换概率和奖励函数可能不是明确的已知条件，那么我们必须从数据中估计这些量。（通常$S,A,\\gamma$都是已知的。）\n",
    "\n",
    "比如在倒置钟摆问题中（见[问题集4](http://cs229.stanford.edu/materials/ps4.pdf)），我对MDP进行了一系列的试验，过程如下：\n",
    "\n",
    "$$\n",
    "\\require{AMScd}\n",
    "\\begin{CD}\n",
    "s_0^{(1)} @>{a_0^{(1)}}>> s_1^{(1)} @>{a_1^{(1)}}>> s_2^{(1)} @>{a_2^{(1)}}>> s_3^{(1)} @>{a_3^{(1)}}>> \\cdots\n",
    "\\end{CD}\\\\\n",
    "\\begin{CD}\n",
    "s_0^{(2)} @>{a_0^{(2)}}>> s_1^{(2)} @>{a_1^{(2)}}>> s_2^{(2)} @>{a_2^{(2)}}>> s_3^{(2)} @>{a_3^{(2)}}>> \\cdots\n",
    "\\end{CD}\n",
    "$$\n",
    "\n",
    "这里$s_i^{(j)}$表示在第$j$此试验中第$i$步的状态，而$a_i^{(j)}$表示在第$j$此试验中第$i$步的状态下执行的动作。在实际操作中，除非MDP过程终止（比如在倒置钟摆问题中，在钟摆倒下时MDP终止），否则每一个试验将一直运行下去；或者可能运行很多但是有限步后停下。有了这些MDP试验的“经验”，我们就可以推导出状态转换概率的最大似然估计：\n",
    "\n",
    "$$\n",
    "P_{sa}\\left(s'\\right)=\\frac{\\textrm{在状态}s\\textrm{下执行动作}a\\textrm{得到状态}s'\\textrm{的次数}}{\\textrm{状态}s\\textrm{下执行动作}a\\textrm{的次数}}\\tag{4}\n",
    "$$\n",
    "\n",
    "如果在应用上式时得到的比值为$\\displaystyle\\frac{0}{0}$时（比如在状态$s$下从未执行过动作$a$时），我们可以简单的将$P_{sa}\\left(s'\\right)$估计为$\\displaystyle\\frac{1}{\\lvert S\\rvert}$。（即将$P_{sa}$估计为在所有状态上的均匀分布。）\n",
    "\n",
    "注意到如果能够在MDP中得到更多的“经验”（即进行更多次试验），我们就可以使用一种更高效的方法，拿新“经验”的数据对状态转换概率的估计进行更新：我们可以记录$(4)$式中分子和分布的计数值，当我们拿到新的试验数据时，不停的累加分子和分母的计数值就可以了。最后只需计算一次比值，就能够得到$P_{sa}$。\n",
    "\n",
    "使用类似的过程也可以处理奖励函数$R$未知的情况，我们可以使用所有状态$s$下的平均奖励来估计状态$s$预期即时奖励$R(s)$。\n",
    "\n",
    "在得到MDP的模型之后，我们可以带入估计出的状态转换概率和奖励函数，使用值迭代或策略迭代解出MDP的最佳策略。举个例子，联合使用模型估计和值迭代，在状态转换概率未知的情况下学习MDP，可以使用下面的算法：\n",
    "1. 随机初始化策略$\\pi$；\n",
    "2. 重复 `{`\n",
    "\n",
    "  * (a) 在MDP中按照策略$\\pi$执行一些试验；\n",
    "  * (b) 使用在MDP的试验中积累的“经验”，更新对$P_{sa}$的估计（如果可以的话也更新$R$）；\n",
    "  * (c) 使用估计的状态转换概率和奖励函数，应用值迭代，估计新的价值函数$V$；\n",
    "  * (d) 使用关于$V$的贪心策略更新$\\pi$；\n",
    "\n",
    "  `}`\n",
    "\n",
    "注意到对于上面这个特定算法，有一个能够大幅优化性能的方法：在内部循环使用值迭代的步骤里，我们不像值迭代定义的那样用$V=0$做初始化，而使用算法上一步迭代得到的解做初始化，这样就给值迭代过程了一个更好的起点，能够使值迭代过程收敛更加迅速。"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}