Reliability-based Design Optimization Considering Correlated Random Variables
WANG Qianrong JIANG Chao FANG Teng
State Key Laboratory of Advanced Design and Manufacturing for Vehicle Body，Hunan University，
Changsha，410082
Abstract： When solving structural reliability⁃based design optimization problems with nonlinear cor⁃ relation，current method for correlated variables might lead to inaccurate optimization results in some sit⁃ uations. Based on Copula function，this paper proposed a method to serve as an effective tool for structur⁃ al reliability ⁃ based design optimization problems where nonlinear correlation existed. The proposed method estimated the parameters of alternative Copula functions according to known samples and select⁃ ed the optimal Copula function by AIC criterion. Therefore，joint distribution function of variables was established to solve structural reliability ⁃ based design optimization problems. Finally，two numerical problems were used to demonstrate the validity of proposed method. Influences on optimization results caused by different Copula were also discussed.
Key words： Copula function；structural reliability ⁃ based design optimization；correlation；AIC cri⁃ terion 0 引言在工程问题中，加工尺寸、外载荷、材料参数等往往存在不确定性，在设计时若不考虑这些不确定性，可能得到不可靠的设计结果 。基于可靠
［1］性的优化设计（ reliability ⁃ based design optimiza⁃ tion，RBDO）充分考虑了各种不确定性，将不确定性变量作为随机变量来处理，从而可以得到满足可靠性要求的设计结果。目前， RBDO已经成为一种十分重要的设计方法，且该领域已发展出一系列求解方法。ZHUANG等 在求解RBDO问
［2］题中不断更新响应面，保证了收敛性和精度； LI⁃
收稿日期： 2017－08－23
基金项目：国家自然科学基金资助项目(51222502，11172096）；湖南省杰出青年基金资助项目(14JJ1016) ANG等 提出了单循环方法，采用K⁃K⁃T条件取
［3］代内层的可靠性分析问题，将两层的优化问题转换成单循环问题； YI等 将序列近似规划策略应用
［4］
到RBDO中，通过构建并求解一系列近似规划子问题来获得RBDO问题的解； DU等 提出了SO⁃
［5］ RA（sequential optimization and reliability assess⁃ ment）方法，通过运用上一循环的可靠性分析信息构造漂移向量，将RBDO问题转换成一系列确定性的优化问题，有效地提高了计算效率，同时具有较高的计算精度并能保证收敛性。
现有的RBDO方法大多假设各随机变量相互独立，然而，在实际工程问题中，很多变量之间具有显著的相关性。现有研究表明，变量之间的相关性对结构的可靠性和优化结果可能产生十分显著
的影响 。在RBDO问题中充分考虑变量间的相
［6］关性，可进一步提高优化模型的准确性，并避免产生不可靠的优化结果。目前， Nataf变换是处理随机变量相关性的主要方法，已被应用于RBDO领域。CHENG 等 利用 Nataf 变换求解RBDO 问
［7］
题。DU 等 提出了一种考虑相关区间变量的
［8］
RBDO方法。Nataf变换仅能描述变量之间的线性相关性，对于一些实际问题中广泛存在的非线性相关性，应用Nataf变换可能会出现比较大的误差，这使得Nataf变换的应用范围受到了限制。
近年来，在不确定性分析领域发展出了一种处理相关性的重要数学工具—— Copula函数。通过 Copula函数，只需确定变量的边缘分布和相关类型，即可方便地构建出联合概率分布函数。此外， Copula函数可描述变量之间不同类型的非线性相关性，并能捕捉到一些重要的信息，如尾部相关性等，这使得Copula函数近年来被广泛用于解决结构可靠性分析问题 。虽然Copula函数在
［9⁃10］可靠性分析方面引起了重视，但将其应用到RB⁃ DO领域的相关工作仍很少。CHOI 等 利用
［ 11⁃12 ］ Copula函数构建联合概率分布函数，并首次将Copula函数用于求解RBDO问题。然而，上述工作只使用了Copula函数的一种类型，即Gaussian Copula函数来描述变量间的相关性。在本质上， Gaussian Copula函数与Nataf变换是等效的，仅能描述变量间的线性相关性 。若能根据样本的实
［13］际分布情况选择出描述变量之间相关类型的最优Copula函数，则可更加准确地建立数学模型，进一步提高不确定性分析精度。
本文提出了一种基于Copula函数的RBDO方法，为存在复杂参数相关性的结构可靠性优化问题提供了一种有效工具。
1 Copula函数基本原理
SKLAR 指出，任一联合概率分布函数都可
［14］
被分解为若干个边缘分布和一个Copula函数，该Copula函数可描述变量之间的相关性。由此看出， Copula 函数实际上是一种将多个变量x1, x2, ⋯, xn的联合概率分布函数F ( x1, x2, ⋯, xn )与它们各自的边缘分布函数F1 ( x1 ) , F2 ( x2 ) , ⋯, Fn ( xn )连接在一起的函数。根据Sklar定理可知，联合概率分布函数与各边缘分布函数之间满足：
F ( x1, x2, ⋯ xn ) = C ( F1 ( x1 ) , F2 ( x2 ) , ⋯, Fn ( xn )) （ 1）
C即为描述 F1 ( x1 ) , F2 ( x2 ) , ⋯, Fn ( xn )之间相关性的 Copula 函数。令 uk = Fk ( xk ), k = 1, 2, ⋯, n，显然， C ( u1, u2, ⋯, un )是一个边缘分布均服从[ 0, 1 ]均匀分布的多元分布函数。
若已知各变量的边缘分布和连接它们的Cop⁃ ula函数，则可根据Sklar定理求解出联合概率分布函数。联合概率分布函数的概率密度
∂ C ( F1 ( x1 ) , F2 ( x2 ) , ⋯, Fn ( xn ))
n f ( x1, x2, ⋯, xn ) = = ∂ x1 ∂ x2⋯∂xn
n c ( u1, u2, ⋯, un ) f i ( xi ) （ 2）
=1
∂ C ( u1, u2, ⋯, un )
n c ( u1, u2, ⋯, un ) =
∂∂ u1 u2⋯∂un
1.1 Copula函数的相关性测度
两个变量之间，若它们的变化趋势一致，则称这两个变量正相关，反之则称负相关。由此可定义一致性的概念： ( x1, y1 )和( x2, y2 )为( X, Y )的两组观测值，若(- x1 x2 ) ( y1 - y2 ) >0 ，称 ( x1, y1 ) 和
( x2, y2 )是一致的；若(- x1 x2 ) ( y1 - y2 ) < 0，则称
( x1, y1 ) 和 ( x2, y2 )是不一致的。HOLLANDER等 基于一致性的概念提出了Kendall 秩相关系
［］ 15
数的定义。假设对连续的随机向量( X, Y )进行观测，得到一个由 N组观测值组成的样本集{( x1, y1 ) ,( x2, y2 ) , ⋯,( xN, yN ) }；将这些观测值两两匹配，得到C 项组合，并将其分为两部分，即2N C =+ c d， c表示一致组合的数量， d表示不一致2N组 合 的 数 量 ，由 此 可 定 义 样 本{( x1, y1 ) ,( x2, y2 ) , ⋯,( xN, yN ) } 的 Kendall 秩相关系数：
c - d τ = （ 3）
c + d
若给定一个Copula函数，则变量之间的相关性也随之确定。Kendall秩相关系数τ可由相应的Cop⁃ ula函数C ( u1, u2 )根据下式计算［ 16⁃17 ］ ：
∫∫ 1 1 τ =4 C ( u1, u2, θ ) dC ( u1, u2, θ )- 1 （ 4）
0 0
式中， θ为该Copula函数的相关性参数。
式（ 4）中不含随机变量边缘分布的表达式，这也表明了对于某个Copula函数，其Kendall秩相关系数不依赖于边缘分布。
在实际工程问题中，变量之间广泛地存在着尾部相关性。由于Kendall秩相关系数τ 和 Pearson相关系数（即线性相关系数ρ）都无法描述尾部相关性，因此引入尾部相关系数。若X1、X2为两个连 续 型 随 机 变 量 ，其 边 缘 分 布 分 别 为F1 ( X1 )、F2 ( X2 )，用 Copula函数表示其联合概率分布函数为C ( F1 ( X1 ) , F2 ( X2 ) )，其尾部相关系数［］ 18 可由下式定义： λU = lim p [> X2 F -1 ( t )> |X1 F -1 ( t )]=
2 1 t → 1-
1 - 2t + C ( t, t ) lim （ 5）
1- t t → 1
C ( t, t ) λL = lim p [< X2 F ( t )< |X1 F ( t ) ] = lim
-1 -1
2 1
t t → 0+ t → 0+
（ 6）式中， λU为上尾相关系数； λL为下尾相关系数； p为事件发生的概率； F 为边缘分布F2的逆函数； t → 1-表示t无限趋-1
2
近于1的左侧。
尾部相关系数可用来表示当一个变量为极值时，另一个变量也出现极值的可能性。尾部相关性可能对结构的可靠性造成显著的影响 ，因
［］ 9⁃10此，在RBDO问题中考虑尾部相关性可提高优化的精度。
1.2 常用的二维Copula函数分类
为方便起见，在本文中，仅考虑二维相关的情形，即多维随机变量中仅有两个变量之间具有相关性。当有两个以上的变量之间两两相关时，可利用 Vine⁃Copula函数将多维Copula函数分解成若干个二维Copula函数。
目前常用的二维Copula 函数见表1，分别为Gaussian Copula（可等效为Nataf变换）、t Copula、Clayton Copula、Gumbel Copula 和 Frank Copula。表1公式中， u1、u2为变量名称， tν和Φ分别表示t分布和高斯分布函数。为了更直观地理解各Copula函数的性质，图1给出了在秩相关系数相同（ τ= 0.5）、边缘分布均为标准正态分布、仿真次数均为2 000的情况下，由5种不同的Copula函数仿真得到的随机变量的散点图。从图1中可以看出，即使具有相同的秩相关系数和相同的边缘分布（标准正态分布），两个随机变量也可呈现不同的相关模式。例如， Gaussian Copula、t Copula 和 Frank Copula的分布具有对称性，且t Copula具有明显的尾部相关性。Gumbel Copula和Clayton Copula能描述非对称的相关模式，其中前者具有明显的上尾相关性，而后者则有明显的下尾相关性。可表1 常用的二维Copula函数信息［］
10
Tab.1 General information of two-dimension
Copula function
图各1 Copula函数仿真散点图
Fig.1 Simulation scatter diagram of Copula functions以看到， Copula函数能够描述变量之间不同的非线性相关类型。通过Copula函数，可更准确地建立随机变量联合概率分布函数。
2 基于Copula函数的RBDO求解算法RBDO问题通常可表述如下［］ 3 ：
式中，为目标函数； f gi为第i个可靠性约束； m为约束总个数；均值向量； d为nd维确定性设计向量； P为nP维随机参数向量； X为nX dd维随机向量； L、 U分别为设计变量μX为X的
dj 的下界和上界； μL Xp、μU Xp分别为μXp的下界和上界； j j Φ为标
准正态分布累积分布函数； βti为第i个约束的目标可靠度。2.1 联合概率分布函数的构建
在求解RBDO的过程中，当变量具有相关性时，需要利用联合概率分布函数进行求解。在实际工程问题中，通常只能获得各变量的边缘分布函数。Sklar定理指出，联合概率分布函数可被分解成边缘分布函数和一个描述变量之间相关性的Copula函数，这提供了一种构建联合概率分布函数的思路：分别获得随机变量的边缘分布及连接两边缘分布的Copula函数，从而得到联合概率分布函数F ( x1, x2, ⋯, xn )。
本文将 Gaussian Copula、t Copula、Clayton Copula、Gumbel Copula 和Frank Copula作为备选Copula函数，基于已知样本，通过极大似然法对这
5种 Copula函数进行参数估计，再由AIC（akaike information criterion）信息准则 选择出拟合效果
［］ 19
最好的Copula函数，确定出变量之间的相关类型。若随机变量X1、X2的边缘累积分布函数分别为F1 ( X1 ) , F2 ( X2 ) ，已知样本集合为 {( x11, x21 ) , ( x12, x22 ) , ⋯,( x1l, x2l ) }，共l组样本。对备选Cop⁃ ula函数建立似然对数函数：
l
L = lnc ( F1 ( x1i ), F2 ( x2i ) |θ ) （ 8）
=1
∂ C ( u1, u2 )
2 c ( u1, u2 ) =
∂∂ u1 u2
式中， θ为备选Copula函数的相关性参数（ t Copula有两个参数： θ和ν）。
由极大似然原理，参数估计值θ̂ = max L （ 9）
θ
对5种备选函数均使用极大似然法求得其参数的估计值后，可用AIC准则或BIC准则择取最优Copula函数：
l
CAIC = -2 lnc ( F1 ( x1i ), F2 ( x2i ) |θ̂ ) + 2r （ 10）
=1
式中，为相应r Copula函数中参数的数目。
CAIC值越小，说明该Copula函数对样本拟合越准确。由Sklar定理可知，当变量之间的Copula函数确定后，其联合概率分布函数也随之确定。
2.2 约束可靠性分析
约束可靠性分析是求解RBDO问题中的重要环节，用于分析当前设计变量d和μX是否满足概率约束。第i个约束gi的失效概率
∫ ⋯∫ p = prob { gi ( dZ , ) ≤ 0 }= fz ( ) dZ （ 11） f i Z
gi ( d, Z )≤0
Z =[ XP ]
式中， f ( )为联合概率密度函数; prob表示事件发生的可能
z Z
性。
一次二阶矩法（ first⁃order reliability method， FORM） 是进行可靠性分析常用的一种高效方
［］ 20
法。该方法需将原功能函数gi从Z空间映射至标准正态空间（即U空间），得到功能函数Gi，并求解如下优化问题： muin u
ï （ 12） s.t. Gi ( d, u )= 0
ï
其中， u表示待优化的随机变量和随机参数。
该优化问题的解为u∗即最可能失效点（ most probable point，MPP），进而可得到第i个约束的可靠度指标βi = u∗ 。若βi ≥ βti ，则说明概率约束 满足要求，反之则不满足。由于该方法用可靠度指标β来度量可靠性，故也称为可靠度指标法（ reli⁃ ability index analysis，RIA）。
此外，在FORM方法中还有另一种可靠性分析方法——功能度量法（ performance measure ap⁃ proach，PMA），该方法同样可以表述为一个优化问题：
min Gi ( d, u )
u ï （ 13） s.t. u = βti
ï该优化问题的解 u∗ 即为 MPP 点。若Gi ( d, u∗ ) ≥ 0，则说明满足概率约束，反之则说明不满足。将u∗ 从U空间映射回Z空间后，记为ZMPPi。对于RBDO问题，在大多数情况下， PMA法相对于RIA法效率更高、稳定性更好，且较少依赖于随机变量的分布类型。
在使用RIA法或PMA法进行可靠性分析时，需要将约束从Z空间映射至U空间，当变量之间存在相关性时，需要用 Rosenblatt 变换进行全概率变换。对于二维问题， Rosenblatt变换和逆变换如下： 时Z2的条件分布函数。
求解条件分布函数需要先得到联合概率分布函数，而在实际工程问题中，通常不能准确地得到联合概率分布函数，这使得Rosenblatt变换的应用受到了限制。通过2.1节中提出的方法可利用Copula函数方便地得到条件分布函数： 显 式 的 表 达 式 ，此 时 可 求 解 非 线 性 方 程h21 ( r1, r2 ) = y2得到数值解。
2.3 设计变量的优化
RBDO本质上是个两层嵌套优化问题，其外层是设计变量的优化，内层为约束可靠性的分析。嵌套优化使得RBDO求解效率较低，尤其是当随机
变量的个数增加时，所需的计算量显著增大。为减小运算量，本文在文献［ 5 ］的基础上，构建了一种解耦算法，对含有非线性相关性的结构可靠性优化设计问题进行高效求解。该方法由若干个迭代步组成，在每个迭代步中，首先对各约束进行可靠性分析，利用MPP点构造一漂移向量，将原问题转换为等效的确定性优化问题，通过求解该问题更新设计变量。由此，将嵌套优化解耦成序列进行的可靠性分析和确定性优化问题。在第k次迭代中，原问题被转化成如下形式： 通过构造漂移向量 S( ，将可靠性约束
k ) i prob { gi ( d, X, P )≤0}≤ p fi 转化为确定性约束gi ( d,μz - S( )≥ 0。由此在求解RBDO问题时，
k ) i只需交替进行可靠性分析和设计变量优化，通常在较少次数的迭代中即可收敛，具有很高的优化效率。
综上所述，本文RBDO方法的计算流程如下。（ 1）根据实际情况构造如式（ 4）形式的RBDO问题。
（ 2）根据已知样本对随机变量进行相关性分析，对存在相关性的变量利用极大似然法求得估计参数值。
（ 3）利用AIC准则求出拟合度最佳的Copula函数，从而求得联合概率分布函数。
（令4） k= 1，给定设计变量的初值d( 、μ( （可0) 0)
X用确定性优化结果作为迭代的初值，从而加快收敛，减小计算量）。
（ 5）在第k次迭代中，利用 Rosenblatt 变换将每一个约束映射至U空间进行可靠性分析，求解式（ 19），得到MPP点Z( 。k - 1)
MPPi
（ 6）构造漂移向量S(。
k ) i
（ 7）求解确定性优化问题式（ 17），得到更新后的设计变量d(、μ( 。k ) k )
X
（ 8）计算εk = | [ d( μ( k ) ]- [ d( k-1 ) μ( k-1 ) ] | /
k )
X X
[ d( k-1 ) μ( k-1 ) ]
X
若 εk > 10- 3，判定为不收敛，则k+ 1→ k，转到步骤 其中， X1 ~ N ( μX1 ,0.52 ), X2 ~ N ( μX2 , 0.52)，优化的初始点( μX1 , μX2 ) = ( 0, 0 )。和的X1 X2 500组样本见图2，从图中可知X1和X2之间具有明显的相关性。按照2.1节的方法对样本进行参数估计和相关性分析，可以得到各备选Copula函数的参数估计值和AIC值，见表2。从表中可知，用Gumbel Copula拟合样本得到的AIC值是各备选Copula函数中最小的，为最优Copula函数。实际上，从样本的散点图亦可知，该组样本具有明显的上尾相关性，而下尾相关性相对不明显。5种备选Copula函数中，只有Gumbel Copula具有这样的性质，因而其拟合效果也最好。根据参数估计和相关性分析所得到的最优Copula函数及其参数值，即可求解该RB⁃ DO问题。 图2 X1和X2的500组样本
Fig.2 500 samples of X1，X2
表2 X1、X2之间备选Copula参数估计值及其AIC值Tab.2 Estimate of Copula parameter and AIC value
between X1、X2
虽然已经得到最优Copula函数，但是为了分析不同类型的Copula函数对于优化结果的影响，本文仍使用3种 Copula函数进行RBDO分析： ① Gumbel Copula， θ = 1.45（最优 Copula 函数）； ② Gaussian Copula， θ = 0.445 ； ③ Clayton Copula， θ = 0.396。如表3所示，优化均经过4个迭代步即达到收敛，说明本文方法具有较好的收敛性和计算效率，但3种相关情形下得到的优化解具有明显的差异。为了解释优化结果产生差异的原因，将3种相关情形下的优化结果和各自的β⁃ circle 画在平面直角坐标系中，见图3。β⁃ circle 是指在U空
表3 算例一优化结果 Tab.3 Optimization results of example 1
间中所有到原点距离为β的点映射至X空间所得点集。在本算例中， X1和X2都为正态分布且标准差σ = 0.5。若X1和 X2相互独立，则β⁃ circle 应该是一个半径R=σβt= 1.5的圆；若X1和X2之间存在线性相关性，则β⁃ circle演化成长轴与X轴成45°角的椭圆。当用Copula函数描述变量间的相关性时， β⁃ circle则为其他类型的形状，该形状与Copu⁃ la函数的种类及其参数有关。在优化过程中， β⁃ circle不与各个约束相交，以保证设计点与约束之间存在安全距离，从而满足可靠性要求。
由图3可知， 3种相关情形下最优解的β⁃ circle与g1、g2均相切，但由于其形状不同，得到的最优解具有较大差异。这说明对于同一组样本，用不同的Copula函数拟合得到的优化结果具有较大的差异，从而说明选择最优Copula的必要性。此外，由于 Gaussian Copula函数可等效为Nataf变换，对比Gaussian Copula 函数和Gumbel Copula函数的优化结果，可知在存在复杂相关性的RBDO问题中，