该书由清华大学李升波教授撰写的,主要面向工业控制领域的研究者和工程师,曾获得2024年度Springer中国新发展奖(China New Development Awards)。全书按照原理剖析、主流算法、典型示例的架构,系统地介绍了用于动态系统决策与控制的强化学习方法。全书共分为11章,内容涵盖了强化学习的基本概念、蒙特卡洛法、时序差分法、动态规划法、函数近似法、策略梯度法、近似动态规划、状态约束的处理和深度强化学习等知识点。书籍及源代码下载网站:书籍及代码链接点这里。
另外,由于本部分篇幅较长,因此分成三篇博客发表。其余两篇博客的地址参考我的系列博客的地址汇总。
文章目录
书籍链接:Reinforcement Learning for Sequential Decision and Optimal Control
本篇博客主要介绍Dynamic Programming(DP),即动态规划,一种Model-Based的Indirect RL方法。DP与动态环境中的随机最优控制问题有紧密的联系。DP通常通过求解Bellman方程(离散时间)或Hamilton-Jacobi-Bellman方程(连续时间)来得到最优策略。对于大部分工程任务来说,以上两个方程是最优解的充要条件。DP的优点在于它的强大的适应性。如线性/非线性系统、离散/连续系统、确定性/随机性系统、可微分/不可微分系统等。但是,DP在实际应用中面临着维度灾难和缺少泛化能力的问题。但是,DP仍然在理论分析最优解的存在以及收敛速度等方面有着重要的作用。
5.1 随机序列决策问题(Stochastic Sequential Decision Problem)
DP的设计包含两个重要的部分: 对于转移的建模、 对于奖励的建模 。因为实际任务的复杂性,对于上述两个问题的回答可能呈现出大相径庭的形式。因此,我们在实际中看到的对于DP的建模可能完全看不出是同一种算法,这对于初学者理解算法是一个很大的困难。因此,书中首先试图给出一种理解DP的统一的视角。
5.1.1 随机环境的建模
有两种方式可以描述随机环境: 状态空间模型(state space model)和概率模型(probabilistic model)。首先来看看状态空间模型。
5.1.1.1 状态空间模型
s t + 1 = f ( s t , a t , ξ t ) s_{t+1}=f(s_t,a_t,\xi_t) st+1=f(st,at,ξt) 这里f(.)表示环境模型,而 ξ t \xi_t ξt表示一个分布已知的随机分布。这里,我们需要做两条假设来保证环境的Markov性质:
- ξ t \xi_t ξt与t时刻之前的噪声 ξ t − 1 , ξ t − 2 , . . . \xi_{t-1},\xi_{t-2},... ξt−1,ξt−2,...及初始状态 s 0 s_0 s0独立。
- 每个动作 a t a_t at的选择是任意的且与历史的状态和随机噪声无关。
有了上述假设后就可证上述状态空间模型具有Markov性质,即下一个状态 s t + 1 s_{t+1} st+1只与当前状态 s t s_t st、当前动作 a t a_t at和随机噪声 ξ t \xi_t ξt有关。
5.1.1.2 概率模型
上述具有Markov性质的状态空间模型可以转化为概率模型。不失一般性地,可以假设噪声
ξ
\xi
ξ为加性噪声。即:
s
t
+
1
=
f
(
s
t
,
a
t
)
+
ξ
t
s_{t+1}=f(s_t,a_t)+\xi_t
st+1=f(st,at)+ξt 再假设噪声
ξ
t
\xi_t
ξt的概率密度函数为
p
ξ
(
ξ
t
;
θ
)
p_{\xi}(\xi_t;\theta)
pξ(ξt;θ),其中
θ
\theta
θ为参数。那么,概率模型可以写为:
P
s
s
′
a
≜
p
ξ
(
ξ
t
;
θ
)
=
p
ξ
(
s
t
+
1
−
f
(
s
t
,
a
t
)
;
θ
)
\mathcal{P}_{ss^{\prime}}^a\triangleq p_\xi(\xi_t;\theta)=p_\xi(s_{t+1}-f(s_t,a_t);\theta)
Pss′a≜pξ(ξt;θ)=pξ(st+1−f(st,at);θ) 这里虽然没有直接给出从当前状态
s
t
s_t
st采取动作
a
t
a_t
at到下一个状态
s
t
+
1
s_{t+1}
st+1的概率,但实际上已经使用了随机噪声的概率密度来表示了这个概率。
为了更形象地说明上述两种模型的关系,我们可以通过一个具体的例子来说明。考虑一个线性的状态空间模型以及加性高斯噪声:
s
t
+
1
=
A
s
t
+
B
a
t
+
ξ
t
s_{t+1}=As_t+Ba_t+\xi_t
st+1=Ast+Bat+ξt 这里假设状态s为n维向量,动作a为m维向量,
A
∈
R
n
×
n
,
B
∈
R
n
×
m
A\in\mathbb{R}^{n\times n},B\in\mathbb{R}^{n\times m}
A∈Rn×n,B∈Rn×m。且原系统能控(回想一下现代控制理论里面的能控性判据,即
r
a
n
k
(
[
B
,
A
B
,
A
2
B
,
.
.
.
,
A
n
−
1
B
]
)
=
n
rank([B,AB,A^2B,...,A^{n-1}B])=n
rank([B,AB,A2B,...,An−1B])=n)。接下来,写出噪声的概率密度函数:
p
ξ
(
ξ
t
)
=
1
(
2
π
)
n
/
2
∣
K
∣
1
/
2
exp
(
−
1
2
(
ξ
t
−
μ
)
T
K
−
1
(
ξ
t
−
μ
)
)
p_\xi(\xi_t)=\frac1{(2\pi)^{n/2}|\mathcal{K}|^{1/2}}\exp\left(-\frac12(\xi_t-\mu)^{\mathrm{T}}\mathcal{K}^{-1}(\xi_t-\mu)\right)
pξ(ξt)=(2π)n/2∣K∣1/21exp(−21(ξt−μ)TK−1(ξt−μ))这是一个n维的高斯噪声,其中
μ
\mu
μ为均值,
K
\mathcal{K}
K为协方差矩阵。那么,接下来把
ξ
t
=
s
t
+
1
−
A
s
t
−
B
a
t
\xi_t=s_{t+1}-As_t-Ba_t
ξt=st+1−Ast−Bat代入上式,就可以得到概率模型:
p
(
s
t
+
1
∣
s
t
,
a
t
)
=
1
(
2
π
)
n
/
2
∣
K
∣
1
/
2
exp
(
−
1
2
z
T
K
−
1
z
)
,
z
=
s
t
+
1
−
(
A
s
t
+
B
a
t
+
μ
)
,
p(s_{t+1}|s_t,a_t)=\frac1{(2\pi)^{n/2}|\mathcal{K}|^{1/2}}\exp\left(-\frac12z^\mathrm{T}\mathcal{K}^{-1}z\right),\\z=s_{t+1}-(As_t+Ba_t+\mu),
p(st+1∣st,at)=(2π)n/2∣K∣1/21exp(−21zTK−1z),z=st+1−(Ast+Bat+μ), 这里
z
∼
N
(
0
,
K
)
z\sim\mathcal{N}(0,\mathcal{K})
z∼N(0,K)。那么,因为
s
t
+
1
=
A
s
t
+
B
a
t
+
μ
+
z
s_{t+1}=As_t+Ba_t+\mu+z
st+1=Ast+Bat+μ+z,所以
s
t
+
1
s_{t+1}
st+1服从高斯分布
s
t
+
1
∼
N
(
(
A
s
t
+
B
a
t
+
μ
)
,
K
)
.
s_{t+1}{\sim}\mathcal{N}\big((As_t+Ba_t+\mu),\mathcal{K}\big).
st+1∼N((Ast+Bat+μ),K).
总结来说,概率模型和状态空间模型都可以描述随机的环境。它们之间的区别如下:
- 如何表示随机性:状态空间模型的随机性是显式表示的,而概率模型的随机性隐含于概率密度之中。
- 是否具有Markov性质:单步的概率模型天然具有Markov性质,而状态空间模型需要通过一些假设来保证Markov性质。
5.1.2 代价函数的建模
序列决策问题的数学建模可以非常多样化。其中一个重要的原因在于对于目标函数/代价函数有着非常不同的建模方式。对于无限时间的序列决策问题,需要精心设计代价函数来保证该函数只能取有限值。两种流行的代价函数是折扣代价函数(discounted cost)和平均代价函数(average cost)。
5.1.2.1 平均代价函数
回顾我们在第二章中介绍的average return:
G
t
a
v
g
=
lim
T
→
∞
1
T
∑
i
=
0
T
−
1
r
t
+
i
G_t^{avg}=\lim_{T\to\infty}\frac{1}{T}\sum_{i=0}^{T-1}r_{t+i}
Gtavg=T→∞limT1i=0∑T−1rt+i 我们可以据此定义average cost函数:
G
a
v
g
(
π
)
=
lim
T
→
∞
1
T
E
{
∑
i
=
0
T
−
1
r
t
+
i
}
G_{\mathrm{avg}}(\pi)=\lim_{T\to\infty}\frac{1}{T}\mathbb{E}\left\{\sum_{i=0}^{T-1}r_{t+i}\right\}
Gavg(π)=T→∞limT1E{i=0∑T−1rt+i} 可以发现,两式的区别在于第二式在第一式的基础上套了一个期望,该期望表示在策略
π
\pi
π下求期望。上式中并没有写出来。通常,我们使用另一种更常见的average cost函数作为优化目标:
J
a
v
g
(
π
)
=
∑
d
π
(
s
)
v
γ
π
(
s
)
,
J_{\mathrm{avg}}(\pi)=\sum d_\pi(s)v_\gamma^\pi(s),
Javg(π)=∑dπ(s)vγπ(s),
v
γ
π
(
s
)
v_\gamma^\pi(s)
vγπ(s)是带衰减因子
γ
\gamma
γ的状态值函数(即状态值函数的标准定义:
v
π
(
s
)
=
d
e
f
E
π
{
G
t
∣
s
}
=
E
π
{
∑
i
=
0
+
∞
γ
i
r
t
+
i
∣
s
t
=
s
}
v^{\pi}(s)\stackrel{\mathrm{def}}{=}\mathbb{E}_{\pi}\left\{G_{t}\left|s\right\}\right.=\mathbb{E}_{\pi}\left\{\sum_{i=0}^{+\infty}\gamma^{i}r_{t+i}\left|s_{t}\right.=s\right\}
vπ(s)=defEπ{Gt∣s}=Eπ{∑i=0+∞γirt+i∣st=s},详见第二单元的博客)。这里状态值函数被使用
d
π
(
s
)
d_\pi(s)
dπ(s)进行加权,因为我们这里的代价函数并不是针对某一个特定的状态来说的,而是针对所有的状态来说的,因此自然要根据哥哥状态出现的概率进行加权。这里的
d
π
(
s
)
d_\pi(s)
dπ(s)表示在策略
π
\pi
π下状态s的分布,也被称为weighting function。但是,这里的
J
a
v
g
(
π
)
J_{\mathrm{avg}}(\pi)
Javg(π)在优化时因为权重
d
π
(
s
)
d_\pi(s)
dπ(s)和状态值函数
v
γ
π
(
s
)
v_\gamma^\pi(s)
vγπ(s)都与当前策略
π
\pi
π有关,所以优化起来非常困难。回顾第二单元讲过的RL问题的标准Formulation:
max
π
/
min
π
J
(
π
)
=
E
s
∼
d
i
n
i
t
(
s
)
{
v
π
(
s
)
}
\max_{\pi} / \min_{\pi} J(\pi)=\mathbb{E}_{s\sim d_{\mathrm{init}} (s)} \{v^{\pi}(s)\}
πmax/πminJ(π)=Es∼dinit(s){vπ(s)} 标准的RL问题的优化目标函数只与初始状态的分布有关,而不是时刻和当前状态相关联。那么,为什么我们这里要使用如此的定义方式呢?这是因为如此定义的
J
a
v
g
(
π
)
J_{\mathrm{avg}}(\pi)
Javg(π)与上述
G
a
v
g
(
π
)
G_{\mathrm{avg}}(\pi)
Gavg(π)有着紧密的联系,而后者有一个很好的性质,就是与策略无关,对于任何策略的形式都是一样的。下面给出这个定理及其证明。
定理1:
J
a
v
g
(
π
)
=
1
1
−
γ
G
a
v
g
(
π
)
\text{定理1:}J_{\mathrm{avg}}(\pi)=\frac1{1-\gamma}G_{\mathrm{avg}}(\pi)
定理1:Javg(π)=1−γ1Gavg(π) 定理的证明如下:
J
a
v
g
(
π
)
=
∑
d
π
(
s
)
v
γ
π
(
s
)
=
∑
d
π
(
s
)
∑
a
∈
A
π
(
a
∣
s
)
{
∑
s
′
∈
S
P
s
s
′
a
(
r
s
s
′
a
+
γ
v
γ
π
(
s
′
)
)
}
=
G
a
v
g
(
π
)
+
γ
∑
d
π
(
s
)
∑
a
∈
A
π
(
a
∣
s
)
∑
s
′
∈
S
P
s
s
′
a
v
γ
π
(
s
′
)
=
G
a
v
g
(
π
)
+
γ
∑
s
′
∈
S
d
π
(
s
′
)
v
γ
π
(
s
′
)
=
G
a
v
g
(
π
)
+
γ
J
a
v
g
(
π
)
=
G
a
v
g
(
π
)
+
γ
G
a
v
g
(
π
)
+
γ
2
J
a
v
g
(
π
)
=
G
a
v
g
(
π
)
+
γ
G
a
v
g
(
π
)
+
γ
2
G
a
v
g
(
π
)
+
⋯
=
1
1
−
γ
G
a
v
g
(
π
)
.
\begin{aligned} J_{\mathrm{avg}}(\pi)& =\sum d_\pi(s)v_\gamma^\pi(s) \\ &=\sum d_\pi(s)\sum_{a\in\mathcal{A}}\pi(a|s)\left\{\sum_{s^{\prime}\in\mathcal{S}}\mathcal{P}_{ss^{\prime}}^a\left(r_{ss^{\prime}}^a+\gamma v_\gamma^\pi(s^{\prime})\right)\right\} \\ &=G_{\mathrm{avg}}(\pi)+\gamma\sum d_\pi(s)\sum_{a\in\mathcal{A}}\pi(a|s)\sum_{s^{\prime}\in\mathcal{S}}\mathcal{P}_{ss^{\prime}}^av_\gamma^\pi(s^{\prime}) \\ &=G_{\mathrm{avg}}(\pi)+\gamma\sum_{s^{\prime}\in\mathcal{S}}d_\pi(s^{\prime})v_\gamma^\pi(s^{\prime}) \\ &=G_{\mathrm{avg}}(\pi)+\gamma J_{\mathrm{avg}}(\pi) \\ &=G_{\mathrm{avg}}(\pi)+\gamma G_{\mathrm{avg}}(\pi)+\gamma^2J_{\mathrm{avg}}(\pi) \\ &=G_{\mathrm{avg}}(\pi)+\gamma G_{\mathrm{avg}}(\pi)+\gamma^2G_{\mathrm{avg}}(\pi)+\cdots \\ &=\frac1{1-\gamma}G_{\mathrm{avg}}(\pi). \end{aligned}
Javg(π)=∑dπ(s)vγπ(s)=∑dπ(s)a∈A∑π(a∣s){s′∈S∑Pss′a(rss′a+γvγπ(s′))}=Gavg(π)+γ∑dπ(s)a∈A∑π(a∣s)s′∈S∑Pss′avγπ(s′)=Gavg(π)+γs′∈S∑dπ(s′)vγπ(s′)=Gavg(π)+γJavg(π)=Gavg(π)+γGavg(π)+γ2Javg(π)=Gavg(π)+γGavg(π)+γ2Gavg(π)+⋯=1−γ1Gavg(π). 其中,从右边第一式到第二式就是利用值函数内部递归关系来代入。详见本系列博客第二单元第2.1.4.3小节。从第二式怎么化简到第三式的呢?其实就是要证明:
G
a
v
g
(
π
)
=
∑
d
π
(
s
)
∑
a
∈
A
π
(
a
∣
s
)
∑
s
′
∈
S
P
s
s
′
a
r
s
s
′
a
G_{\mathrm{avg}}(\pi)=\sum d_\pi(s)\sum_{a\in\mathcal{A}}\pi(a|s)\sum_{s^{\prime}\in\mathcal{S}}\mathcal{P}_{ss^{\prime}}^ar_{ss^{\prime}}^a
Gavg(π)=∑dπ(s)a∈A∑π(a∣s)s′∈S∑Pss′arss′a 那么上式右边就是对于reward的期望。而左边的
G
a
v
g
(
π
)
G_{\mathrm{avg}}(\pi)
Gavg(π)是对于一个长度趋近于无限的序列的reward的平均值再求期望。而当序列的长度趋向于无穷时,就也是对于reward的期望。得证。下面再来看看怎么从第三式到第四式。要证三到四式,就是要证明:
∑
d
π
(
s
)
∑
a
∈
A
π
(
a
∣
s
)
∑
s
′
∈
S
P
s
s
′
a
v
γ
π
(
s
′
)
=
∑
s
′
∈
S
d
π
(
s
′
)
v
γ
π
(
s
′
)
\sum d_\pi(s)\sum_{a\in\mathcal{A}}\pi(a|s)\sum_{s^{\prime}\in\mathcal{S}}\mathcal{P}_{ss^{\prime}}^av_\gamma^\pi(s^{\prime})=\sum_{s^{\prime}\in\mathcal{S}}d_\pi(s^{\prime})v_\gamma^\pi(s^{\prime})
∑dπ(s)a∈A∑π(a∣s)s′∈S∑Pss′avγπ(s′)=s′∈S∑dπ(s′)vγπ(s′) 左边的式子可通过调整求和顺序来等价变形为:
∑
s
′
∈
S
(
∑
d
π
(
s
)
∑
a
∈
A
π
(
a
∣
s
)
)
P
s
s
′
a
)
v
γ
π
(
s
′
)
\sum_{s^{\prime}\in\mathcal{S}}\Big(\sum d_\pi(s)\sum_{a\in\mathcal{A}}\pi(a|s))\mathcal{P}_{ss^{\prime}}^a\Big)v_\gamma^\pi(s^{\prime})
s′∈S∑(∑dπ(s)a∈A∑π(a∣s))Pss′a)vγπ(s′) 则此时括号里的那一大堆就是
d
π
(
s
′
)
d_\pi(s^{\prime})
dπ(s′),因为我们先是枚举了所有状态s得分布再通过
P
s
s
′
a
\mathcal{P}_{ss^{\prime}}^a
Pss′a得到了
s
′
s^{\prime}
s′的分布。第四式到第五式就是利用
J
a
v
g
(
π
)
J_{\mathrm{avg}}(\pi)
Javg(π)的定义。接下来就是递归的套用上述过程并在最后使用等比数列的求和公式即可。得证。
有了上述定理,我们就证明了优化
J
a
v
g
(
π
)
J_{\mathrm{avg}}(\pi)
Javg(π)与优化
G
a
v
g
(
π
)
G_{\mathrm{avg}}(\pi)
Gavg(π)是等价的:
max
π
J
a
v
g
(
π
)
⇔
max
π
G
a
v
g
(
π
)
\max_\pi J_{\mathrm{avg}}(\pi)\Leftrightarrow\max_\pi G_{\mathrm{avg}}(\pi)
πmaxJavg(π)⇔πmaxGavg(π) 也即优化以上两种average cost函数得到的最优策略是相同的:
π
a
v
g
∗
=
arg
max
π
J
a
v
g
(
π
)
=
arg
max
π
G
a
v
g
(
π
)
\pi_{\mathrm{avg}}^*=\arg\max_\pi J_{\mathrm{avg}}(\pi)=\arg\max_\pi G_{\mathrm{avg}}(\pi)
πavg∗=argπmaxJavg(π)=argπmaxGavg(π) 整理一下上面的式子,可得以下结论:
J
a
v
g
(
π
a
v
g
∗
)
=
1
1
−
γ
G
a
v
g
(
π
a
v
g
∗
)
=
∑
d
a
v
g
∗
(
s
)
v
γ
π
a
v
g
∗
(
s
)
.
J_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)=\frac1{1-\gamma}G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)=\sum d_{\mathrm{avg}}^*(s)v_\gamma^{\pi_{\mathrm{avg}}^*}(s).
Javg(πavg∗)=1−γ1Gavg(πavg∗)=∑davg∗(s)vγπavg∗(s). 两种average cost的等价性也为我们在实际中根据它们各自的优化难度来选择合适的优化目标提供了依据。
5.1.2.2 折扣代价函数
折扣代价函数的定义如下:
J
γ
(
π
)
=
∑
d
i
n
i
t
(
s
)
v
γ
π
(
s
)
J_\gamma(\pi)=\sum d_{\mathrm{init}}(s)v_\gamma^\pi(s)
Jγ(π)=∑dinit(s)vγπ(s) 这里
d
i
n
i
t
(
s
)
d_{\mathrm{init}}(s)
dinit(s)表示初始状态的分布,可以取任意一个分布,没有影响。为了后面方便与average cost函数进行比较,我们可以使用average cost下的最优策略对应的分布来作为初始分布,即令
d
i
n
i
t
(
s
)
=
d
a
v
g
∗
(
s
)
d_{\mathrm{init}}(s)=d_{\mathrm{avg}}^*(s)
dinit(s)=davg∗(s)。那么,discounted cost函数可以具体化的写为:
J
γ
(
π
)
=
∑
d
a
v
g
∗
(
s
)
v
γ
π
(
s
)
J_\gamma(\pi)=\sum d_{\mathrm{avg}}^*(s)v_\gamma^\pi(s)
Jγ(π)=∑davg∗(s)vγπ(s) 定义
π
γ
∗
\pi_\gamma^*
πγ∗为discounted cost下的最优策略:
π
γ
∗
=
arg
max
π
J
γ
(
π
)
\pi_\gamma^*=\arg\max_\pi J_\gamma(\pi)
πγ∗=argπmaxJγ(π)
下面来考虑一个有趣的问题:在各自的cost函数下,代入各自的最优策略
π
a
v
g
∗
\pi_{\mathrm{avg}}^*
πavg∗和
π
γ
∗
\pi_\gamma^*
πγ∗,哪个cost值更优(更大,因为我们优化的目标是max)呢?给出下述结论:
定理2:
J
a
v
g
(
π
a
v
g
∗
)
≤
J
γ
(
π
γ
∗
)
\text{定理2:}J_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)\leq J_\gamma(\pi_\gamma^*)
定理2:Javg(πavg∗)≤Jγ(πγ∗) 给出这个定理的证明:
首先,观察discounted cost函数的定义式,因为其对于任何的初始分布得到的最优策略都是一样的,因此就等价于把每个状态s对应的值函数
v
γ
π
γ
∗
(
s
)
v_\gamma^{\pi_\gamma^*}(s)
vγπγ∗(s)都优化到最大值。即:
v
γ
π
(
s
)
≤
v
γ
π
γ
∗
(
s
)
,
for all
s
∈
S
.
v_\gamma^{\pi}(s)\leq v_\gamma^{\pi_\gamma^*}(s),\text{for all }s\in\mathcal{S}.
vγπ(s)≤vγπγ∗(s),for all s∈S. 那么显然有:
v
γ
π
a
v
g
∗
(
s
)
≤
v
γ
π
γ
∗
(
s
)
,
for all
s
∈
S
.
v_\gamma^{\pi_{\mathrm{avg}}^*}(s)\leq v_\gamma^{\pi_\gamma^*}(s),\text{for all }s\in\mathcal{S}.
vγπavg∗(s)≤vγπγ∗(s),for all s∈S. 那么就有:
J
a
v
g
(
π
a
v
g
∗
)
=
∑
d
π
a
v
g
∗
(
s
)
v
γ
π
a
v
g
∗
(
s
)
≤
∑
d
π
a
v
g
∗
(
s
)
v
γ
π
γ
∗
(
s
)
=
J
γ
(
π
γ
∗
)
.
\begin{aligned}J_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)&=\sum d_{\pi_{\mathrm{avg}}^*}(s)v_\gamma^{\pi_{\mathrm{avg}}^*}(s)\\&\leq\sum d_{\pi_{\mathrm{avg}}^*}(s)v_\gamma^{\pi_\gamma^*}(s)\\ & =J_\gamma(\pi_\gamma^*). \end{aligned}
Javg(πavg∗)=∑dπavg∗(s)vγπavg∗(s)≤∑dπavg∗(s)vγπγ∗(s)=Jγ(πγ∗).
下面来看看对于两种cost函数及其对应的两种最优策略,共有四种组合的情况:
G
a
v
g
(
π
a
v
g
∗
)
G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)
Gavg(πavg∗)、
G
a
v
g
(
π
γ
∗
)
G_{\mathrm{avg}}(\pi_\gamma^*)
Gavg(πγ∗)、
J
a
v
g
(
π
a
v
g
∗
)
J_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)
Javg(πavg∗)、
J
γ
(
π
γ
∗
)
J_\gamma(\pi_\gamma^*)
Jγ(πγ∗)。它们之间的大小关系如下:
1
1
−
γ
G
a
v
g
(
π
γ
∗
)
≤
1
1
−
γ
G
a
v
g
(
π
a
v
g
∗
)
=
J
γ
(
π
a
v
g
∗
)
≤
J
γ
(
π
γ
∗
)
.
\frac1{1-\gamma}G_{\mathrm{avg}}(\pi_\gamma^*)\leq\frac1{1-\gamma}G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)=J_\gamma(\pi_{\mathrm{avg}}^*)\leq J_\gamma(\pi_\gamma^*).
1−γ1Gavg(πγ∗)≤1−γ1Gavg(πavg∗)=Jγ(πavg∗)≤Jγ(πγ∗). 证明如下。首先,由
π
a
v
g
∗
\pi_{\mathrm{avg}}^*
πavg∗和
π
γ
∗
\pi_\gamma^*
πγ∗分别是各自cost函数下的最优策略,有:
G
a
v
g
(
π
a
v
g
∗
)
≥
G
a
v
g
(
π
γ
∗
)
,
J
γ
(
π
a
v
g
∗
)
≤
J
γ
(
π
γ
∗
)
.
G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)\geq G_{\mathrm{avg}}(\pi_\gamma^*),\\J_{\gamma}(\pi_{\mathrm{avg}}^{*})\leq J_{\gamma}(\pi_{\gamma}^{*}).
Gavg(πavg∗)≥Gavg(πγ∗),Jγ(πavg∗)≤Jγ(πγ∗). 由这两式的的第一式可证得
1
1
−
γ
G
a
v
g
(
π
γ
∗
)
≤
1
1
−
γ
G
a
v
g
(
π
a
v
g
∗
)
\frac1{1-\gamma}G_{\mathrm{avg}}(\pi_\gamma^*)\leq\frac1{1-\gamma}G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)
1−γ1Gavg(πγ∗)≤1−γ1Gavg(πavg∗)。而
1
1
−
γ
G
a
v
g
(
π
a
v
g
∗
)
=
J
a
v
g
(
π
a
v
g
∗
)
=
∑
d
a
v
g
∗
(
s
)
v
γ
π
a
v
g
∗
(
s
)
=
J
γ
(
π
a
v
g
∗
)
\frac1{1-\gamma}G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)=J_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)=\sum d_{\mathrm{avg}}^*(s)v_\gamma^{\pi_{\mathrm{avg}}^*}(s)=J_\gamma(\pi_{\mathrm{avg}}^*)
1−γ1Gavg(πavg∗)=Javg(πavg∗)=∑davg∗(s)vγπavg∗(s)=Jγ(πavg∗) 最后,有这两式的第二式可得
J
γ
(
π
a
v
g
∗
)
≤
J
γ
(
π
γ
∗
)
J_\gamma(\pi_{\mathrm{avg}}^*)\leq J_\gamma(\pi_\gamma^*)
Jγ(πavg∗)≤Jγ(πγ∗)。综上,得证。
在实际中,average cost更受青睐。但是,其分布和值函数的耦合问题导致不易求解。故而常常使用discounted cost函数代替。
5.1.2.3 使用Stochastic Linear Quadratic控制器来说明两种cost函数之间的关系
这部分与最优控制理论里面的LQR控制器有关。有关这部分的基础知识可以参考如下链接:LQR控制器。那么我们就follow5.1.1那里对于随机环境的线性状态空间模型的描述:
s
t
+
1
=
A
s
t
+
B
a
t
+
ξ
t
s_{t+1}=As_t+Ba_t+\xi_t
st+1=Ast+Bat+ξt 其中的
ξ
t
\xi_t
ξt是一个服从高斯分布的随机噪声,即
ξ
t
∼
N
(
0
,
σ
2
I
)
\xi_t\sim\mathcal{N}(0,\sigma^2I)
ξt∼N(0,σ2I)。那么,奖励信号也可写成一个关于当前状态
s
t
s_t
st和动作
a
t
a_t
at的二次型函数:
r
(
s
t
,
a
t
)
=
s
t
T
Q
s
t
+
a
t
T
R
a
t
r(s_t,a_t)=s_t^\mathrm{T}Qs_t+a_t^\mathrm{T}Ra_t
r(st,at)=stTQst+atTRat 这里的Q和R分别是半负定和负定的矩阵。这里与一般的LQR里面的对于Q、R的规定可能恰好相反,这里其实就是考虑一个最大化还是最小化的问题。一般的LQR那里需要最小化cost function,而我们这里需要最大化奖励信号。
- Average cost函数
结合针对于average cost函数的differential Bellman方程以及上述对于环境的状态空间建模,可得下述针对于LQR控制器的differential Bellman方程:
h ∗ ( s ) + ρ ∗ = max a { s T Q s + a T R a + E ξ { h ∗ ( A s + B a + ξ ) } } , h^*(s)+\rho^*=\max_a\left\{s^\mathrm{T}Qs+a^\mathrm{T}Ra+\mathbb{E}_\xi\{h^*(As+Ba+\xi)\}\right\}, h∗(s)+ρ∗=amax{sTQs+aTRa+Eξ{h∗(As+Ba+ξ)}}, 关于 h ∗ ( s ) h^*(s) h∗(s)和 ρ ∗ \rho^* ρ∗请参看第5.3.2.2节那里的叙述。可以证明,对于上述的二次型形式的代价函数(即reward),h函数也是一个关于状态s的二次型函数:
h ( s ) = s T P A v g s h(s)=s^\mathrm{T}P_\mathrm{Avg}s h(s)=sTPAvgs 其中, P A v g P_\mathrm{Avg} PAvg可由下式确定:
P A v g = A T P A v g A − A T P A v g B ( B T P A v g B + R ) − 1 B T P A v g A + Q . P_{\mathrm{Avg}}=A^\mathrm{T}P_{\mathrm{Avg}}A-A^\mathrm{T}P_{\mathrm{Avg}}B\left(B^\mathrm{T}P_{\mathrm{Avg}}B+R\right)^{-1}B^\mathrm{T}P_{\mathrm{Avg}}A+Q. PAvg=ATPAvgA−ATPAvgB(BTPAvgB+R)−1BTPAvgA+Q. 那么最优的average cost如下:
G a v g ( π a v g ∗ ) = σ 2 t r ( P A v g ) , G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)=\sigma^2\mathrm{tr}(P_{\mathrm{Avg}}), Gavg(πavg∗)=σ2tr(PAvg), 最佳的策略(控制率)如下:
π a v g ∗ ( s ) = − K A v g s , w h e r e K A v g = ( B T P A v g B + R ) − 1 B T P A v g A . \pi_{\mathrm{avg}}^{*}(s)=-K_{\mathrm{Avg}}s,\\\mathrm{where}\quad K_{\mathrm{Avg}}=\left(B^{\mathrm{T}}P_{\mathrm{Avg}}B+R\right)^{-1}B^{\mathrm{T}}P_{\mathrm{Avg}}A. πavg∗(s)=−KAvgs,whereKAvg=(BTPAvgB+R)−1BTPAvgA. - Discounted cost函数
写出此时的Bellman方程:
v γ ∗ ( s ) = max a E ξ { s T Q s + a T R a + γ v γ ∗ ( A s + B a + ξ ) } , v_\gamma^*(s)=\max_a\mathbb{E}_\xi\big\{s^\mathrm{T}Qs+a^\mathrm{T}Ra+\gamma v_\gamma^*(As+Ba+\xi)\big\}, vγ∗(s)=amaxEξ{sTQs+aTRa+γvγ∗(As+Ba+ξ)}, 最优的值函数也可写成一个关于状态s的二次型函数:
v γ ∗ ( s ) = s T P γ s + M , P γ = γ A T P γ A − γ 2 A T P γ B ( γ B T P γ B + R ) − 1 B T P γ A + Q , M = γ 1 − γ σ 2 t r ( P γ ) . \begin{gathered}v_\gamma^*(s)=s^\mathrm{T}P_\gamma s+M,\\ P_{\gamma}=\gamma A^{\mathrm{T}}P_{\gamma}A-\gamma^{2}A^{\mathrm{T}}P_{\gamma}B\big(\gamma B^{\mathrm{T}}P_{\gamma}B+R\big)^{-1}B^{\mathrm{T}}P_{\gamma}A+Q,\\M=\frac{\gamma}{1-\gamma}\sigma^{2}\mathrm{tr}\big(P_{\gamma}\big). \end{gathered} vγ∗(s)=sTPγs+M,Pγ=γATPγA−γ2ATPγB(γBTPγB+R)−1BTPγA+Q,M=1−γγσ2tr(Pγ). 最优的策略仍然是一个线性的反馈控制律:
π γ ∗ ( s ) = − K γ s , where K γ = γ ( γ B T P γ B + R ) − 1 B T P γ A . \pi_{\gamma}^{*}(s)=-K_{\gamma}s,\\\text{where}\quad K_{\gamma}=\gamma\left(\gamma B^{\mathrm{T}}P_{\gamma}B+R\right)^{-1}B^{\mathrm{T}}P_{\gamma}A. πγ∗(s)=−Kγs,whereKγ=γ(γBTPγB+R)−1BTPγA. 还需要说明的是,Linear Quadratic控制问题的一个有意思的性质是,最优的控制律由Ricaatti方程(Riccati equation)给出。尽管针对于不同的cost function,Ricatti方程的解可能不同,但是彼此之间都有联系。比如,我们上述的分别使用average cost和discounted cost函数得到的最优cost函数的值有如下关系:
lim γ → 1 1 ∑ i = 0 ∞ γ i v γ ∗ ( s ) = lim γ → 1 s T P γ s + γ 1 − γ σ 2 t r ( P γ ) ∑ i = 0 ∞ γ i = lim γ → 1 s T P γ s + γ 1 − γ σ 2 t r ( P γ ) 1 1 − γ = lim γ → 1 { ( 1 − γ ) s T P γ s + γ σ 2 t r ( P γ ) } = σ 2 t r ( P γ ) = G a v g ( π a v g ∗ ) . \begin{aligned} \operatorname*{lim}_{\gamma\to1}\frac{1}{\sum_{i=0}^{\infty}\gamma^{i}}v_{\gamma}^{*}(s)& =\lim_{\gamma\to1}\frac{s^\mathrm{T}P_\gamma s+\frac\gamma{1-\gamma}\sigma^2\mathrm{tr}(P_\gamma)}{\sum_{i=0}^\infty\gamma^i} \\ & =\lim_{\gamma\to1}\frac{s^\mathrm{T}P_\gamma s+\frac\gamma{1-\gamma}\sigma^2\mathrm{tr}(P_\gamma)}{\frac{1}{1-\gamma}} \\ &=\lim_{\gamma\to1}\{(1-\gamma)s^\mathrm{T}P_\gamma s+\gamma\sigma^2\mathrm{tr}(P_\gamma)\} \\ &=\sigma^2\mathrm{tr}(P_\gamma) \\ &=G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*). \end{aligned} γ→1lim∑i=0∞γi1vγ∗(s)=γ→1lim∑i=0∞γisTPγs+1−γγσ2tr(Pγ)=γ→1lim1−γ1sTPγs+1−γγσ2tr(Pγ)=γ→1lim{(1−γ)sTPγs+γσ2tr(Pγ)}=σ2tr(Pγ)=Gavg(πavg∗). 这里比较有趣的是 v γ ∗ ( s ) v_{\gamma}^{*}(s) vγ∗(s)是一个与具体状态有关的量,而 G a v g ( π a v g ∗ ) G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*) Gavg(πavg∗)与具体状态无关。但是在 γ → 1 \gamma\to1 γ→1的前提下,两者是相等的。这也再次告诉我们,尽管在实际建模中average cost函数更受青睐,但是其不易计算。而当我们采用 γ \gamma γ较大的discounted cost函数时,是可以得到较好的近似的。
下面我们用一个实际的例子来比较一下上述两种cost函数下的Linear Quadratic控制器。例子的各个参数如下表所示:
参数 | 值 |
---|---|
状态维数n | 1 |
A | 2 |
B | 1 |
ξ t \xi_t ξt | N ( μ = 0 , σ 2 = 1 ) \mathcal{N}(\mu=0,\sigma^2=1) N(μ=0,σ2=1) |
Q | -1 |
R | -2 |
γ \gamma γ | 0.7 |
控制律均采用 π ( s ) = − K s \pi(s)=-Ks π(s)=−Ks的形式。初始状态均置为0。我们按以下几个方面来比较:
- 达到稳态时的状态分布
将反馈控制律代入,可以将状态空间模型化简为下面的形式:
s t + 1 = ( A − B K ) s t + ξ t , ξ t ∼ N ( 0 , σ 2 ) s_{t+1}=(A-BK)s_t+\xi_t,\quad \xi_t{\sim}\mathcal{N}(0,\sigma^2) st+1=(A−BK)st+ξt,ξt∼N(0,σ2)
当达到稳态之后,对上式左右两侧取方差,并注意到此时 s t + 1 = s t = s s_{t+1}=s_t=s st+1=st=s,因此移项整理之后可得状态s的分布:
s ∼ N ( 0 , σ 2 1 − ( A − B K ) 2 ) . s\sim\mathcal{N}\left(0,\frac{\sigma^2}{1-(A-BK)^2}\right). s∼N(0,1−(A−BK)2σ2). 带入相应的参数值,我们可以得到两种cost函数下的稳态状态分布如下:
d π a v g ∗ ( s ) = 1 − ( A − B K A v g ) 2 1 2 π σ exp ( − s 2 2 σ 2 ( 1 − ( A − B K A v g ) 2 ) ) , d π γ ∗ ( s ) = 1 − ( A − B K γ ) 2 1 2 π σ exp ( − s 2 2 σ 2 ( 1 − ( A − B K γ ) 2 ) ) . \begin{gathered} d_{\pi_{\mathrm{avg}}^{*}}(s)=\sqrt{1-\left(A-BK_{\mathrm{Avg}}\right)^{2}}\frac{1}{\sqrt{2\pi} \sigma}\operatorname{exp}\left(-\frac{s^{2}}{2\sigma^{2}}\Big(1-\Big(A-BK_{\mathrm{Avg}}\Big)^{2}\Big)\right), \\ d_{\pi_{\gamma}^{*}}(s)=\sqrt{1-\left(A-BK_{\gamma}\right)^{2}}\frac{1}{\sqrt{2\pi} \sigma}\exp\biggl(-\frac{s^{2}}{2\sigma^{2}}\biggl(1-\bigl(A-BK_{\gamma}\bigr)^{2}\biggr)\biggr). \end{gathered} dπavg∗(s)=1−(A−BKAvg)22πσ1exp(−2σ2s2(1−(A−BKAvg)2)),dπγ∗(s)=1−(A−BKγ)22πσ1exp(−2σ2s2(1−(A−BKγ)2)). - 最优cost函数的值
我们在第5.1.2.2节中讨论了四种cost函数的值之间的关系。根据那里的公式以及Linear Quadratic控制器的具体情况,可得:
J γ ( π γ ∗ ) = ∑ d π a v g ∗ ( s ) v γ ∗ ( s ) = σ 2 P γ ( γ 1 − γ + 1 1 − e 2 ) , G a v g ( π γ ∗ ) = ( 1 − γ ) ∑ d π γ ∗ ( s ) v γ ∗ ( s ) = σ 2 P γ ( γ + 1 − γ 1 − e 2 ) , J γ ( π a v g ∗ ) = σ 2 P A v g 1 1 − γ , G a v g ( π a v g ∗ ) = σ 2 P A v g . \begin{gathered} J_\gamma\big(\pi_\gamma^*\big)=\sum d_{\pi_{\mathrm{avg}}^*}(s)v_\gamma^*(s)=\sigma^2P_\gamma\left(\frac{\gamma}{1-\gamma}+\frac{1}{1-e^2}\right), \\ G_{\mathrm{avg}}(\pi_{\gamma}^{*})=(1-\gamma)\sum d_{\pi_{\mathrm{\gamma}}^{*}}(s)v_{\gamma}^{*}(s)=\sigma^{2}P_{\gamma}\left(\gamma+\frac{1-\gamma}{1-e^{2}}\right), \\ J_{\gamma}\big(\pi_{\mathrm{avg}}^{*}\big)=\sigma^{2}P_{\mathrm{Avg}}\frac{1}{1-\gamma} , \\ G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)=\sigma^2P_{\mathrm{Avg}} . \end{gathered} Jγ(πγ∗)=∑dπavg∗(s)vγ∗(s)=σ2Pγ(1−γγ+1−e21),Gavg(πγ∗)=(1−γ)∑dπγ∗(s)vγ∗(s)=σ2Pγ(γ+1−e21−γ),Jγ(πavg∗)=σ2PAvg1−γ1,Gavg(πavg∗)=σ2PAvg. 其中,前两个式子直接带入公式即可。而第四式原本是 G a v g ( π a v g ∗ ) = σ 2 t r ( P A v g ) G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)=\sigma^2\mathrm{tr}(P_{\mathrm{Avg}}) Gavg(πavg∗)=σ2tr(PAvg),这里因为状态维数为1,所以直接就是 P A v g P_{\mathrm{Avg}} PAvg。第三式是因为 J γ ( π a v g ∗ ) = 1 1 − γ G a v g ( π a v g ∗ ) J_{\gamma}(\pi_{\mathrm{avg}}^{*})=\frac{1}{1-\gamma}G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^{*}) Jγ(πavg∗)=1−γ1Gavg(πavg∗)由第四式得到的。那么,经过仿真可得到四个值分别为:
Performance measure Real Simulation Error J γ ( π γ ∗ ) − 14.64 − 14.38 1.9 % J γ ( π a v g ∗ ) − 16.17 − 16.00 1.1 % ( 1 − γ ) − 1 G a v g ( π a v g ∗ ) − 16.17 − 16.41 1.5 % ( 1 − γ ) − 1 G a v g ( π γ ∗ ) − 21.28 − 21.22 0.2 % \begin{array}{|c|c|c|c|}\hline\text{Performance measure}&\text{Real}&\text{Simulation}&\text{Error}\\\hline J_\gamma(\pi_\gamma^*)&-14.64&-14.38&1.9\%\\\hline J_\gamma(\pi_{\mathrm{avg}}^*)&-16.17&-16.00&1.1\%\\\hline(1-\gamma)^{-1}G_{\mathrm{avg}}(\pi_{\mathrm{avg}}^*)&-16.17&-16.41&1.5\%\\\hline(1-\gamma)^{-1}G_{\mathrm{avg}}(\pi_\gamma^*)&-21.28&-21.22&0.2\%\\\hline\end{array} Performance measureJγ(πγ∗)Jγ(πavg∗)(1−γ)−1Gavg(πavg∗)(1−γ)−1Gavg(πγ∗)Real−14.64−16.17−16.17−21.28Simulation−14.38−16.00−16.41−21.22Error1.9%1.1%1.5%0.2% 与我们在第5.1.2.2节中的结论一致。
5.1.3 策略迭代vs.值迭代
就像之前第二单元博客介绍过的,策略迭代是通过PEV和PIM两个步骤交替的优化最终取得最优解的。而值迭代则是首先通过fixed-point iteration求解Bellman方程,然后再通过贪心策略来计算最优策略。两种方法的收敛性分别由以下原因保证:策略迭代通过PEV和PIM两个步骤的交替优化,每一轮经过这两步后都可以保证更新后的策略 π ′ \pi' π′优于原策略 π \pi π。而值迭代则是通过Bellman算子的收缩性(contractive)来保证收敛性。
5.1.3.1 策略迭代
- Policy Evaluation(PEV):PEV适用于model-based和model-free的RL。尽管各种PEV表面上看起来可能有很大差别,但是其核心是一样的,就是试图为当前的策略找到一个准确的估计。
- Policy Improvement(PIM):PIM需要使得更新后的策略比更新前的策略更优。这可由Policy Improvement Theorem保证。详见第二单元博客。
5.1.3.2 值迭代
值迭代把Bellman方程看作一个contractive的operator,这个operator定义了一个在Banach空间里从一个点到另一个点的非线性映射。之后值迭代采用fixed-point iteration来求最优解。