基于因子分析的地震动特征提取及潜在破坏势评估

胡进军，刘巴黎，谢礼立

(1. 中国地震局工程力学研究所，黑龙江，哈尔滨 150080；2. 中国地震局地震工程与工程振动重点实验室，黑龙江，哈尔滨 150080)

地震动破坏势反映了结构响应和地震动的综合特征，是衡量地震动破坏强度的重要物理量。地震动的随机性以及结构响应的复杂性使得评估及量化地震动破坏势成为工程抗震领域的热点和难点。自从1994 年Northridge 地震以后，基于性能的地震工程(Performance-based earthquake engineering, PBEE)和基于性能抗震设计(Performancebased seismic design, PBSD)方法得到土木工程和地震工程领域众多学者的广泛共识，并在结构抗震设计和研究中得到广泛应用[1]。在PBEE 框架中，地震动强度指标(Intensity measure,IM)是联系地震危险性和结构地震响应的桥梁。合理的选择IM不仅能够减小结构损伤预测的不确定性[2]，也有助于全面准确评估地震动的潜在破坏势。

近年来，国内外学者针对单一的IM与结构损伤指标(Damage measure,DM)开展了一系列研究。RIDDELL[3]将IM分为3 类，即与加速度相关、与速度相关以及与位移相关，针对单自由度(Single degree of freedom, SDOF)体系，对三类IM和4 个DM之间进行相关性分析，结果表明：没有一个IM能在全部周期段内能与DM指标保持良好的相关性。ELENAS[4]及ELENAS 和MESKOURIS[5]研究了地震动谱指标、峰值指标及能量指标和结构整体DM的相关性，从而评估了这些IM的潜在破坏势。韩建平等[6]选取18 组汶川地震动记录，针对双线性SDOF 体系对12 个IM和3 个DM之间的相关性进行了研究。叶列平等[7]针对SDOF和多自由度(Multi-degree of freedom, MDOF)系统研究了32 个IM和不同响应指标的相关性，分析了不同IM的优缺点和适用范围。李爽等[8]针对近场和远场地震动，对IM和结构整体DM的相关性进行了研究，评估了近场和远场地震动的潜在破坏势。李雪红等[9]针对SDOF 体系和减震桥梁结构，通过相关性分析研究了不同IM的敏感性。杨参天等[10]基于高层隔震结构，评估了25 个IM与高层隔震结构DM的相关性。张艺欣等[11]针对高层钢筋混凝土结构提出适用于高层建筑结构的考虑多个周期的均值谱加速度指标，从而评估了地震动对高层建筑结构的潜在破坏势。卢啸等[12]从相关性和离散度等角度研究了超高层建筑结构IM与DM变化规律以及不同IM的适用性。张成明等[13]针对圆形隧道，对20 个常用IM在表征地震动破坏势的适用性进行了评估。

尽管上述关于单一IM与DM之间相互关系研究取得了一定的进展，但是单一IM在反映地震动特征和特性方面始终显得不足。针对单一IM在评估地震动破坏势时的片面性，有些学者提出了标量型IM和向量型IM。标量型IM通常由多个单一的IM组成[14-16]，而向量型IM通过引入多个自振周期的谱加速度而构成[17]，或者采用谱形参数来预测结构的概率地震需求[18]。此外，还有些学者通过考虑多个IM来构造复合型IM。邱意坤等[19]提出了一种综合考虑周期延长和高阶振型效应的复合型IM，基于高耸结构在三维地震动作用下对该指标的充分性和有效性进行了检验。LIU 等[20]基于偏最小二乘回归提出一种复合型IM的构造方法，分别从充分性、有效性、有益性和鲁棒性4 个角度对构造的复合型IM开展统计性分析。还有一些学者提出了基于人工神经网络和遗传算法[21-27]和支持向量机[28-29]的复合IM构造方法。然而，人工神经网络的一个潜在缺点是耗时。与人工神经网络相比，支持向量机对核函数的类型有较强的依赖性。

近年来，一些学者引入多元统计学中的主成分分析[30]和典型相关分析[31-32]，将多个IM线性组合成为一个新的IM，即：地震动多元强度指标(Multivariate intensity measure,MIM)，并对MIM与DM之间的相关性进行全面研究。多元统计学有一个重要的共同点就是多变量数据处理，这使得基于多元统计方法构造MIM成为可能。如果能够从更加广泛的角度应用多元统计分析方法来构造MIM，一方面可以通过引入多个IM使得构造的MIM与DM保持较好的相关性，另一方面可以使组成MIM的每一个IM分量都有着较好的可解释性。此外，在评估地震动潜在破坏势过程中，影响同一DM的同一类IM之间通常都具有较强的相关性，它们是表征地震动潜在破坏势且相互关联的IM。如图1 所示，与加速度相关的IM地震动峰值加速度PGA和有效峰值加速度EPA之间有很强的相关性。这些构造MIM的基本IM在表征地震动潜在破坏势时贡献不一，甚至有些贡献是相似或相同的，这样就容易使IM产生相关重叠性，而因子分析(Factor analysis, FA)正是解决上述问题的有效方法。它在损失最少信息的前提下，将原始变量(IM)综合成较少的互不相关的因子变量(同类IM)，消除了不同因子变量间的相关重叠性。同时，还可以对不同因子变量进行解释，使得FA在MIM的构建过程中具有更加广泛的分析意义。

图1 PGA 与EPA 散点图Fig. 1 Scatter plots of PGA and EPA

鉴于此，本文引入多元统计学中的FA 进行地震动特征提取，通过对公共因子进行旋转使得旋转后的公共因子具备较好的可解释性，同时将保留的互不相关的主因子(Principal factor,PF)以方差贡献率作为权重系数进行线性组合，构造地震动多元强度主因子。以SDOF 和MDOF 体系为研究对象，选取2766 组地震动记录为输入，同时选取9 个基本IM及不同类型的结构响应指标为DM，对构造的地震动多元强度主因子与不同DM进行相关性分析，以期为地震动潜在破坏势的准确评估提供新的思路。

1 因子分析

1.1 基本原理

为了消除代表不同类IM因子变量间的相关重叠性，本文引入统计学中常用的FA 方法，该方法利用降维思想，在损失最少信息的前提下用少数几个互不相关的变量来解释原始变量[33]。假设有p个成分的观测随机向量X=[X1,X2, ···,Xp]，有均值μ和协方差矩阵Σ。因子模型要求X是线性依赖于几个不能观测的称之为公共因子的随机变量F=[F1,F2, ···,Fm]和p个附加的称之为特殊因子的变差源ε=[ε1, ε2, ···, εp]。FA 模型为：

式中：lij为因子载荷系数；L为因子载荷矩阵；F为公共因子；ε 为特殊因子。由于模型中lij、F以 及ε 为 不 能 观 测 量，不 能 直 接 从X1,X2,···,Xp的观测值来直接确认此因子模型。为了得到因子模型，需要对随机向量F和ε 作附加假设：

式中，Sij为随机变量Xi的方差。

保留少数几个公共因子的目的是以少数几个PF来尽可能充分的反映原样本的信息。公共因子的保留通常采用两大标准：① 累计方差贡献率CVC(Cumulative variance contribution, CVC)准则(CVC 大于85%)；② 特征值准则(特征值大于1)。CVC 表达式如下：

通过上述方法求得的公共因子的解释意义有可能不明确。本文采用最大方差旋转法进行因子旋转，其目的是使得每个因子上的载荷一部分趋于0 而另一部分±1。

为了反映原始变量的相关关系且同时能够描述样本的特征，需要用公共因子来代表原始变量，公共因子表示为原始变量线性组合的方程称为因子得分方程：

式中，βjp为因子得分系数。将原始变量取值和因子得分系数代入式(9)便可求得因子得分。

1.2 数据前处理和计算

由于随机变量X1,X2, ···,Xp的量纲往往有所差别，为剔除量纲带来的影响，需要对原始数据进行0-1 标准化。随机变量标准化后会使每个变量Xi中的数值平均变为0、标准差变为1，此方法被广泛的使用在机器学习算法(例如：支持向量机、逻辑回归和类神经网络)和多元统计分析(例如：主成分分析和典型相关分析)中。此外，FA 模型中因子载荷矩阵L、公共因子F、特殊因子ε 以及协方差矩阵Σ 等的求解参考文献[33]。

1.3 因子分析逻辑框图

以上介绍了FA 的基本思想和理论方法，为了更加清楚地展示FA 各分析步骤的脉络关系，更好地运用FA 进行地震动特征提取，本文给出了FA 的逻辑框图，如图2 所示。

图2 因子分析逻辑框图Fig. 2 Logical structure diagram of factor analysis

2 地震动记录及强度指标的选取

2.1 地震动记录

本文从美国太平洋地震工程研究中心(PEER NGA-West2)强震数据库中挑选2766 条水平地震动记录，挑选原则简述如下：① 自由场地震动记录；② 矩震级大于5.0；③ 震中距介于10 km～60 km；④PGA大于50 cm/s2；⑤ 不包含脉冲型地震动。选取的地震动记录的震级-震中距分布图如图3 所示。可见，选取的地震动记录的震级和震中距范围分布较广，能较为充分地体现各因素对地震动特性的影响。参考文献[34]场地类别的划分标准，根据30 m 覆盖层平均剪切波速Vs30将挑选的2766 条水平地震动记录分为B、C、D、E 4 类。地震动记录详细分类信息如表1 所示。

表1 本文选取的地震动记录按照场地类别分类结果Table 1 Classification of ground motion records according to site class in this study

图3 震级-震中距分布图Fig. 3 Magnitude-fault distance distribution of the selected ground motions

2.2 地震动强度指标

本文首选取了23 个代表性的IM[3]作为初选指标，可将IM分为与加速度相关的参数、与速度相关的参数及与位移相关的参数3 类。由于23 个IM之间有些表达式得相近，需要对这些IM进行有效筛选。有效筛选的原则主要有：1) 两类IM在表达式上相差一个指数的情况，通常只选择其一；2) 某一类IM的表达式在对数坐标系下是另一类或者几类IM(表达式较为简单)的线性组合，通常只选择表达式较为简单的IM。

通过上述有效筛选后最终选取9 个IM为基本指标，具体信息如表2 所示。其中，u¨g(t) 、u˙g(t)和ug(t)分别为地震动加速度、速度和位移时程；Td和tf分别为地震动显著持时和总持时；PSV为伪速度谱值；IA为Arias 强度，t5和t95分别表示5%和95%的Arias 强度所对应的时间，v0为地震动加速度时程曲线单位时间内通过时间轴的次数；arms为Pa的平方根。值得注意的是，文献[3]将修正的阿里亚斯烈度PD划分为与与速度相关的参数，本文通过分别计算2766 条地震动记录对应的PD与PGA、PGV及PGD的相关性发现，PD与PGA、PGV及PGD的相关系数分别为0.694、0.830 及0.708，PD与PGV的相关性略好于PD与PGV和PGD的相关性。这说明，在相关性方面，PD与加速度相关参数、速度相关参数及位移相关参数都保持了相近的表现。2766 条地震动记录对应的PD与PGA、PGV及PGD的散点图及相关性如图4 所示。

图4 PD 与PGA、PGV 及PGD 散点图Fig. 4 Scatter plots of PD and PGA, PGV and PGD, respectively

表2 本文选取的基本地震动强度指标Table 2 The selected basic ground motion intensity measures

3 结构建模及损伤指标

3.1 单自由度体系

SDOF 体系建模以及非线性时程分析时，SDOF体系恢复力模型分别选择理想弹塑性模型(EP)、双线性模型(BP)以及修正的Clough 模型(MC)。上述三种恢复力模型分别如图5 所示。其中，BP 模型和MC 的屈服后刚度系数α 取0.05。采用Newmark-β 法求解SDOF 体系运动方程。积分采用文献[35]建议的步长，即：地震动记录采样时间间隔，自振周期T的1/25 以及0.01 s 三者的最小值。SDOF 体系T取为0.1 s～6 s，周期间隔为0.1 s。强度折减系数R取值为2、3、4 和5，其中R的表达式为：

图5 单自由度体系恢复力模型Fig. 5 Hysteretic models of SDOF systems

式中：Fe为地震作用下结构保持弹性所需的最小强度；Fy为结构屈服强度。

3.2 多自由度体系

选择OpenSees 平台进行MDOF 体系的建模和非线性时程分析。计算模型为3、6、8、10、12、15、20、30、40 和50 层的集中质量悬臂剪切模型(串杆模型)，模型自振周期分别为0.35 s、0.64 s、0.83 s、1.02 s、1.23 s、1.52 s、2.10 s、3.00 s、4.10 s、5.01 s。这种简化的MDOF 计算模型，是由在各层集中质量的悬臂剪切模型，来进行计算模拟在地震作用下可能在柱间形成塑性铰的多层多跨框架结构。计算模型的基本滞回规则采用文献[36]提出的双折线模型。模型参数信息：每层质量m为600 t，每层的层间弹性剪切刚度k为1×109N/m，每层的层高h为3.0 m。层间屈服位移角dy为1/150，强化段刚度系数αs为0.05，软化段刚度系数αc为-0.05，名义延性系数为6。模型参数如图6 所示。

图6 层剪切模型参数Fig. 6 Lumped mass shear model parameters

3.3 结构损伤指标

结构损伤指标DM是反映结构地震响应以及和表征结构地震损伤的参数。针对SDOF 体系，采用文献[37]中所使用的DM。其中，最大位移响应umax为表征结构和非结构构件损伤的参数；最大加速度amax响应为表征非结构构件损伤的参数；等效滞回耗能速度响应veq为表征累积损伤的参数。针对MDOF 体系，分别采用最大层间位移角(Maximum interstory drift,ISD)、最大楼层加速度(Maximum floor acceleratiom,Amax)以及最大楼层速度(Maximum floor velocity,Vmax)作为DM。滞回耗能Eh及等效滞回耗能速度veq的表达式如下：

式中：f(t)为恢复力时程；u为加速度时程；k为结构刚度；m为结构质量。

通过计算选取地震动记录的IM构造数据矩阵，建立FA 模型并计算公共因子，采用最大方差旋转法进行因子旋转来获取具备较好可解释性的PF。PF可以表示为原始变量IM的线性组合，即：因子得分方程。因子得分方程中的的系数可以通过加权最小二乘法回归获得。将保留的PF以方差贡献率作为权重系数进行线性组合，构造地震动多元强度主因子。

4 地震动特征提取

4.1 多指标数据矩阵构建

为了采用FA 方法进行地震动特征提取，本文首先定义选取的9 个IM作为基本指标，即原始变量为lnIM=[lnPGA, lnEa, lnPa, lnPGV, lnEv, lnSI,lnPD, lnPGD, lnEd]T。

记地震动数据矩阵为X，包括m个记录、p个地震动强度指标以及S个场地类别。则数据矩阵可表示为：

将数据矩阵X按照场地类别展开，得到地震动记录和IM的二维矩阵XS，可表示为：

4.2 因子分析及主因子的选取原则

将数据矩阵按照场地类别展开，对各场地类别的二维数据矩阵分别进行FA，公共因子和因子载荷矩阵的求解采用主成分法[33]。同时，采用CVC 准则和特征值准则相结合的方法(即：同时满足CVC 大于85%和特征值大于1)对公共因子进行保留。由图7 可知：通过对B、C、D 和E 四类场地的二维数据矩阵分别进行FA 得到的前两个公共因子的CVC 超过85%，且前两个公共因子的特征值都大于1，其他公共因子的特征值都小于1.0。这说明选取前两个公共因子作为PF，即主因子PF1和主因子PF2，能较为全面地解释IM的方差。

图7 不同场地数据因子分析对应的公共因子特征值及累积贡献率Fig. 7 Eigenvalues and cumulative contribution rates corresponding to different site categories

4.3 因子旋转

为了使保留的PF具有较好的可解释性，需要对因子模型进行旋转。本文采用正交旋转法[33]中的最大方差旋转法对因子模型进行旋转，并求取旋转后的因子载荷系数FL(Factor loading,FL)。FL反映了原始变量对PF 的依赖程度，FL越大说明该原始变量与PF的相关系数越大。图8 给出了B、C、D 和E 四类场地FA 对应的FL。由图可知，Ed和PGD在PF1上的载荷系数超过0.9，Ev在PF1上的载荷系数超过0.8，SI和PGV在PF1上的载荷系数超过0.7，其他IM在PF1上的载荷系数都小于0.7。这说明第一主因子PF1主要由与速度相关的IM(PGV、Ev及SI)及与位移相关的IM(PGD、Ed)决定，因此主因子PF1可以解释为速度和位移联合作用因子；Pa和PGA在PF2上的载荷系数超过0.85，Ea在PF2上的载荷系数超过0.7，其他IM在PF2上的载荷系数都小于0.7。这说明第二主因子PF2主要由与加速度相关的IM(PGA、Ea及Pa)3 个IM决定，因此主因子PF2可以解释为加速度作用因子。值得注意得是，修正的阿里亚斯烈度PD在PF2上的载荷系数大于0.7，同时在PF1上的载荷系数大于0.6，这与2.2 节中PD与加速度相关参数PGA、速度相关参数PGV及位移相关参数PGD相关系数较为接近的结论一致。

图8 不同场地数据因子分析对应的因子载荷系数Fig. 8 Factor loading corresponding to different site categories

综合上述可知，所有基本IM在不同场地类别下都稳定地属于某个PF，即各场地类别各指标所属的PF趋于稳定。各场地类别FA 中，9 个基本IM都转化为2 个PF，根据PF与基本IM的关系，将其分别命名为速度和位移联合作用因子及加速度作用因子。

通过加权最小二乘法[33]回归获得因子得分系数，将PF表示成原始变量IM的线性方程，并将两个PF以方差贡献率占比作为权重系数进行线性组合，构造地震动多元强度主因子。其中，权重系数ω的计算如式(15)所示。根据FA 的原理可知，保留的两个PF之间互不相关。采用FA 既消除了PF之间的相关重叠性，又使得保留的PF具备较好的可解释性(速度和位移联合作用因子及加速度作用因子)。构造的地震动多元强度主因子lnIMFA的表达如式(16)所示。

4.4 地震动多元强度主因子的构造

5 SDOF 体系地震动多元强度主因子与结构损伤指标相关性分析

5.1 相关系数

本文选取Pearson 相关系数[30,38-39]来考察IM和DM之间的相关性，相关系数越大，则IM表征地震动破坏势的能力越强。Pearson 相关系数的表达式如下式所示：

5.2 基本地震动强度指标的潜在破坏势

为了衡量IM与DM之间的相关性，同时为了量化IM表征地震动潜在破坏势的能力，本文以基本地震动强度指标PGA、PGV及PGD为例，SDOF体系恢复力模型为EP 模型，且强度折减系数R=4 的条件下，给出了单一IM 分别与amax、veq及umax之间的Pearson 相关系数分析结果，如图9 所示。由图可知，PGA与不同DM之间的相关系数随T增大而逐渐减小，PGD与不同DM之间的相关系数随T而逐渐增大，而PGV与不同DM之间的相关系数随T增大先增大而后缓慢减小。PGA、PGV及PGD与不同DM之间的相关系数分别在短周期(0.0 s～0.5 s)、中等周期(0.5 s～3 s)以及长周期段(T>3 s)内取得峰值，这表明采用单一IM来表征地震动潜在破坏势时需要考虑结构周期的影响。

图9 SDOF 体系基本地震动强度指标与结构损伤指标的相关系数Fig. 9 Correlation coefficient between the basic intensity measure and damage measures of SDOF systems

5.3 地震动多元强度主因子的潜在破坏势

以强度折减系数R=4 的EP 模型、BP 模型以及MC 模型分析结果为例，对IMFA与DM进行相关性分析，如图10 所示。由图可知，针对同一恢复力模型，IMFA与不同DM之间的相关系数在T<1.0 s 时缓慢上升，当T>1.0 s 时相关系数保持不变。不同DM与IMFA的相关系数除在T<1.0 s 时有较小差别外，其他周期范围内差别很小。

图10 SDOF 体系地震动多元强度主因子与结构损伤指标的相关系数Fig. 10 Correlation coefficient between principal factor of multivariate earthquake intensity measures and damage measures of SDOF systems

通过图9 和图10 的对比可以发现，IMFA既保持了与速度相关及与位移相关的单一IM在中长周期段内表征地震动潜在破坏势的良好表现，又较大改善了与加速度相关的单一IM在短周期内表征地震动潜在破坏势的不足。主要原因在于，IMFA考虑了更多的地震动信息，从而保持了其与不同DM之间稳定且较好的相关性。

为了研究场地类别和恢复力模型参数对IMFA和DM之间相关性的影响，图11(a)给出了EP 模型R=4 对应的B、C、D 和E 类场地的相关性分析结果，图11(b)给出了2766 条地震动记录EP 模型对应的R=2、3、4 和5 的相关性分析结果，图11(c)给出了2766 条地震动记录R=4 对应的EP 模型、BP 模型和MC 模型的相关性分析结果，同时以DM为veq为例。由图可知，针对不同场地类别，相关系数在整个周期段内差别很小；针对强度折减系数R，相关系数在T<1.0 s 时有较小差别，T>1.0 s 时差别很小；针对恢复力模型，相关性分析结果在整个周期段内差别很小。总体而言，场地类别和恢复力模型参数对IMFA和DM之间相关性影响较小。

图11 SDOF 体系不同参数对相关系数的影响Fig. 11 Influence of different parameters on correlation coefficient between principal factor of multivariate earthquake intensity measures and damage measures of SDOF systems

6 MDOF 体系地震动多元强度主因子与结构损伤指标相关性分析

图12 为MDOF 层剪切模型相关性分析结果。以基本地震动强度指标PGA、PGV及PGD为例，给出了这三个IM分别与Amax、Vmax及ISD之间的Pearson 相关系数分析结果。由图可知，PGV与体系不同DM之间的相关系数随T增大先增大而后减小。PGD与不同DM之间的相关系数随T增大而逐渐增大。MDOF 体系PGV和PGD与不同DM之间的相关系数的变化趋势与对应SDOF 体系相关系数的变化趋势一致。PGA与表征非结构构件损伤的指标Amax以及表征结构累积损伤的指标Vmax之间的相关系数随T增大而先增大后趋于平稳，这一点与SDOF 体系对应PGA与相应DM之间相关系数的变化趋势不同。MDOF体系IMFA与不同DM之间的相关系数在T<1.0 s时上升，当T>1.0 s 时相关系数保持不变，这与SDOF 体系IMFA与相应DM之间相关系数的变化趋势一致。

综合对比图12(a)～图12(c)可以发现，单一IM与DM之间的相关系数受结构周期的影响较大，且在相同周期区间内与不同DM之间的相关系数差别也较大。单一IM的适用性无法在各个周期段内以及针对不同DM达到最佳[7]。IMFA和DM之间的相关系数除在短周期段略小于0.9 外，在其他周期范围内都超过0.9 且保持稳定状态。说明IMFA通过考虑更多的地震动信息，保持了其与不同DM之间稳定且良好的相关性。

图12 MDOF 体系地震动强度指标与结构损伤指标的相关系数Fig. 12 Correlation coefficient between intensity measure and damage measures of MDOF systems

7 结论

本文基于统计学中的因子分析方法进行地震动特征提取并构造地震动多元强度主因子，针对SDOF 和MDOF 体系，以2766 条地震动记录为输入，通过分析地震动多元强度主因子和损伤指标之间相关性，以此来评估构造的地震动多元强度主因子表征地震动潜在破坏势的能力。主要结论如下：

(1) 所选基本地震动强度指标在不同场地类别下都稳定地属于某个主因子，即不同场地类别下各地震动强度指标所属的主因子趋于稳定。各场地类别因子分析中，9 个基本地震动强度指标都转化为2 个主因子。

(2) 基于因子分析提取的互不相关的两个主因子，根据主因子和基本地震动强度指标的关系，将其分别命名为速度和位移联合作用因子和加速度作用因子。第一主因子主要由与速度相关以及与位移相关的基本地震动强度指标决定，第二主因子主要由与加速度相关的基本地震动强度指标决定。基于因子分析提取的地震动强度主因子具备较好的可解释性。

(3) 针对SDOF 体系，场地类别、强度折减系数以及恢复力模型对地震动多元强度主因子与结构损伤指标相关性的影响较小。

(4) 针对SDOF 体系和MDOF 体系，构造的地震动多元强度主因子既保持了与速度相关及与位移相关的基本地震动强度指标在中长周期段内表征地震动潜在破坏势时的良好表现，又能较大改善与加速度相关的基本地震动强度指标在短周期内表征地震动潜在破坏势的不足。构造的地震动多元强度主因子通过考虑更多的地震动信息，保持了其与不同结构损伤指标之间稳定且良好的相关性。