相比于常规能源，风电、光伏等可再生能源出力具有明显的随机性、波动性以及间歇性的特征，随着风电、光伏渗透率的不断增加，电力系统中长期规划运行与短期调度中都必须考虑风电、光伏出力的不确定性问题。场景分析（scenario analysis）是一种通过构建确定性场景来分析电力系统不确定性问题的方式，它是解决含可再生能源的电力系统优化规划运行问题的一种有效途径。本文将对现有的场景生成方法做一个梳理。首先，介绍随机优化和场景分析的基本概念和数学模型，其次对电力系统中常用的场景生成方法进行分类，并分别介绍方法和思路。

数学优化是应用数学的重要分支，在当今工程、商业、计算机科学等领域有着极其广泛的应用。数学优化的目的是在一定的约束条件下使得某个目标函数达到最优，其数学表达式如下：

min ⁡ f ( x )  s.t.  g i ( x ) ≤ 0 , i = 1 , … , m h j ( x ) = 0 , j = 1 , … , k x ∈ X ⊂ R n (1) \begin{array}{ll} \min & f(x) \\ \text { s.t. } & g_{i}(x) \leq 0, i=1, \ldots, m \\ & h_j(x)=0, j=1, \ldots, k\\ & x \in \mathcal{X} \subset \mathbb{R}^{n} \end{array} \tag{1} min s.t. ​f(x)gi​(x)≤0,i=1,…,mhj​(x)=0,j=1,…,kx∈X⊂Rn​(1) 式中， x x x 表示决策变量； X ⊂ R n \mathcal{X} \subset \mathbb{R}^{n} X⊂Rn 为决策变量的可行域； f ( x ) f(x) f(x) 为数学优化问题的目标函数； g ( x ) g(x) g(x) 和 h ( x ) h(x) h(x)分别表示数学优化问题的不等式约束和等式约束。

通常的理论研究均假设公式 (1) 中的参数都是确定性的已知数，然而在现实世界的实际应用中很多最优化决策问题都存在着大量的不确定性，即在 f ( x ) f(x) f(x) ， g ( x ) g(x) g(x) 或 h ( x ) h(x) h(x)中存在随机变量。因此，根据决策问题中是否存在随机变量，数学优化被分类为“确定性优化”和“不确定性优化”。随机优化(Stochastic optimization, SO)是一种对涉及不确定性优化问题进行建模的重要工具，其目标函数或约束条件中的部分参数是不确定的。SO的数学表达式为： min ⁡ f ( ω , x )  s.t.  g i ( ω , x ) ≤ 0 , i = 1 , … , m h j ( ω , x ) = 0 , j = 1 , … , k x ∈ X ⊂ R n (2) \begin{array}{ll} \min & f(\omega, x) \\ \text { s.t. } & g_{i}(\omega, x) \leq 0, \quad i=1, \ldots, m \\ & h_j(\omega, x)=0, \quad j=1, \ldots, k\\ & x \in \mathcal{X} \subset \mathbb{R}^{n} \end{array} \tag{2} min s.t. ​f(ω,x)gi​(ω,x)≤0,i=1,…,mhj​(ω,x)=0,j=1,…,kx∈X⊂Rn​(2)式中 ω \omega ω 表示随机变量。

尽管由于随机变量的存在导致SO无法直接求得确定性最优解，但是随机变量在概率的角度存在对应的概率分布，则使期望收益/损失达到最大/最小更加具有现实意义。

场景分析方法是在观测到随机变量的实现后再做决策的方法，依赖于对不确定性的充分观测。场景分析模型通常是求解SO问题的期望值模型，一般化数学形式如下： min ⁡ E [ f ( ω , x ) ]  s.t.  E [ g i ( ω , x ) ] ≤ 0 , i = 1 , … , m E [ h j ( ω , x ) ] = 0 , j = 1 , … , k x ∈ X ⊂ R n (3) \begin{aligned} &\min \mathbb{E}[f(\omega, x)] \\ \text { s.t. } &\mathbb{E}\left[g_{i}(\omega, x)\right] \leq 0, \quad i=1, \ldots, m \\ & \mathbb{E}\left[h_j(\omega, x)\right]=0, \quad j=1, \ldots, k\\ &\quad x \in \mathcal{X} \subset \mathbb{R}^{n} \end{aligned}\tag{3}  s.t. ​minE[f(ω,x)]E[gi​(ω,x)]≤0,i=1,…,mE[hj​(ω,x)]=0,j=1,…,kx∈X⊂Rn​(3)式中 E [ ⋅ ] \mathbb{E}[\cdot] E[⋅] 为期望算子。当确定了随机变量 ω \omega ω 的概率分布后，可以利用场景分析方法将该问题转换为确定性优化进行求解。基于场景分析（Scenario-based）的SO广泛应用于含风电场最优潮流问题和电力系统经济调度问题。

在实际问题中，公式 (3) 表述的基于场景分析的期望值模型可以由以下凸SO的期望值模型所表示： min ⁡ x ∈ X E P f ( ω , x ) = ∫ Ω f ( ω , x ) P ( d ω ) (4) \min _{x \in \mathcal{X}} \mathbb{E}_{\mathbf{P}} f(\omega, x)=\int_{\Omega} f(\omega, x) \mathbf{P}(d \omega) \tag{4} x∈Xmin​EP​f(ω,x)=∫Ω​f(ω,x)P(dω)(4)式中 P \mathbf{P} P 是 ω \omega ω 在空间 Ω \Omega Ω 上的概率测度。

为了求解公式(4)，现有研究设计出了许多针对不同实际问题的近似模型。其中 ζ \zeta ζ 结构概率测度是一种主流的近似方法，数学表达式如下： d F ( P , Q ) = sup ⁡ f ∈ F ∣ ∫ Ω f ( ω ) P ( d ω ) − ∫ Ω f ( ω ) Q ( d ω ) ∣ (5) d_{\mathrm{F}}(\mathbf{P}, \mathbf{Q})=\sup _{f \in \mathbf{F}}\left|\int_{\Omega} f(\omega) \mathbf{P}(d \omega)-\int_{\Omega} f(\omega) \mathbf{Q}(d \omega)\right| \tag{5} dF​(P,Q)=f∈Fsup​ ​∫Ω​f(ω)P(dω)−∫Ω​f(ω)Q(dω) ​(5)式中 sup ⁡ \sup sup 表示上确界。

通过 ζ \zeta ζ 结构概率测度，可以知道求解 ∫ Ω f ( ω ) P ( d ω ) \int_{\Omega} f(\omega) \mathbf{P}(d \omega) ∫Ω​f(ω)P(dω) 的期望最优值等价于求解 ∫ Ω f ( ω ) Q ( d ω ) \int_{\Omega} f(\omega) \mathbf{Q}(d \omega) ∫Ω​f(ω)Q(dω) 的期望最优值。因此，当公式(4)规模较大且难以完全描述的时，可以釆用 P \mathbf{P} P 的近似简化概率测度 Q \mathbf{Q} Q 求解模型。为了进一步求得简化概率测度 ，将对 ζ \zeta ζ 结构概率测度进行进一步分析。公式(5)其满足如下性质： ζ c ( P , Q ) : = d F c ( P , Q ) ≤ μ ^ c ( P , Q ) (6) \zeta_{c}(\mathbf{P}, \mathbf{Q}):=d_{\mathrm{F}_{c}}(\mathbf{P}, \mathbf{Q}) \leq \hat{\mu}_{c}(\mathbf{P}, \mathbf{Q}) \tag{6} ζc​(P,Q):=dFc​​(P,Q)≤μ^​c​(P,Q)(6)式中 μ ^ c \hat{\mu}_{c} μ^​c​ 为康托洛维奇泛函 (Kantorovich functional)，其表达式为： μ ^ c ( P , Q ) : = inf ⁡ { ∫ Ω × Ω c ( ω , ω ~ ) η ( d ( ω , ω ~ ) ) } (7) \hat{\mu}_{c}(\mathbf{P}, \mathbf{Q}):=\inf \left\{\int_{\Omega \times \Omega} c(\omega, \tilde{\omega}) \eta(d(\omega, \tilde{\omega}))\right\} \tag{7} μ^​c​(P,Q):=inf{∫Ω×Ω​c(ω,ω~)η(d(ω,ω~))}(7) 总何公式(6)与(7)可以得出结论：若 P \mathbf{P} P 与 Q \mathbf{Q} Q 的关系满足公式(7)，则“简化概率测度 Q \mathbf{Q} Q 求得SO模型的最优值式(8)”和“原始概率侧度 P \mathbf{P} P 求得SO模型的最优值公式(4)”是最接近的。换句话说，一旦原模型由于其概率测度 P \mathbf{P} P 规模庞大而无法求解时，可通过寻求满足公式(7)的简化概率测度 Q \mathbf{Q} Q 将原凸SO问题进一步化简为公式(8)，且简化后的模型经过计算得到的近似结果是与原模型是十分接近。 min ⁡ x ∈ X E Q f ( ω , x ) = ∫ Ω f ( ω , x ) Q ( d ω ) (8) \min _{x \in \mathcal{X}} \mathrm{E}_{\mathbf{Q}} f(\omega, x)=\int_{\Omega} f(\omega, x) \mathbf{Q}(d \omega) \tag{8} x∈Xmin​EQ​f(ω,x)=∫Ω​f(ω,x)Q(dω)(8)

在诸如电力系统随机优化等现实问题中，原始概率测度 P \mathbf{P} P 往往是满足某种概率分布的连续随机变量，无法直接对公式(4)进行求解。因此需要将概率测度 P \mathbf{P} P 离散化，即用数量庞大的有限个离散的确定性样本近似表示 ，这些离散样本 P ~ \mathbf{\widetilde{P}} P 即为“场景”(scenario)，而如何根据连续且解析的概率分布函数获得有限个数量庞大的样本近似表示概率测度 P \mathbf{P} P 的问题即为本文中所谓的 场景生成”。可以描述为： ( ξ s , ρ s ) s = 1 , 2 , … , S (\xi_s, \rho_s) \quad s = 1, 2, \ldots, S (ξs​,ρs​)s=1,2,…,S其中 ξ s \xi_s ξs​ 为一个场景，而 ρ s \rho_s ρs​ 为对应场景的概率。

例如：下图为采用均匀分布将原概率测度 P \mathbf{P} P 进行离散化的示意图，图中，曲线代表原概率密度分布，矩形条代表场景，矩形条的高度代表场景对应的概率，图片来自参考文献[3]。

由于计算复杂度的限制，等价于概率测度 P \mathbf{P} P 的场景 P ~ \mathbf{\widetilde{P}} P 往往十分庞大。因此需要采用一个仅含有少数典型场景所组成的概率测度 Q \mathbf{Q} Q 来进一步近似概率测度 。如何获得最佳的典型场景集的问题即为“场景削减”。

近年来，已经提出了许多方法实现场景生成。可以分为：基于抽样的方法、基于人工智能的方法、基于时序模型的方法等

基于抽样的方法通过对概率分布进行抽样并将输出样本作为生成的场景来得到离散场景。这类方法包括：蒙特卡罗（Monte Carlo, MC）方法、拉丁超立方体采样（Latin Hypercube Sampling, LHS）方法和Copula函数采样，Gibbs采样、马尔可夫链蒙特卡洛（Markov chain Monte Carlo, MCMC），样本平均近似（Sample Average Approximation, SAA）等。下面主要介绍三种方法：

基于MC的方法是最常用、最基本、也是最简单的场景生成方法，其过程可以总结为：

构建一个概率分布模型：

f ( v ; c , k ) = ( k c ) ( v c ) k − 1 exp ⁡ [ − ( v c ) k ] (9) f(v ; c, k)= \left(\frac{k}{c}\right)\left(\frac{v}{c}\right)^{k-1} \exp \left[-\left(\frac{v}{c}\right)^{k}\right]\tag{9} f(v;c,k)=(ck​)(cv​)k−1exp[−(cv​)k](9)式中： v v v 为实际风速， c c c 是尺度参数（scale parameter）， k k k是形状参数（shape parameter）

f ( s ; α , β ) = Γ ( α + β ) Γ ( α ) Γ ( β ) s ( α − 1 ) ( 1 − s ) ( β − 1 ) (10) f(s;\alpha,\beta)=\frac{\Gamma(\alpha+\beta)}{\Gamma(\alpha) \Gamma(\beta)}s^{(\alpha-1)}(1-s)^{(\beta-1)} \tag{10} f(s;α,β)=Γ(α)Γ(β)Γ(α+β)​s(α−1)(1−s)(β−1)(10)式中：其中 s s s 是太阳辐照度； α , β \alpha, \beta α,β 表示Beta分布的形状参数。

从概率分布模型中采样

生成满足概率分布模型的场景

根据上述步骤得到风速的场景后，可以利用“风速-风功率转换公式"得到风电功率场景 P W = { 0 , v < v in  or  v > v out  P ∗ , v n ≤ v ≤ v out  ( v − v in ) P ∗ / ( v n − v in ) , v in ≤ v ≤ v n (11) {{P}_{\text{W}}}=\left\{ \begin{matrix} 0, & v{v}_{\text{in}}}\text{ or }v{{v}_{\text{out }}} \\ {{P}_{*}}, & {{v}_{\text{n}}}\le v\le {{v}_{\text{out }}} \\ \left( v-{{v}_{\text{in}}} \right){{P}_{*}}/\left( {{v}_{\text{n}}}-{{v}_{\text{in}}} \right), & {{v}_{\text{in}}}\le v\le {{v}_{\text{n}}} \\ \end{matrix} \right. \tag{11} PW​=⎩ ⎨ ⎧​0,P∗​,(v−vin​)P∗​/(vn​−vin​),​vvout ​vn​≤v≤vout ​vin​≤v≤vn​​(11)式中： v n v_\text{n} vn​ 表示额定风速； v in {v}_{\text{in}} vin​ 和 v out  {v}_{\text{out }} vout ​ 分别表示切入风速和切出风速； P ∗ {P}_{*} P∗​ 和 P W {P}_{\text{W}} PW​ 分别表示额定功率和输出功率。

而辐照度则可以利用公式，转换为光伏输出功率场景： P PV = s A pv η pv (12) P_{\text{PV}}=sA_{\text{pv}}\eta_{\text{pv}} \tag{12} PPV​=sApv​ηpv​(12)式中：其中 P PV P_{\text{PV}} PPV​ 为光伏功率输出， A pv A_{\text{pv}} Apv​ 为光伏组件辐射面积， η pv \eta_{\text{pv}} ηpv​ 为光电转换效率。

至此我们就得到了风功率和光伏的出力场景，是不是很简单。如图所示:

那么引出三个问题：

如何得到Weibull分布和Beta分布的参数呢? 可以参考我的另一篇文章：风速Weibull分布和光伏Beta分布的参数拟合方法。

为什么不直接对风功率建模并采样呢？ 这是因为目前为止还没有一个被广泛接受的风功率精确数学解析表达式。或许你说风速的Weibull分布模型具有普遍的适用性，可以利用“风速-风功率转换公式"得到风功率的PDF，但是从风速转化到风功率的转换却具有较大的误差，不能代表风功率的概率模型。光伏功率就没那么麻烦，因为 s s s 与 P PV P_{\text{PV}} PPV​是线性关系。

假设每个时刻的参数相互独立，合理吗？ 不太合理，风功率/光伏功率作为一个关于时间连续的随机变量，在不同时刻的彼此存在相关性。这就要考虑到时间序列的相关性了，多数研究采用多元高斯分布构建时间相关关系。也可以利用马尔科夫链蒙特卡洛(MCMC)等方法，得到考虑时序自相关性的场景。若考虑时序相关性，场景示意图如下，采用的模型参数和上边两张图一样：

总结：基于MC的方法简单、高效，但是必须知道随机变量（风速和辐照度）的概率分布。当随机变量不服从一般分布时，PDF 通常难以获得。而且，MC方法不仅依赖于这种方法依赖于统计假设和 ，还需要大量的样本来匹配观察场景的统计特征，导致计算量大，而且不能很好地捕捉时间序列的自相关特征。

MC是抽样中的简单随机抽样法，有一个缺点是如果随机抽样的样本数小， 样本则难以反映随机变量的理论分布； 简单随机抽样法需要大量的循环，抽样数必须足够大，因而效率不高。为了解决简单随机抽样效率不高且精度较差的不足，拉丁超立方抽样LHS技术被广泛采用。 LHS属于分层抽样， 样本数量小且避免重复，仍可反映随机变量的理论分布，保证所有的抽样区域都被抽样点均匀覆盖。

假设 X 1 , X 2 , … , X K X_1,X_2,\ldots,X_K X1​,X2​,…,XK​ 为 K K K 个相互独立的随机变量，其中第 k k k 个随机变量的CDF F k ( X k ) F_k(X_k) Fk​(Xk​) 及其反函数 F k − 1 ( Y k ) F_k^{-1}(Y_k) Fk−1​(Yk​) 存在，数学表达式如下： Y k = F k ( X k ) , k = 1 , 2 , … , K X k = F k − 1 ( Y k ) (13) Y_k = F_k(X_k), k=1,2,\ldots,K\\ X_k =F_k^{-1}(Y_k) \tag{13} Yk​=Fk​(Xk​),k=1,2,…,KXk​=Fk−1​(Yk​)(13)随机变量 X k X_k Xk​ 的CDF如图所示。 设 N N N 代表生成的场景总数，LHS的抽样方法是，将图中的CDF的纵轴等分成 N N N 个连续递增的数值区间，区间宽度为 1 / N 1/N 1/N，将 Y Y Y 轴每个区间的重点 y k n y_k^n ykn​ 抽样，与中点 y k n y_k^n ykn​ 相对应的反函数 x k n x_k^n xkn​ 即为第 n n n 个样本，其数值等于： x k n = F k − 1 ( y k n ) = F k − 1 ( n − 0.5 N ) , n = 1 , 2 , … , N x_k^n=F_k^{-1}(y_k^n)=F_k^{-1}(\frac{n-0.5}{N}), n=1,2,\ldots,N xkn​=Fk−1​(ykn​)=Fk−1​(Nn−0.5​),n=1,2,…,N LHS和MC的区别在于： 具有抽样的“记忆”功能， 可以避免逆变换抽样法数据点集小而导致的重复抽样，它强制抽样过程中抽样点必须离散分布于整个抽样空间，正因为此 LHS的抽样点的数量不必太多即可大致刻画随机变量的概率分布。但同样的，LHS多用于静态场景的生成，当要生成动态场景时没有考虑各时段的相关性。

上述抽样方法多是对单一风电场/风电机组生成场景。当我们要考虑到一个区域内的多个风电场/风电机组时，如何构建具有相关性的多风场联合出力场景是一个值得研究的问题。但是多元联合分布函数的场景不易构造和模拟，且通常要求所有的边缘分布服从同一分布。因此，引入了一种Copula函数(连接函数)描述空间相邻风电场间的相关性，基于 Copula 函数生成风电场出力场景的方法对边缘分布没有限制，能捕捉变量之间非线性、非对称性以及尾部相关关系。

Copula 函数可将多个变量的联合分布和相应的边缘分布连接到一起，该函数可灵活获取联合分布函数。

多元分布的 Sklar 定理

令 F ( x 1 , x 2 , … , x N ) F(x_1,x_2,\ldots,x_N) F(x1​,x2​,…,xN​) 为具有边缘分布为 F 1 ( x 1 ) , F 2 ( x 2 ) , … , F N ( x N ) F_1(x_1),F_2(_x2),\ldots,F_N(x_N) F1​(x1​),F2​(x​2),…,FN​(xN​) 的联合分布函数，则必定存在一个 Copula 函数 C ( F 1 ( x 1 ) , F 2 ( x 2 ) , … , F N ( x N ) ) C(F_1(x_1),F_2(x_2),\ldots,F_N(x_N)) C(F1​(x1​),F2​(x2​),…,FN​(xN​))，满足 F ( x 1 , x 2 , … , x N ) = C ( F 1 ( x 1 ) , F 2 ( x 2 ) , … , F N ( x N ) ) F(x_1,x_2,\ldots,x_N)=C(F_1(x_1),F_2(x_2),\ldots,F_N(x_N)) F(x1​,x2​,…,xN​)=C(F1​(x1​),F2​(x2​),…,FN​(xN​)) 式中 F ( x 1 , x 2 , … , x N ) F(x_1,x_2,\ldots,x_N) F(x1​,x2​,…,xN​) 代表随机变量 x 1 , x 2 , … , x N x_1,x_2,\ldots,x_N x1​,x2​,…,xN​ 的联合分布函数。

Copula函数总体包括 Elliptical Copulas函数族和 Archimedean Copulas 函数簇，前者包括正态Copula、t-Copula，后者包括Clayton copula、Frank Copula和Gumbel Copula。例如：

二元正态 Copula 分布函数的数学表达式为： C G a ( u , v , ρ ) = ∫ − ∞ Φ − 1 ( u ) ∫ − ∞ Φ − 1 ( v ) 1 2 π 1 − ρ 2 ⋅ exp ⁡ [ − s 2 − 2 ρ s t + t 2 2 ( 1 − ρ 2 ) ] d s   d t (14) \begin{gathered} C_{G a}(u, v, \rho)=\int_{-\infty}^{\Phi^{-1}(u)} \int_{-\infty}^{\Phi^{-1}(v)} \frac{1}{2 \pi \sqrt{1-\rho^{2}}} \cdot\exp \left[-\frac{s^{2}-2 \rho s t+t^{2}}{2\left(1-\rho^{2}\right)}\right] \mathrm{d} s \mathrm{~d} t \end{gathered} \tag{14} CGa​(u,v,ρ)=∫−∞Φ−1(u)​∫−∞Φ−1(v)​2π1−ρ2 ​1​⋅exp[−2(1−ρ2)s2−2ρst+t2​]ds dt​(14)式中 Φ − 1 \Phi^{-1} Φ−1 表示标准正态分布的分布函数逆函数； ρ \rho ρ 为变量间线性相关系数。

二元t-Copula分布函数的数学表达式为： C t ( u , v ; ρ , k ) = ∫ − ∞ t k − 1 ( u ) ∫ − ∞ t k − 1 ( v ) 1 2 π 1 − ρ 2 ⋅ [ 1 + s 2 − 2 ρ s t + t 2 k ( 1 − ρ 2 ) ] − k ( k + 2 ) 2   d s   d t (15) \begin{aligned} &C_{t}(u, v ; \rho, k)=\int_{-\infty}^{t_{k}^{-1}(u)} \int_{-\infty}^{t_{k}^{-1}(v)} \frac{1}{2 \pi \sqrt{1-\rho^{2}}} \cdot\left[1+\frac{s^{2}-2 \rho s t+t^{2}}{k\left(1-\rho^{2}\right)}\right]^{-\frac{k(k+2)}{2}} \mathrm{~d} s \mathrm{~d} t \tag{15} \end{aligned} ​Ct​(u,v;ρ,k)=∫−∞tk−1​(u)​∫−∞tk−1​(v)​2π1−ρ2 ​1​⋅[1+k(1−ρ2)s2−2ρst+t2​]−2k(k+2)​ ds dt​(15)式中： ρ \rho ρ 为变量间线性相关系数； k k k 为自由度； t k − 1 t_{k}^{-1} tk−1​ 表示自由度为 k k k 的一元 t 分布函数的逆函数。

不同的Copula 用于捕捉随机变量间具有差异的尾部相关关联。如：

基于Copula的场景生成方法主要步骤如下：

对于二元Copula，针对模型的评价在此引入经验 Copula 函数的概念。定义如下：设 ( x i , y i ) ( i = 1 , 2 , … , n ) (x_i,y_i)(i=1,2,\ldots,n) (xi​,yi​)(i=1,2,…,n) 为变量 ( X , Y ) (X,Y) (X,Y) 的样本， X , Y X,Y X,Y 的经验分布函数定义为 F ( x ) F(x) F(x) 和 G ( y ) G(y) G(y)，则经验 Copula 为： C N ( u , v ) = 1 n ∑ i = 1 n I [ F ( x i ≤ u ) ] I [ G ( y i ) ≤ v ] (16) C_N(u,v)=\frac{1}{n}\sum_{i=1}^nI_{[F(x_i\leq u)]}I_{[G(y_i)\leq v]} \tag{16} CN​(u,v)=n1​i=1∑n​I[F(xi​≤u)]​I[G(yi​)≤v]​(16)其中， I I I 为示性函数。

欧式距离是由经验分布与理论分布函数值的差值平方和计算而来，如下式所示 d = ∑ i = 1 n ∣ C n ( u i , v i ) − C ( u i , v i ) ∣ 2 (17) d=\sum_{i=1}^n|C_n(u_i,v_i)-C(u_i,v_i)|^2 \tag{17} d=i=1∑n​∣Cn​(ui​,vi​)−C(ui​,vi​)∣2(17)其中： C n C_n Cn​ 和 C C C 分别为经验 Copula 函数和拟合的 Copula 函数的累积分布函数。 欧式距离越小，则说明该 Copula 函数的拟合精确度越高，由此可选择最优Copula 函数。

算例展示：

人工智能方法一般是数据驱动的方法，生成场景的质量取决于历史观察样本。它能够捕捉到可再生能源的多样性和时间序列的复杂非线性关系。人工智能方法不考虑有关分布函数的任何假设，也不依赖于任何采样技术。主要方法由：径向基函数神经网络（RBFNN）、人工神经网络（ANN）、长短时记忆（LSTM）网络、变分自编码器(VAE)和生成对抗网络（GANs）等。本文主要介绍基于GAN的可再生能源场景生成方法：

GANs是由Goodfellow等人于2014年提出的一种无监督生成式模型，其包含有生成器G和判别器D两个结构独立的深度神经网络，通过生成器和判别器之间的博弈学习来提高生成器生成样本的真实。GANs场景生成方法的基本结构如图所示： 生成器 G ( ⋅ ; θ g ) G(\cdot;\theta_g) G(⋅;θg​) 用于拟合输入噪声分布 P Z P_Z PZ​ 与真实样本分布 P d P_d Pd​ 之间的映射函数，其目标是使其生成的样本能够符合真实样本分布。而判别器 D ( ⋅ ; θ d ) D(\cdot;\theta_d) D(⋅;θd​) 通常是一个分类器，用于拟合高维样本映射至低维度分类值0~1的映射函数，其目标是区分真实样本和生成样本。生成器与判别器之间形成一个动态的“博弈”过程，最终实现纳什均衡：生成器生成的样本足够真实，而判别器无法区分真实样本与生成样本的真假，对于，真实样本与生成样本的判别器输出值都基本接近1/2。

在训练过程中，将一批噪声 z ∼ P Z z\sim P_Z z∼PZ​送到生成器之中，输出的生成样本为 G ( z ; θ g ) G(z;\theta_g) G(z;θg​) 。同时，来自真实数据 x ∼ P d x\sim P_d x∼Pd​ 或生成数据 G ( z ; θ g ) G(z;\theta_g) G(z;θg​) 的一批样本被送到判别器，该判别器目标是识别数据来自何处。根据不同的数据源，判别器的输出可以表示为： { p real  = D ( x ) p fake  = D ( G ( z ) ) (18) \left\{\begin{array}{l} p_{\text {real }}=D(x) \\ p_{\text {fake }}=D(G(z)) \end{array}\right. \tag{18} {preal ​=D(x)pfake ​=D(G(z))​(18)式中 p r e a l p_{real} preal​ 为判别器对真实样本的判别概率，而 p f a k e p_{fake} pfake​为 判别器对生成样本的判别概率，为方便描述省略生成器和判别器神经网络的权重 θ g , θ d \theta_g,\theta_d θg​,θd​。对于给定的判别器，生成器希望增大判别器对 G ( z ) G(z) G(z) 的概率输出 P f a k e P_{fake} Pfake​，因为较大的判别器输出意味着样本更加真实。对于给定的生成器，判别器则要在最小化 p f a k e p_{fake} pfake​ 的同时寻求的 p r e a l p_{real} preal​ 最大化。因此，这两个神经网络的损失函数可以分别定义为： L G = − E z ∼ P Z [ D ( G ( z ) ) ] L D = − E x ∼ P d [ D ( x ) ] + E z ∼ P Z [ D ( G ( z ) ) ] (19) \begin{gathered} \quad L_{G}=-\mathrm{E}_{z \sim P_{Z}}[D(G(z))] \quad \\ L_{D}=-\mathrm{E}_{x \sim P_{d}}[D(x)]+\mathrm{E}_{z \sim P_{Z}}[D(G(z))] \tag{19} \end{gathered} LG​=−Ez∼PZ​​[D(G(z))]LD​=−Ex∼Pd​​[D(x)]+Ez∼PZ​​[D(G(z))]​(19)在生成器和判别器之间建立一个博弈关系，使两个网络可以同时训练，因此GAN的训练目标是二人极小极大值博弈： min ⁡ G max ⁡ D V ( G , D ) = E x ∼ P d [ D ( x ) ] − E z ∼ P Z [ D ( G ( z ) ) ] (20) \begin{aligned} \min _{G} \max _{D} V(G, D) &=\mathrm{E}_{x \sim P_{d}}[D(x)]-\mathrm{E}_{z \sim P_{Z}}[D(G(z))] \tag{20} \end{aligned} Gmin​Dmax​V(G,D)​=Ex∼Pd​​[D(x)]−Ez∼PZ​​[D(G(z))]​(20)在训练结束后利用生成器即可生成可再生能源的场景。使用GAN生成可再生能源随机场景己经得到了初步的研宄，其可以基于历史数据充分挖掘可再生能源的不确定性生成符合可再生能源出力特性的随机场景。并且，基于GAN的方法是一种数据驱动方法，不需要依赖于统计假设和抽样技术。在GAN的基础上，引申出了许多改进模型，如：DCGAN，WGAN，WGAN-GP，与条件生成对抗网络结合，还有基于联邦学习的分布式场景生成方法等等。

下面是简单的复现了一些论文：

基于WGAN的光伏场景生成

基于条件生成对抗网络的风功率日前场景：

ARMA、ARIMA方法、状态空间方法。

ARMA 模型又称为自回归-滑动平均模型，假设一个时间序列在某时刻的值可以用 p p p 个历史观测值的线性组合加上一个白噪声序列的 q q q 项滑动平均来表示，则该时间序列即为 ARMA(p，q)过程，模型具体如下表示 未完待续

矩匹配

聚类方法 常见的聚类划分算法一般有经典的 K-means 算法、K-medoids 算法、模糊 C 均值聚类(FCM)、DBSCAN等。其中 K-means 算法是最常见并被广泛应用的聚类算法，它的基本步骤是首先指定 K 个数据为算法的初始聚类中心，随后计算其余数据点到这 K 个中心的距离，并按距离大小依次指定各数据点分配到最近的 K 个质心并形成 K 个类，接着重新计算 各类别中的质心，重复以上步骤，直至各类中质心不再发生变化，算法结束。然而，如果直接对历史观测值进行聚类生成场景，则得到的仅仅是历史已经发生过的某些场景，缺乏预见性和多样性。 因此，聚类方法多用于场景削减。

其他组合方法

未完待续

随机抽样过程将得到组等概率的随机场景，在这些场景中，大量场景是相似的，从场景所能提供信息量大小的角度讲，对这些场景的精确分辨意义不大，并且会显著影响计算效率。因而，在场景抽样的基础上，合并部分场景形成不等概率的场景集合是十 分有必要的。此时，一种场景可能替代了抽样形成的多种相近场景，从而提高了场景的描述效率。

现有场景缩减问题的常用算法包括：后向缩减法(Backward Reduction， BR)、 前向选择法(Forward Selection, FS)、场景树构建法(Scenario Tree Construction, STC )、聚类法等。

未完待续

参考文献： [1] Li J, Zhou J, Chen B. Review of wind power scenario generation methods for optimal operation of renewable energy systems[J]. Applied Energy, 2020, 280: 115992. [2] 黎静华, 韦化, 莫东. 含风电场最优潮流的Wait-and-See模型与最优渐近场景分析[J]. 中国电机工程学报, 2012, 30(22): 15-24. [3] Chen Y, Wang Y, Kirschen D, et al. Model-free renewable scenario generation using generative adversarial networks[J]. IEEE Transactions on Power Systems, 2018, 33(3): 3265-3275. [4] 黎静华, 文劲宇, 程时杰, 等. 考虑多风电场出力 Copula 相关关系的场景生成方法[J]. 中国电机工程学报, 2013, 33(16): 30-36. [5] Li Y, Li J, Wang Y. Privacy-preserving spatiotemporal scenario generation of renewable energies: A federated deep generative learning approach[J]. IEEE Transactions on Industrial Informatics, 2022, 18(4): 2310-2320. [6] Dupačová J, Gröwe-Kuska N, Römisch W. Scenario reduction in stochastic programming[J]. Mathematical Programming, 2003, 95(3): 493-511. [7] 董骁翀, 孙英云, 蒲天骄. 基于条件生成对抗网络的可再生能源日前场景生成方法[J]. 中国电机工程学报, 2020, 40(17): 5527-5536.