使用凸优化方法确定随机过程的不安全概率和风险轮廓

《Automatica》：Unsafe probabilities and risk contours for stochastic processes using convex optimization

　　贾里德·米勒（Jared Miller）|马泰奥·塔基（Matteo Tacchi）|迪迪埃·亨里翁（Didier Henrion）|马里奥·斯奈尔（Mario Sznaier）
德国斯图加特大学数学系数学系统理论系主任

**摘要**
在评估动态系统轨迹的安全规范时，能够限定最坏情况下不安全（违反约束）的概率至关重要。随机安全性和最坏情况下不安全概率的认证可以表示为无限维线性程序（例如随机屏障函数、占用度量问题）。本文证明了在随机动力学的紧凑性和规则性条件下，这些无限维线性程序及其有限维矩和平方截断（Moment-Sum-of-Squares truncations）对真实的不安全概率是不一致的（非保守的）。例如，可以为随机过程生成不安全概率估计和风险轮廓。

**1. 引言**
如果动态系统的轨迹在从初始集合X0开始执行的过程中永远不会接触到危险集合Xu，则该轨迹被认为是安全的。确定性系统的安全规范实例可能包括飞机不会坠毁到地面、自动驾驶车辆不会撞到行人，或者房间内的温度不会违反舒适度约束。在随机环境中，这样的安全声明可以表述为几乎必然成立的安全规范（Clark, 2021）。然而，几乎必然的安全保证可能过于保守，或者在某些情况下无法实现（例如非消失的高斯噪声等遍历过程）。因此，可以根据（1?p）安全性来认证安全性，其中以最坏概率p从X0接触到危险集合Xu。之前关于使用（控制）屏障函数认证（1?p）安全性的工作包括（Prajna等人，2004年；Salamati等人，2021年；Santoyo等人，2021年）。这些方法通过利用超鞅概率不等式（Kushner, 1967年中的引理1）提供了安全的充分条件，但它们没有提供存在性、非保守性或收敛性的条件。其他随机验证方法包括可达性（Digailova和Kurzhanski，2004年；Prandini等人，2006年）、避达分析（Bansal等人，2017年；Mohajerin Esfahani等人，2016年）、毁坏（Asmussen & Albrecher，2010年）以及组合抽象（Lavaei等人，2020年；Nejati等人，2022年）。

本文的重点是精确计算从X0开始达到危险集合Xu的最坏情况下不安全概率，并提供在非保守条件下可以使用凸优化方法计算这一最坏情况下不安全概率的条件。最坏情况下的不安全问题通常是一个非凸但有限维的优化问题。利用基于占用度量的最优控制提升方法（Lewis和Vinter，1980年；Rubio，1975年），非凸不安全问题将被提升为无限维线性程序（LP）。本研究的主要贡献是提出LP公式和条件，使得：
(1) 凸LP目标等于最坏情况下的不安全概率；
(2) 从LP导出的风险轮廓（具有初始时间和初始条件参数的连续函数）将从上方以L1意义收敛到真实的最坏情况下概率。

为了获得数值解，本文开发的无限维LP必须截断为有限维程序。用于通用随机过程的截断算法包括网格化（Cho & Stockbridge，2002年）、径向基函数选择（Kariotoglou等人，2017年）、随机采样（Huynh等人，2012年）以及带有反例引导验证的神经网络参数化（Abate等人，2021年；Abate等人，2018年）。在随机过程和所有问题数据都有有理表示的特定情况下，可以使用矩和平方（SOS）层次结构（Lasserre，2009年）将LP截断为按多项式程度参数化的半定规划（SDPs）层次结构。矩和平方层次结构在随机分析和控制中的应用包括区域验证（Steinhardt & Tedrake，2012年）、避达估计（Xue等人，2023年）、退出时间估计（Henrion等人，2021年）以及（条件）风险价值上界（Miller等人，2025年）。

与本文最相关的两项先前工作是Sloth和Wisniewski（2015年）以及Jasour和Williams（2019年）。Sloth和Wisniewski（2015年）的工作提出了无限维LP，用于生成固定初始分布μ0和概率水平p的随机微分方程（SDE）的概率（1?p）安全集的超水平外逼近。我们的问题评估了相对于自由初始分布的特定初始集合X0的安全性，进行了风险轮廓分析，并将框架扩展到更一般的随机过程。我们在提出这个LP时使用的论点与之前工作中引理2中的步骤类似（Sloth & Wisniewski，2015年）。Jasour和Williams（2019年）的工作为具有分布不确定性的优化问题合成了风险轮廓（重点在于障碍物和路径规划），而我们的论文通过分析随机动态系统的安全性扩展了这项工作。

本文的结构如下：第2节介绍符号、随机过程和非负性的SOS证明等预备知识。第3节提出无限维凸LP，用于计算相对于初始集合X0或初始分布μ0的随机过程的最坏情况下不安全概率。第4节使用矩和平方层次结构将不安全LP截断为收敛的有限维SDP序列。第5节通过具有多项式动态的示例随机过程演示了不安全概率和风险轮廓的计算。第6节总结了本文。

**2. 预备知识**
2.1. 符号
介于a和b之间的自然数子集表示为a?B。两个量的最小值表示为a∧b=min(a,b)。在Rn中具有实值系数的多项式集合表示为R[x]。每个多项式c∈R[x]可以通过对多索引α∈Nn的有限求和唯一描述，即c(x)=∑αcαxα。多项式的次数degc等于最大的指数和∑i=1nαi，使得cα≠0。x中次数最多为2k的多项式集合表示为R[x]≤2k。

2.2. 分析和测度理论
集合S上的连续函数集合表示为C(S)。其非负连续函数的子锥表示为C+(S)。给定一个乘积集S×H，在第一个变量s上一次连续可微且在第二个变量h上二次连续可微的函数集合表示为C1,2(S×H)。在S上支持的带符号的Borel测度集合表示为M(S)，其在S上支持的非负Borel测度的子集表示为M+(S)。测度μ∈M+(Rn)的支持集是最小的闭集s，使得每个开邻域N?(s)都有正测度μ(N?(s))>0。当S是紧凑的时，C+(S)和M+(S)在拓扑上是对偶的，并且它们通过Lebesgue积分定义了对偶配对〈?,?〉：?f∈C+(S),μ∈M+(S):〈f,μ〉=∫Sf(s)dμ(s)。这种对偶配对将扩展为C(S)和M+(S)之间的双线性配对。集合A?S的μ-测度也可以表示为μ与指示函数IA的配对：μ(A)=〈IA,μ〉=inf{〈w,μ〉∣w∈C(S),w≥IA}。μ的质量表示为μ(S)=〈1,μ〉，如果这个质量为1，则μ是一个概率分布。一个重要的概率分布例子是在s′处仅支持的Dirac delta δs′，遵循点评估配对规则?f∈C(S):〈f,δs′〉=f(s′)。当μ是一个概率分布时，事件x∈A的概率表示为Probμ[x∈A]=μ(A)。

2.3. 随机过程
随机过程是一组通过系统动态（Stroock & Varadhan，1997）相关联的时间索引随机变量{μt}（沿流图的前向传播）。随机过程的属性将根据适当的集合C（如C=C([t0,T]×X）中的测试函数进行分析。根据分布μt，在时间t处测试函数v(s,x)的期望表示为E[v(t,x)∣μt]。设Δt>0是一个时间步长，随机过程的生成器LΔt满足（对于所有适当的测试函数v∈domLΔt）LΔtv(t,x)=limΔt′→ΔtE[v(t+Δt′,x)∣μt+Δt′]?v(t,x)Δt′。我们将LΔt的定义域表示为C=domLΔt，使得C是LΔt下连续函数的原像的子集。离散时间随机过程的生成器是LΔt，其中Δt>0，它是相对于测试函数类C=C([t0,T]×X)定义的。对于遵循xt+Δt=f(t,xt,λt)的离散时间法则，其中时变参数λt∈Λ具有概率分布ξ(λt)，相关的生成器满足LΔtv(t,x)=1Δt∫λ∈Λvt+Δt,f(t,x,λ)?v(t,x)dξ(λ)。

2.4. 占用度量
随机过程的生成器可以使用占用度量理论来描述。给定一个初始条件（点）x0∈X0?X和一个时间t∈[t0,T]，我们用x(t∣x0)表示从初始分布δx=x0开始的随机过程{μt}的轨迹，其生成器为L。随机过程{μt}从开始时间t0到停止时间t?的占用度量μ∈M+([t0,T]×X)和终止度量μτ∈M+([t0,T]×X)满足?A?[tt0,T],B?X: μ(A×B)=∫X0∫t=t0t?IA×Bt,x(t∣x0)dtdμ0(x0)（4a）μτ(A×B)=∫X0IA×Bt?,x(t?∣x0)dμ0(x0）。（4）中的测度也可以相对于停止时间t?∈M+([t0,T])的概率分布来定义。一组度量（μ0,μτ,μ）∈M+(X0)×M+([t0,T]×X)2满足具有生成器L的随机过程的鞅关系（3）?v∈C:〈v,μτ〉=〈v(t0,?),μ0〉+〈Lv,μ〉。关系（5）被称为常微分方程（ODE）的Liouville方程，在随机微分方程（SDE）的背景下称为Dynkin方程。我们也可以用（6）μτ=δt0?μ0+L?μ的隐式（?v∈C）量化等价地表示（5）。满足（6）的度量元组（μ0,μτ,μ）是放松的占用度量。并非每个放松的占用度量都是真正的占用度量，Cho和Stockbridge（2002）的定理3.1概述了放松占用度量在随机轨迹图上得到支持的条件。

**3. 不安全问题和线性程序**
本节将提出LP来计算L生成的随机过程的最坏情况下不安全概率。

3.1. 问题陈述
本文的重点是通过计算轨迹达到不安全集合Xu的概率来进行随机过程的风险分析。轨迹在状态空间X?Rn中的最大时间范围T内演变，从初始集合X0?X开始，可能会进入不安全集合Xu?X。我们考虑由生成器L描述的随机过程（例如，SDE），其中μt是从X0开始的过程的时变状态概率分布。从X到达Xu的撞击时间分别为τX=inf{t:[x∈?X∣μt]}和τu=inf{t:[x∈Xu∣μt]}，它们共同构成了停止时间τm=min(τu,τX,T)。从X0开始的L轨迹在时间范围[t0,T]内的最坏情况下不安全概率为：
(7a)P?(t0,x0)=ProbμT[x∈Xu]
(7b)s.t.x(0)=x0
(7c)x(t)遵循L?t∈[t0,τm]。遵循停止时间τm的轨迹将在到达不安全集合Xu、集合?X的边界或经过T时间单位时停止（以先发生的情况为准）。停止时间τm还确保轨迹x(t)在所有相关t期间始终位于X内。从初始集合X0开始的L轨迹的最坏情况下不安全概率为：
(8a)P?(t0,X0)=supx0∈X0P?(t0,x0)，
而从初始条件概率分布μ0∈M+(X0)开始的平均最坏情况下不安全概率为：
(8b)P?(t0,μ0)=∫X0P?(t0,x0)dμ0(x0)。对于单个初始条件x0∈X，单一和不平均的不安全概率P(t0,x0)和P(t0,δx0)是相等的。如果P?(t0,X0)=0，则初始集合X0是安全的；如果P?(t0,X0)=1，则初始集合X0是不安全的。类似地，如果P?(t0,μ0)=0，则初始条件分布μ0是安全的；如果P?(t0,μ0)=1，则初始条件分布μ0是不安全的。任何其他值P?(t0,X0)∈(0,1)都表示随机过程的最大不安全概率。本文的目标是形成一个收敛的凸优化问题，其目标等于（8a）和（8b），并找到一个连续的风险轮廓v∈C，该轮廓将从上方以L1意义收敛到最优映射P?(t,x0)。因此，任何这样的风险函数v的水平集都可以用来上界限制从任何(t,x0)∈[t0,T]×X开始的L的随机执行风险。**假设**
在本文中，我们将提出以下假设，以便通过凸规划获得非保守的安全性分析条件：
A1：状态集X0、Xu、X都是紧凑的。
A2：时间区间[t0,T]是紧凑的。
A3：轨迹在首次离开X时停止（τX∧T）。
A4：测试函数集C=dom(L)满足C?C([t0,T]×X)，且1∈C以及L1=0。
A5：C能够区分不同点，并且是乘法封闭的。
A6：存在一个可数集{vk}∈C，使得对于所有v∈C，(v,Lv)包含在{(vk,Lvk)}的线性张成的有界点态闭包中。

假设A1–A2确保了空间的紧凑性。停止假设A3确保轨迹永远不会离开空间X并在之后返回。假设A4–A6与Cho和Stockbridge（2002）中的条件1相同。假设A4–A6非常通用，因为它们适用于马尔可夫过程、随机微分方程（SDE）和Lévy过程。

**3.3. 测度规划**
我们将提出一个凸的但无限维的线性规划（LP），以对程序（8）进行上界限制。该LP将基于以下变量定义：
(9a) μ0∈M+(X0)（初始分布）
(9b) μ∈M+([t0,T]×X)（占用分布）
(9c) μc∈M+([t0,T]×X)（互补分布）
(9d) μu∈M+([t0,T]×Xu)（不安全分布）
总和μu+μc将作为一个具有自由终止时间的终端测度。这些测度变量将被设计为随机过程的松弛占用测度。在解决X0问题时（8a），初始分布μ0是一个优化变量；而在解决平均不安全问题时（8b），它被固定为一个给定的初始分布。

我们论文的核心贡献是以下关于无松弛差距的声明：
**定理1**
在假设A3下，以下LP将对（8a）中的值进行上界限制，即p?(t0,X0)≥P?(t0,X0)：
(10a) p?(t0,X0)=sup〈1,μu〉
(10b) s.t. μu+μc=δt0?μ0+L?μ
(10c) 〈1,μ0〉=1
(10d) μ0∈M+(X0)
(10e) μ,μc∈M+([t0,T]×X)
(10f) μu∈M+([t0,T]×Xu)
在假设A1–A6下，目标（7c）是相等的（p?(t0,X0)=P?(t0,X0)）。

**证明**
见附录A。

**推论2**
在μ0预先指定的情况下，以下具有目标p?(t0,μ0)的LP将在假设A3下对（8b）中的P?(t0,μ0)进行上界限制，并且在假设A1–A6下没有松弛差距，即p?(t0,μ0)=P?(t0,μ0)：
(11a) p?(t0,μ0)=sup〈1,μu〉
(11b) s.t. μu+μc=δt0?μ0+L?μ
(11c) μ,μc∈M+([t0,T]×X)
(11d) μu∈M+([t0,T]×Xu)
证明**
上界和无松弛差距的结果是通过回顾附录A中用于证明定理1的步骤得出的。对于上界，设（μ[x0],μc[x0],μu[x0]）是由附录A.1中的构造过程给定初始点x0和停止时间t?得到的测度。因此，可以通过对所有x0∈μ0的平均值来收集（μ[x0],μc[x0],μu[x0]），在此过程中成本得到保持。对于无松弛差距，元组（μ0,μu+μc,μ）仍然支持在随机轨迹的图上，因此可以无需修改地使用定理1的证明步骤。

**3.4. 函数规划**
本小节将定义在连续函数变量v∈C上的LP，这些LP与之前的LP是对偶的。我们首先介绍这些对偶性声明。然后，我们将分析这个解v的属性（注意它将在逐点意义上对P?(t0,x0)进行上界限制），并概述其收敛性质。由于这种逐点满足性，这个解v也被称为“风险轮廓”。这些风险轮廓将属于一个受限的函数类：
(12) R=Lv(t,x)≤0 对所有(t,x)∈[t0,T]×X
v∈C ∣ v(t,x)≥0 对所有(t,x)∈[t0,T]×X
v(t,x)≥1 对所有(t,x)∈[t0,T]×Xu。

集合R取决于集合C、[t0,T]、X和Xu。

**定理3**
关于辅助函数v的LP形成了一个弱对偶（p?(t0,X0)≥d?(t0,X0)）：
(13a) d?(t0,X0)=infγ∈R,v∈R
(13b) s.t. γ≥v(t0,x) 对所有x∈X0。在假设A1–A4下，将成立强对偶（p?(t0,X0)=d?(t0,X0）。

**证明**
见附录B。

**注释4**
约束（13b）–（13b）的形式与固定概率水平γ下的随机（时间依赖）屏障函数相同（来自Prajna等人，2004年）。我们的工作涉及将γ作为优化变量使用，以及在假设A1–A6下证明无松弛差距。

**推论5**
以下LP与（11）弱对偶，并且在假设A1–A4下具有强对偶性：
(14a) d?(t0,μ0)=infv∈R∫X0v(t0,x)d(x0)
(14b) s.t. Lv(t,x)≤0 对所有(t,x)∈[t0,T]×X。

**证明**
这种强对偶性通过对定理3和附录B中使用的方法进行轻微修改而成立。

我们现在概述通过解决问题找到的风险轮廓v函数的属性。首先分析v的属性是其有界性和对P?的逐点上界限制。

**定理6**
在假设A1–A6下，每个v∈R都满足v(t0,x0)≥P?(t0,x0)对于所有初始条件x0∈X。

**证明**
见附录C。

**命题7**
d?(t0,X0)和d?(t0,μ0)都被1上界限制。

**证明**
函数v(t,x)=1是R的成员。点(v=1,γ=1)对于（10）是可行的，其值为γ=1，而(v=1)对于（11）是可行的，其值为〈1,μ0〉=1。因此，任何v∈R的级别集都可以用来对从任何(t,x0)∈[t0,T]×X开始的L的随机执行风险进行上界限制。根据以下定理，解决（14）得到的解v将在L1(μ0)意义上最接近真实的不安全概率映射：
**定理8**
设{vk}k≥1是一系列解，每个vk∈R，并且∫X0vk(t0,x0)dμ0(x0)→d?(t0,μ0)（在假设A1–A6下）。那么vk在L1(μ0)中收敛到P?，如（15）所示：
∫Xvk(t0,x0)?P?(t0,x0)dμ0(x0)?k→∞0。

**证明**
见附录D。

**注释9**
选择参考分布μ0为紧凑集X上的均匀分布可以为实践者提供关于P?结构的最可能信息。

**4. 不安全半定规划**
在中的无限维线性规划必须被截断为有限维优化问题，以便进行可行的数值计算。执行这种截断的一种方法是矩-SOS层次结构。本节将首先回顾用于截断具有多项式结构数据的无限维LP的矩-SOS层次结构，然后将该层次结构应用于的截断。

**4.1. 平方和方法**
设c∈R[x]≤2k是一个多项式。如果对于所有x∈Rn：c(x)≥0，则多项式c是非负的。如果存在N个多项式qj∈R[x]使得c(x)=∑j=1Nqj2(x)，则多项式c是SOS的。SOS多项式的锥体是Σ[x]，其度数≤2k的SOS多项式的子集是Σ[x]≤2k。只有在一元多项式（n=1）、二次多项式（2k=2）或二元四次多项式（n=2,2k=4）的情况下，SOS多项式集合才等于非负多项式集合（Hilbert, 1888）。验证多项式的非负性通常是一个NP难题，而且随着n和d的增加，SOS多项式在非负多项式集合中会消失（Blekherman, 2006）。

**4.2. 平方和假设**
为了将SOS方法应用于问题的收敛近似，我们需要额外的多项式结构假设：
A7：集合X0、Xu、X都是阿基米德基本半代数（BSA）集。
A8：生成器L在多项式v(t,x)∈R[t,x]?Lv(t,x)∈R[t,x]下是封闭的。

给定一个度数k和满足假设A8的生成器L，我们可以定义动态度数k?为k?=maxv∈R≤2k[t,x]deg(Lv(t,x))/2。

**4.3. 不安全概率SDP**
程序（13）的度数k的SOS紧化是（18a）dk?(t0,X0)=infγ∈Rγ（18b）s.t. γ?v(t0,x)∈Σ[X0]≤2k
（18c）?Lv(t,x)∈Σ[[0,T]×X]≤2k
（18d）v(t,x)∈Σ[[0,T]×X]≤2k
（18e）v(t,x)?1∈Σ[[0,T]×Xu]≤2k
程序（14）的度数k的SOS紧化是（19a）dk?(t0,μ0)=inf∫Xv(t0,x)dμ0(x）
（19b）s.t. ?Lv(t,x)∈Σ[[0,T]×X]≤2k
（19c）v(t,x)∈Σ[[0,T]×X]≤2k
（19d）v(t,x)?1∈Σ[[0,T]×Xu]≤2k

**注释10**
目标总是被dk?(t0,X0)≤1和dk?(t0,μ0)≤1上界限制的，因为对于它们的约束，v(t,x)=1总是可行的（来自命题7）。

**定理11**
在假设A1–A8下，程序（18）将收敛到（13），即limk→∞dk?(t0,X0)=d?(t0,X0），程序（19）将收敛到（14），即limk→∞dk?(t0,μ0)=d?(t0,μ0）。

**证明**
见附录E。

设vk(t,x)∈Σ[[0,T]×X]≤2k是程序（19）在度数k下的解，设Iu是在不安全集Xu上的0/1指示函数。那么从点x0∈X开始时的不安全概率被（20）上界限制：
q1:k(x)=min(1,mink′∈1..kvk′(t0,x))。

**推论12**
函数序列q1:k(x)在k递增的情况下将在测度μ0中收敛到v(t0,x)。

**证明**
这个推论来自定理8，其中序列{vk′}k′∈1…k用于近似v?。

**注释13**
设v1,v2是通过解决程序（19）在度数k1∈k2∈1..k（且k1≤k2）且目标dk1(t0,μ0)≥dk2(t0,μ0）得到的多项式。在（20）中沿所有度数取最小值可以比仅使用vk作为风险轮廓得到更精确的P?估计。

**4.4. 计算复杂性**
我们将通过它们SOS程序中涉及的最大Gram矩阵的大小来量化度数k紧化的计算复杂性。在典型的情况下，当k?>k（仅在L将每个多项式映射为常数时违反A1–A8）时，最大的Gram矩阵将出现在Lie约束（18c）、（19b）中。Lie约束每个都有n+1个变量（t,x），因此使用单项式基础时的Gram矩阵大小为（n+1+k?k?）。所有其他约束的度数都较低（k而不是k?），或者仅基于n个变量x提出。因此，使用内点方法解决SOS程序的复杂性将基于O((n+1)6k)（对于固定的k）或O(k(n+1)（对于固定的n）（Lasserre, 2009）。

**5. 数值示例**
MATLAB（2021a）代码可用于重现所有示例，可在https://github.com/jarmill/prob_unsafe获取。所有程序都使用YALMIP（Lofberg, 2004）建模，并使用MOsek 10（MOSEK ApS, 2023）求解。所有示例都将初始时间设为t0=0。实验在配备16.0 GB RAM的Intel Core i9 2.30 GHz上运行。

**5.1. 二维三次SDE**
我们的第一个演示分析了Prajna等人（2004）示例1中的三次多项式SDE的安全性：
(21) dx=x2?x1?x2?1
2x1
3dt+0
0
0.1dW。
在状态空间X=[?2,2]2内评估（21）的安全性，直到时间范围T=5。不安全集Xu是一个月牙形区域，定义为Xu={x∈R2∣0.90832≤(x1+0.5006)2+(x2+0.2902)2,0.52≥(x1?0.2)2+x22}。初始集X0是一个半径为R0、中心为[0.85;?0.75]的圆。图1绘制了从X0（洋红色区域，R0=0.2）开始的（21）的轨迹。注意一些采样轨迹触及红色半圆Xu的右角（因此是不安全的）。表1报告了程序（18）计算出的半径为R0的圆形初始集的不安全概率。图2(a)绘制了在T=5时通过解决（18）（阶数k=6，R0=0.2）获得的不安全上界函数min(1,v(0,x))。不安全集Xu是红色的月牙形区域。洋红色圆圈是初始集X0的边界。注意，在X0周围的区域，概率估计比远离X0的Prob=1饱和区域更精确。

**下载：** 下载高分辨率图像（192KB）
**下载：** 下载全尺寸图像

**图1.** SDE（21）的轨迹（青色）在圆盘X0（紫色边界）中初始化，击中不安全区域Xu（红色月牙形）。
（关于此图例中颜色的解释，请参阅本文的网页版本。）

**表1.** 系统（21）的不安全概率上界。

**顺序**
1 2 3 4 5 6
R0=0 10.944
2 0.607
2 0.480
5 0.406
8 0.369
6 R0=0.2 10.973
6 0.794
3 0.709
10 0.639
0 0.613

**相比之下，** 图2(b)解决了程序（14）的六阶SOS紧化，以生成在X=[?2,2]2中有效的风险图v。该风险图得出v(0,[0.85;?0.75])=0.4366，这比表1中的六阶概率界限0.3696要宽松。然而，图2(b)提供了整个集合X上的风险的可解释可视化。

**下载：** 下载高分辨率图像（60KB）
**下载：** 下载全尺寸图像

**图2.** 在t=T=5且阶数k=6时，SDE（21）的不安全上界函数v(0,x)的风险水平，初始化在圆盘X0（紫色边界）中，不安全区域Xu（红色月牙形）。
（关于此图例中图3(a)展示了根据动力学方程(22)从初始条件（洋红色点）出发的300个随机实现的两个平台在相空间中的位置和速度。左侧（蓝色）轨迹是平台1（x1,x?1）的运动情况，而右侧（红色）轨迹是平台2（x2,x?2）的响应。不安全区域Xu={x∈R4∣(x1?x2)2≤0.152}，描述了平台相互碰撞的潜在风险。图3(b)绘制了平台位置差x2(t)?x1(t)随时间t的变化情况，其中t∈[0,T]。图3(b)中的红色虚线表示不安全状态的开始，即x2(t)?x1(t)=0.15的水平集。图3(b)中的一些经验采样轨迹穿过了这条红线，表明初始点x0有可能进入不安全区域。表2列出了从x0出发进入Xu的概率上限，这些上限是通过对方程(18)进行有限阶截断计算得出的。

图3. 根据方程(22)得出的浮动平台的位置。

表2. 系统(22)的不安全概率上限。

6. 结论
本文提出了一种通过构建占位测度的线性规划（LPs）来分析随机过程不安全概率的方法。这些LPs受到随机微分方程（SDPs）的矩-SOS层次结构的约束，从而得到一系列收敛的不安全概率上限。通过修改目标函数，可以生成用于分析和运动规划的可视化风险轮廓。所提出的方法适用于满足假设A1–A8的任何随机过程。未来的工作将包括寻找冗余约束以优化和改进不安全概率的上限。其他工作还包括开发对风险敏感的随机控制方法的凸重构形式。