【动手学强化学习】第三章 马尔可夫决策过程 知识点总结

本章知识点：

马尔可夫性质

• 随机过程：随时间连续变化的随机现象（随机现象是状态的变化过程）
• 在某时刻$t$的状态$S_t$通常取决于时刻之前的状态，我们将已知历史信息$(S_1,\dots ,S_t)$时，下一个时刻状态为$S_{t+1}$的概率表示成$P(S_{t+1}|S_1,S_2,\dots ,S_t)$
• 马尔科夫性质：当且仅当某时刻的状态只取决于上一时刻的状态时，一个随机过程被称为具有马尔可夫性质（Markov property）$P(S_{t+1}|S_t)=P(S_{t+1}|S_1,S_2,\dots ,S_t)$
• 当前状态是未来的充分统计量，即下一个状态只取决于当前状态，而不会受到过去状态的影响
• 具有马尔可夫性质并不代表这个随机过程就和历史完全没有关系。因为虽然$t+1$时刻的状态只与t时刻的状态有关，但是$t$时刻的状态其实包含了之前时刻的状态的信息，通过这种链式的关系，历史的信息被传递到了现在。马尔可夫性可以大大简化运算，因为只要当前状态可知，所有的历史信息都不再需要了，利用当前状态信息就可以决定未来。

马尔可夫过程(也称为马尔可夫链）

• 马尔可夫过程（Markov process）指具有马尔可夫性质的随机过程，也被称为马尔可夫链（Markov chain）
• 用元组$<S,P>$描述一个马尔可夫过程，其中是$S$有限数量的状态集合，$P$是状态转移矩阵（state transition matrix)
  
• 举个栗子🌰
  
  
• 采样的定义：给定一个马尔可夫过程，从某个状态出发，根据状态转移矩阵生成一个状态序列（episode），这个步骤也被叫做采样（sampling）
• “episode”（回合）是指智能体从一个特定的初始状态开始，到达到某个终止状态或满足特定条件结束的一系列连续的时间步（timesteps），episode强调从开始状态到终止状态的整个过程。

马尔可夫奖励过程

• 在马尔可夫过程的基础上加入奖励函数$r$和折扣因子$\gamma$ --》得到马尔科夫奖励过程$<S,P,r,\gamma >$
  

• 回报:在一个马尔可夫奖励过程中，从第$t$时刻状态$S_t$开始，直到终止状态时，所有奖励的衰减之和称为回报（Return）

• 回报计算公式： $G_t=R_t+\gamma R_{t+1}+{\gamma }^2{R}_{t+2}+\dots ={\sum }_{k=0}^{\infty }{\gamma }^k{R}_{t+k}$
  $$\begin{array}{rl}{G}_t& =R_t+\gamma R_{t+1}+{\gamma }^2{R}_{t+2}+\dots \\ & =R_t+\gamma (R_{t+1}+\gamma R_{t+2}+\dots )\\ & =R_t+\gamma G_{t+1}\end{array}$$
  

# 以书上的3-2马尔科夫奖励过程 定义计算回报的函数
import numpy as np
np.random.seed(0)
# 定义状态转移概率矩阵P
P = [
    [0.9, 0.1, 0.0, 0.0, 0.0, 0.0],
    [0.5, 0.0, 0.5, 0.0, 0.0, 0.0],
    [0.0, 0.0, 0.0, 0.6, 0.0, 0.4],
    [0.0, 0.0, 0.0, 0.0, 0.3, 0.7],
    [0.0, 0.2, 0.3, 0.5, 0.0, 0.0],
    [0.0, 0.0, 0.0, 0.0, 0.0, 1.0],
]
P = np.array(P)

rewards = [-1, -2, -2, 10, 1, 0] #定义奖励函数
gamma = 0.5 #定义折扣因子
# 给定一条序列,计算从某个索引（起始状态）开始到序列最后（终止状态）得到的回报
def compute_return(start_index, chain, gamma):
    G = 0
    for i in reversed(range(start_index,len(chain))): # reversed函数来反转上面生成的序列。例如，如果上面的序列是[3, 4]，那么反转后就是[4, 3] 
        G = gamma * G + rewards[chain[i]-1] #倒着计算的 从后往前 因reversed
    return G

# 一个状态序列,s1-s2-s3-s6
chain = [1, 2, 3, 6] 
start_index = 0   #因这里从S1开始，因此start_index = 0
G = compute_return(start_index, chain, gamma)
print("根据本序列计算得到回报为：%s。" % G)

• 价值value：一个状态的期望回报（即从这个状态出发的未来累积奖励的期望）被称为这个状态的价值

• 所有状态的价值就组成了状态价值函数（value function），输入为某个状态，输出为这个状态的价值
  $$\begin{array}{rl}V(s)& =E[G_t|]S_t=s]\\ & =E[R_t+\gamma R_{t+1}+{\gamma }^2{R}_{t+2}+\dots |S_t=s]\\ & =E[R_t+\gamma (R_{t+1}+\gamma R_{t+2}+\dots )|S_t=s]\\ & =E[R_t+\gamma G_{t+1}|S_t=s]]\\ & =E[R_t+\gamma V(S_{t+1})|S_t=s]]\end{array}$$

• 在上式的最后一个等号中，一方面，即时奖励的期望正是奖励函数的输出, 即$E[R_t|S_t=s]=r(s)$；一方面，等式中剩余部分$E[\gamma V(S_{t+1})|S_t=s]]$可以根据从状态$s$出发的转移概率得到

• 马尔可夫奖励过程的贝尔曼方程 $V(s)=r(s)+\gamma {\sum }_{s^{\prime }\in S}P(s^{\prime }|s)V(s^{\prime })$

• 请注意这是在马尔可夫奖励过程下，只涉及状态，状态转移矩阵和奖励函数$<S,P,r,\gamma >$，还没有引入动作$A$和策略$\pi$

• 通过矩阵运算求解析解，用于求解马尔可夫奖励过程的价值函数（仅适用于小规模）
  

• 除了使用动态规划算法，马尔可夫奖励过程中的价值函数也可以通过蒙特卡洛方法估计得到，将在本篇博文的最后介绍该方法

求解价值函数的解析解方法

def compute(P, rewards, gamma, states_num):
    ''' 利用贝尔曼方程的矩阵形式计算解析解,states_num是MRP的状态数 '''
    rewards = np.array(rewards).reshape((-1, 1))  #将rewards写成列向量形式
    value = np.dot(np.linalg.inv(np.eye(states_num, states_num) - gamma * P),
                   rewards)
    return value
V = compute(P, rewards, gamma, 6)
print("MRP中每个状态价值分别为\n", V)

马尔可夫决策过程

• 前两章所说的强化学习中的环境一般就是一个马尔可夫决策过程

• 与多臂老虎机问题不同，马尔可夫决策过程包含状态信息以及状态之间的转移机制

• 如果要用强化学习去解决一个实际问题，第一步要做的事情就是把这个实际问题抽象为一个马尔可夫决策过程，也就是明确马尔可夫决策过程的各个组成要素

• 马尔可夫过程和马尔可夫奖励过程都是自发改变的随机过程；而如果有一个外界的“刺激”来共同改变这个随机过程，就有了马尔可夫决策过程（Markov decision process，MDP），将这个来自外界的刺激称为智能体（agent）的动作

• 马尔科夫决策过程：在马尔可夫奖励过程（MRP）的基础上加入动作，就得到了马尔可夫决策过程（MDP），由五元组组成$<S,A,P,r,\gamma >$
  

• MDP 的定义中，我们不再使用类似 MRP 定义中的状态转移矩阵方式，而是直接表示成了状态转移函数。这样做一是因为此时状态转移与动作也有关，变成了一个三维数组，而不再是一个矩阵（二维数组）；二是因为状态转移函数更具有一般意义，例如，如果状态集合不是有限的，就无法用数组表示，但仍然可以用状态转移函数表示
  

• 智能体根据当前状态$S_t$选择动作$A_t$；对于状态$S_t$和动作$A_t$，MDP 根据奖励函数和状态转移函数得到$S_{S_{t+1}}$和$R_t$并反馈给智能体。

• 智能体的目标是最大化得到的累计奖励。智能体根据当前状态从动作的集合中选择一个动作的函数，被称为策略。

策略

• 智能体根据当前状态从动作的集合A中选择一个动作的函数，被称为策略 通常使用$\pi$表示
• 策略是一个函数 $\pi (a|s)=P(A_t=a|S_t=s)$，表示在输入状态 $s$ 情况下采取动作$a$的概率
• 当一个策略是确定性策略（deterministic policy）时，它在每个状态时只输出一个确定性的动作，即只有该动作的概率为 1，其他动作的概率为 0；
• 当一个策略是随机性策略（stochastic policy）时，它在每个状态时输出的是关于动作的概率分布，然后根据该分布进行采样就可以得到一个动作
• 在 MDP 中，由于马尔可夫性质的存在，策略只需要与当前状态有关，不需要考虑历史状态

马尔可夫决策过程的状态价值函数

• $V^{\pi }(s)$表示在 MDP 中基于策略的状态价值函数（state-value function），定义为从状态出发遵循策略能获得的期望回报
  $$\begin{array}{r}{V}^{\pi }(s)=E_{\pi }[G_t|S_t=s]\end{array}$$
• 用于评估当前局势的好坏，靠期望回报衡量 能够体现当前局势好不好，好-期望大，不好-期望小
• 价值函数与策略有关

马尔可夫决策过程的动作价值函数

• 不同于 MRP，在 MDP 中，由于动作的存在，我们额外定义一个动作价值函数（action-value function），

• 用$Q^{\pi }(s,a)$表示在 MDP 遵循策略$\pi$时，对当前状态$s$执行动作$a$得到的期望回报
  $$\begin{array}{r}{Q}^{\pi }(s,a)=E_{\pi }[G_t|S_t=s,A_t=a]\end{array}$$

• 用于评估在状态$s$下执行动作$a$的好坏

状态价值函数与动作价值函数的关系

• 在使用策略$\pi$中，状态$s$的价值等于在该状态下基于策略$\pi$采取所有动作的概率与相应的价值相乘再求和的结果（期望的定义）：
  $$\begin{array}{r}{V}^{\pi }(s)=\sum _{a\in A}\pi (a|s)Q^{\pi }(s,a)\end{array}$$

• 推导过程
  $$\begin{array}{rl}{V}^{\pi }(s)& =E_{\pi }[G_t|S_t=s]\\ & =\sum _{a\in A}\pi (a|s)E_{\pi }[G_t|S_t=s,A_t=a]\\ & =\sum _{a\in A}\pi (a|s)Q^{\pi }(s,a)\end{array}$$

• 使用策略$\pi$时，状态$s$下采取动作$a$的价值等于即时奖励加上经过衰减后的所有可能的下一个状态的状态转移概率与相应的价值的乘积
  $$\begin{array}{r}{Q}^{\pi }(s,a)=r(s,a)+\gamma \sum _{s^{\prime }\in S}P(s^{\prime }|s,a)V^{\pi }(s^{\prime })\end{array}$$

• 二者关系的理解建议结合下面的图理解，白色代表状态，黑色代表动作

• 从$s\to a$，由$\pi (a|s)$决定

• 从$s,a\to s^{\prime }$，由$P(s^{\prime }|s,a)$决定

贝尔曼期望方程

• 贝尔曼期望方程：对动作价值函数和状态价值函数分别进行推导当前时刻$t$与下一时刻$t+1$的关系，类似推导回报$G_t$与$G_{t+1}$的关系

• 贝尔曼方程（Bellman Equation），也称为贝尔曼期望方程

• $V^{\pi }(s)$推导过程
  $$\begin{array}{rl}{V}^{\pi }(s)& =E_{\pi }[G_t|]S_t=s]\\ & =E_{\pi }[R_t+\gamma R_{t+1}+{\gamma }^2{R}_{t+2}+\dots |S_t=s]\\ & =E_{\pi }[R_t+\gamma (R_{t+1}+\gamma R_{t+2}+\dots )|S_t=s]\\ & =E_{\pi }[R_t+\gamma G_{t+1}|S_t=s]\\ & =E_{\pi }[R_t+\gamma V_{\pi }(S_{t+1})|S_t=s]]\\ & =\sum _{a_t\in A}\pi (a_t|s_t)[R(s_t,a_t)+\gamma \sum _{s_{t+1}\in S}P(s_{t+1}|s_t,a_t)V_{\pi }(s_{t+1})]\end{array}$$

• $Q^{\pi }(s,a)$推导过程
  $$\begin{array}{rl}{Q}^{\pi }(s,a)& =E_{\pi }[G_t|S_t=s,A_t=a]\\ & =E_{\pi }[R_t+\gamma R_{t+1}+{\gamma }^2{R}_{t+2}+\dots |S_t=s,A_t=a]\\ & =E_{\pi }[R_t+\gamma (R_{t+1}+\gamma R_{t+2}+\dots )|S_t=s,A_t=a]\\ & =E_{\pi }[R_t+\gamma G_{t+1}|S_t=s,A_t=a]]\\ & =R(s,a)+\gamma E_{\pi }[G_{t+1}|S_t=s,A_t=a]\\ & =R(s,a)+\gamma E_{\pi }[\sum _{s_{t+1}\in S}P(s_{t+1}|s,a)G_{t+1}|S_t=s,A_t=a]这里考虑状态转移函数\\ & =R(s,a)+\gamma \sum _{s_{t+1}\in S}P(s_{t+1}|s,a)E_{\pi }[G_{t+1}|S_t=s,A_t=a]\\ & =R(s,a)+\gamma \sum _{s_{t+1}\in S}P(s_{t+1}|s,a)\sum _{a_{t+1}\in A}\pi (a_{t+1}|s_{t+1})Q^{\pi }(s_{t+1},a_{t+1})这里考虑策略\end{array}$$

• 上述推导比较复杂，另一种理解方式：根据状态价值函数和动作价值函数的关系，互相带入，具体如下图所示
  
  

• 编写代码来表示图 3-4 中的马尔可夫决策过程，并定义两个策略

  第一个策略是一个完全随机策略，即在每个状态下，智能体会以同样的概率选取它可能采取的动作。
  第二个策略是一个提前设定的一个策略。

S = ["s1", "s2", "s3", "s4", "s5"]  # 状态集合
A = ["保持s1", "前往s1", "前往s2", "前往s3", "前往s4", "前往s5", "概率前往"]  # 动作集合
# 状态转移函数
P = {
    "s1-保持s1-s1": 1.0,
    "s1-前往s2-s2": 1.0,
    "s2-前往s1-s1": 1.0,
    "s2-前往s3-s3": 1.0,
    "s3-前往s4-s4": 1.0,
    "s3-前往s5-s5": 1.0,
    "s4-前往s5-s5": 1.0,
    "s4-概率前往-s2": 0.2,
    "s4-概率前往-s3": 0.4,
    "s4-概率前往-s4": 0.4,
}
# 奖励函数
R = {
    "s1-保持s1": -1,
    "s1-前往s2": 0,
    "s2-前往s1": -1,
    "s2-前往s3": -2,
    "s3-前往s4": -2,
    "s3-前往s5": 0,
    "s4-前往s5": 10,
    "s4-概率前往": 1,
}
gamma = 0.5  # 折扣因子
MDP = (S, A, P, R, gamma)

# 策略1,随机策略
Pi_1 = {
    "s1-保持s1": 0.5,
    "s1-前往s2": 0.5,
    "s2-前往s1": 0.5,
    "s2-前往s3": 0.5,
    "s3-前往s4": 0.5,
    "s3-前往s5": 0.5,
    "s4-前往s5": 0.5,
    "s4-概率前往": 0.5,
}
# 策略2
Pi_2 = {
    "s1-保持s1": 0.6,
    "s1-前往s2": 0.4,
    "s2-前往s1": 0.3,
    "s2-前往s3": 0.7,
    "s3-前往s4": 0.5,
    "s3-前往s5": 0.5,
    "s4-前往s5": 0.1,
    "s4-概率前往": 0.9,
}


# 把输入的两个字符串通过“-”连接,便于使用上述定义的P、R变量
def join(str1, str2):
    return str1 + '-' + str2

MDP马尔可夫决策过程转换为MRP马尔可夫奖励过程

• 给定一个 MDP 和一个策略，可以将策略的动作选择进行边缘化（marginalization)，就可以得到没有动作的 MRP

• 对于某一个状态，我们根据策略所有动作的概率进行加权，得到的奖励和就可以认为是一个 MRP 在该状态下的奖励，以加权的奖励和作为MRP的奖励
  $$\begin{array}{r}{r}^{\prime }(s)=\sum _{a\in A}\pi (a|s)r(s,a)\end{array}$$

• 同理，计算采取动作的概率与使$s$转移到$s^{\prime }$的概率的乘积，再将这些乘积相加，其和就是一个 MRP 的状态$s$从转移至$s^{\prime }$的概率
  $$\begin{array}{r}{P}^{\prime }(s^{\prime }|s)=\sum _{a\in A}\pi (a|s)P(s^{\prime }|s,a)\end{array}$$

• 根据上面的式子，将MDP转化为了MRP $<S,P^{\prime },r^{\prime },\gamma >$

根据价值函数的定义可以发现,转化前的 MDP 的状态价值函数和转化后的 MRP 的价值函数是一样的

• 将MDP转化为MRP后，可以用 MRP 中计算价值函数的解析解来计算这个 MDP 中该策略的状态价值函数
  计算用随机策略（也就是代码中的Pi_1）时的状态价值函数

gamma = 0.5
# 转化后的MRP的状态转移矩阵
P_from_mdp_to_mrp = [
    [0.5, 0.5, 0.0, 0.0, 0.0],
    [0.5, 0.0, 0.5, 0.0, 0.0],
    [0.0, 0.0, 0.0, 0.5, 0.5],
    [0.0, 0.1, 0.2, 0.2, 0.5],
    [0.0, 0.0, 0.0, 0.0, 1.0],
]
P_from_mdp_to_mrp = np.array(P_from_mdp_to_mrp)
R_from_mdp_to_mrp = [-0.5, -1.5, -1.0, 5.5, 0]

V = compute(P_from_mdp_to_mrp, R_from_mdp_to_mrp, gamma, 5)
print("MDP中每个状态价值分别为\n", V)

蒙特卡洛方法

• 用蒙特卡洛方法的好处在于我们不需要知道 MDP 的状态转移函数和奖励函数，它可以得到一个近似值，并且采样数越多越准确

• 蒙特卡洛方法（Monte-Carlo methods）也被称为统计模拟方法，是一种基于概率统计的数值计算方法

• 运用蒙特卡洛方法时，我们通常使用重复随机抽样，然后运用概率统计方法来从抽样结果中归纳出我们想求的目标的数值估计

• 使用蒙特卡洛估计圆面积的例子：
  

• 使用蒙特卡洛方法来估计一个策略在一个马尔可夫决策过程中的状态价值函数：

• 一个状态的价值是它的期望回报，那么一个很直观的想法就是用策略在 MDP 上采样很多条序列，计算从这个状态出发的回报再求其期望就可以了
  $$\begin{array}{r}{V}^{\pi }(s)=E_{\pi }[G_t|S_t=s]\approx \frac{1}{N}\sum _{i=1}^N{G}_t^{(i)}\end{array}$$

• 个人补充内容：首次访问蒙特卡洛FVMC方法 和 每次访问蒙特卡EVMC方法
  

• 首次访问蒙特卡洛FVMC方法主要包含两个步骤，①产生一个回合的完整轨迹，②遍历轨迹计算每个状态的回报。只在第一次遍历到某个状态时会记录并计算对应的回报，而每次访问蒙特卡EVMC方法中不会忽略同一状态的多个回报

• 
  

• 蒙特卡洛方法有一定的局限性，即只适用于有终止状态的马尔可夫决策过程

• 假设我们现在用策略$\pi$从状态$s$开始采样序列，据此来计算状态价值。我们为每一个状态维护一个计数器和总回报，计算状态价值的具体过程如下所示
  

• 使用增量式更新 估计状态价值的推导过程, $V(s)$为状态$s$的价值， $M(s)$为状态$s$的总回报，$N(s)$为采样到状态$s$的总计数，根据上述可得，$V(s)=\frac{M(s)}{N(s)}$
  $$\begin{array}{rl}V(s)& =\frac{M(s)}{N(s)}\\ & =\frac{1}{N(s)}M(s)\\ & =\frac{1}{N(s)}[G+M(s)]这里的M(s)是减去最新一次遇到状态s的总回报\\ & =\frac{1}{N(s)}[G+(N(s)-1)V(s)]右边的V(s)也是去掉最新一次状态s价值的估计\\ & =V(s)+\frac{1}{N(s)}(G-V(s))\end{array}$$

• 定义一个采样函数，采样函数需要遵守状态转移矩阵和相应的策略，每次将(s,a,r,s_next)元组放入序列中，直到到达终止序列

def sample(MDP, Pi, timestep_max, number):
    ''' 采样函数,策略Pi,限制最长时间步timestep_max,总共采样序列数number '''
    S, A, P, R, gamma = MDP
    episodes = []
    for _ in range(number):
        episode = []
        timestep = 0
        s = S[np.random.randint(4)]  # 随机选择一个除s5以外的状态s作为起点
        # 当前状态为终止状态或者时间步太长时,一次采样结束
        while s != "s5" and timestep <= timestep_max:
            timestep += 1
            rand, temp = np.random.rand(), 0
            # 在状态s下根据策略选择动作
            for a_opt in A:
                temp += Pi.get(join(s, a_opt), 0)
                if temp > rand:
                    a = a_opt
                    r = R.get(join(s, a), 0)
                    break
            rand, temp = np.random.rand(), 0
            # 根据状态转移概率得到下一个状态s_next
            for s_opt in S:
                temp += P.get(join(join(s, a), s_opt), 0) # 使用 get 方法获取键为 "s,a" 的值，如果键不存在，返回默认值 0
                if temp > rand:
                    s_next = s_opt
                    break
            episode.append((s, a, r, s_next))  # 把（s,a,r,s_next）元组放入序列中
            s = s_next  # s_next变成当前状态,开始接下来的循环
        episodes.append(episode)
    return episodes


# 采样5次,每个序列最长不超过20步
episodes = sample(MDP, Pi_1, 20, 5)
print('第一条序列\n', episodes[0])
print('第二条序列\n', episodes[1])
print('第五条序列\n', episodes[4])

• 对所有采样序列计算所有状态的价值

def MC(episodes, V, N, gamma):
    for episode in episodes:
        G = 0
        for i in range(len(episode) - 1, -1, -1):  #一个序列从后往前计算   从len(episode)-1开始递减到-1 步长为-1 range 函数的参数分别是起始值、终止值和步长
            (s, a, r, s_next) = episode[i]
            G = r + gamma * G
            N[s] = N[s] + 1
            V[s] = V[s] + (G - V[s]) / N[s]  #增量式更新


timestep_max = 20
# 采样1000次,可以自行修改
episodes = sample(MDP, Pi_1, timestep_max, 1000)
gamma = 0.5
V = {"s1": 0, "s2": 0, "s3": 0, "s4": 0, "s5": 0}
N = {"s1": 0, "s2": 0, "s3": 0, "s4": 0, "s5": 0}
MC(episodes, V, N, gamma)
print("使用蒙特卡洛方法计算MDP的状态价值为\n", V)

占用度量

• 占用度量（occupancy measure）：是关于数据分布的概念，归一化的的占用度量用于衡量一个智能体决策与一个动态环境的交互过程中，采样到一个具体的状态动作对state-action pair的概率分布
  
• 对上面的式子我也不太理解，注意上面的$v^{\pi }(s)$是状态访问分布
• 状态访问概率表示一个策略和 MDP 交互会访问到的状态的分布
• $P_t^{\pi }(s)$表示采取策略$\pi$使得智能体在$t$时刻状态为$s$的概率
• 占用度量有一个很重要的性质：给定两个策略及其与一个动态环境交互得到的两个占用度量，那么当且仅当这两个占用度量相同时，这两个策略相同。也就是说，如果一个智能体的策略有所改变，那么它和环境交互得到的占用度量也会相应改变
• 策略的占用度量的定义：${\rho }^{\pi }(s,a)=(1-\gamma ){\sum }_{t=0}^{\infty }{\gamma }^t{P}_t^{\pi }(s)\pi (a|s)$，表示动作状态对$(s,a)$被访问到的概率
• 占用度量${\rho }^{\pi }(s,a)$是在策略$\pi$下，智能体在状态$s$下采取动作$a$的概率。它是状态和动作对$(s,a)$的访问频率的度量
• 占用度量${\rho }^{\pi }(s,a)$与状态访问分布$v^{\pi }(s)$的关系：
  $$\begin{array}{r}{\rho }^{\pi }(s,a)=v^{\pi }(s)\pi (a|s)\end{array}$$
  
• “合法”占用度量是指存在一个策略使智能体与 MDP 交互产生的状态动作对被访问到的概率
  近似估计占用度量。这里我们采用近似估计，即设置一个较大的采样轨迹长度的最大值，然后采样很多次，用状态动作对出现的频率估计实际概率

def occupancy(episodes, s, a, timestep_max , gamma):
    ''' 计算状态动作对（s,a）出现的频率,以此来估算策略的占用度量 '''
    rho = 0
    total_times = np.zeros(timestep_max) # 记录每个时间步t各被经历过几次
    occur_times = np.zeros(timestep_max) # 记录(s_t,a_t)=(s,a)的次数
    for episode in episodes:
        for i in range(len(episode)):
            (s_opt, a_opt, r, s_next) = episode[i]
            total_times[i] += 1
            if s == s_opt and a == a_opt:
                occur_times[i] += 1
    for i in reversed(range(timestep_max)):
        if total_times[i]:
            rho += gamma**i* occur_times[i] / total_times[i]
    return (1- gamma) * rho

gamma = 0.5
timestep_max = 1000

episodes_1 = sample(MDP, Pi_1, timestep_max, 1000)
episodes_2 = sample(MDP, Pi_2, timestep_max, 1000)
rho_1 = occupancy(episodes_1, "s4", "概率前往", timestep_max, gamma)
rho_2 = occupancy(episodes_2, "s4", "概率前往", timestep_max, gamma)
print(rho_1, rho_2)

最优策略

• 定义策略之间的偏序关系：当且仅当对于任意的状态$s$都有$V^{\pi }(s)\ge V^{{\pi }^{\prime }}(s)$ ，记$\pi \ge {\pi }^{\prime }$

• 于是在有限状态和动作集合的 MDP 中，至少存在一个策略比其他所有策略都好或者至少存在一个策略不差于其他所有策略，这个策略就是最优策略（optimal policy）

• 最优策略可能有很多个，我们都将其表示为${\pi }^{\ast }(s)$

• 最优策略都有相同的状态价值函数，我们称之为最优状态价值函数
  $$\begin{array}{rl}{V}^{\ast }(s)& =\underset{\pi }{\max }{V}^{\pi }(s)，\forall s\in S\end{array}$$

• 最优动作价值函数
  $$\begin{array}{rl}{Q}^{\ast }(s,a)& =\underset{\pi }{\max }{Q}^{\pi }(s,a)，\forall s\in S,a\in A\end{array}$$

• 对比状态价值函数与动作价值函数的关系

• 最优状态价值函数与最优动作价值函数的关系
• 最优状态价值是选择此时使最优动作价值最大的那一个动作时的状态价值
  $$\begin{array}{rl}{V}^{\ast }(s)& =\underset{\pi }{\max }{Q}^{\ast }(s,a)\\ Q^{\ast }(s,a)& =r(s,a)+\gamma \sum _{s^{\prime }\in S}P(s’|s,a)V^{\ast }(s^{\prime })\end{array}$$

贝尔曼最优方程（Bellman optimality equation）

• 对照贝尔曼期望方程理解
  

• 贝尔曼最优方程
  V ∗ ( s ) = max ⁡ a ∈ A { r ( s , a ) + γ ∑ s t + 1 ∈ S P ( s t + 1 ∣ s t , a t ) V ∗ ( s t + 1 ) } Q ∗ ( s , a ) = r ( s , a ) + γ ∑ s t + 1 ∈ S P ( s t + 1 ∣ s , a ) max ⁡ a t + 1 ∈ A Q ∗ ( s t + 1 , a t + 1 ) \begin{aligned} V^*(s) &= \max_{a \in A} \{r(s,a)+ \gamma \sum_{s_{t+1} \in S}P(s_{t+1}|s_t,a_t)V^*(s_{t+1})\} \\ Q^*(s,a) &= r(s,a) +\gamma \sum_{s_{t+1} \in S} P(s_{t+1}|s,a) \max_{a_{t+1} \in A} Q^*(s_{t+1},a_{t+1}) \end{aligned} V∗(s)Q∗(s,a)​=a∈Amax​{r(s,a)+γst+1​∈S∑​P(st+1​∣st​,at​)V∗(st+1​)}=r(s,a)+γst+1​∈S∑​P(st+1​∣s,a)at+1​∈Amax​Q∗(st+1​,at+1​)​

• 强化学习中的环境就是一个马尔可夫决策过程,接下来将要介绍的强化学习算法通常都是在求解马尔可夫决策过程中的最优策略