地震偿付能力资本要求的重新评估：一种开源的时空建模方法

《Annals of Actuarial Science》：Earthquake solvency capital requirements’ re-assessment: an open-source spatio-temporal modeling approach

　　摘要：由自然灾害（如地震）引起的保险风险是联邦监管的财产保险和意外伤害保险公司最低资本测试的关键组成部分。本文提出了一种集成的、开源的、基于模拟的精算框架，用于评估地震保险风险和偿付能力资本要求。该框架结合了时空地震发生建模、基于加拿大地震危险地图的物理信息地面震动估计、建筑物暴露和脆弱性建模，以及在统一流程中的详细保险损失和索赔计算。通过基于核的时空点过程建模来捕捉地震风险的空间异质性，而Voronoi基偏离残差被用作验证模型适用性的局部诊断工具。模拟的保险损失用于估计区域和全国范围内的可能最大损失（PML），并且提出了一种新的资本聚合公式，该公式明确考虑了地震损失跨省的依赖性，与当前的基于区域的监管聚合方法不同。所提出的框架能够在精细的地理尺度上进行空间解析的损失和资本评估，并在一个完全可复制的开源环境中实现。还提供了一个交互式Web应用程序，允许用户模拟在用户指定的震中位置发生的地震及其造成的财务损失和保险索赔。

1. 引言
地震对财产保险和意外伤害（P&C）保险公司及其再保险公司的财务稳定性构成重大风险，尤其是在人口密集地区发生时。从精算角度来看，地震的随机性质——其特征是时间和地点不可预测——给地震相关保险产品的建模和定价带来了很大的不确定性（美国地质调查局，2022a）。与气象灾害不同，灾难性事件的长期回归期和低频率使得历史地震损失经验成为预测未来风险暴露的较差指标。全球每年记录约50万次地震，其中大约100次强烈到足以造成结构损坏（美国地质调查局，2021）。破坏规模在很大程度上取决于建筑基础设施、设计标准和材料。例如，2010年海地的地震造成了近25万人死亡，而同年新西兰发生的同样强烈的地震却没有造成人员死亡，这主要是由于更严格的建筑规范（Wallemacq和House，参考Wallemacq和House2017）。在加拿大，地震频繁发生，但很少强烈到足以造成重大破坏（NRC，2019）。西海岸属于环太平洋地震带，也被称为“火环”，地震活动特别活跃，容易受到高震级事件的影响（美国地质调查局，2022b）。相比之下，加拿大东部的地震活动受到区域应力场的驱动，发生频率较低。然而，北美东部的地质特性使得地震波能够传播更远的距离，即使是从中等强度的地震中也能放大风险（美国地质调查局，2018）。尽管风险增加，但在加拿大的许多地区，尤其是魁北克等东部地区，地震保险的普及率仍然很低。魁北克超过70%的人口居住在地震活跃地区，但在蒙特利尔和魁北克市，只有3.4%的房屋所有者购买了地震保险，而在温哥华和维多利亚这一比例为65%（Swiss Re，2017）。一个主要障碍是地震保险政策的结构，其中保险覆盖通常只在高免赔额之后才生效，这个免赔额通常是以财产保险价值的百分比计算的。这意味着在发生赔付之前，损失必须非常巨大，从而降低了许多保单持有人对保险价值的感知，导致高风险地区的保险不足。加拿大联邦保险和银行监管机构——金融机构监管办公室（OSFI）要求保险公司使用地震模型来估计其可能最大损失（PML），即500年一遇事件的地震损失阈值美元价值，这是最低资本测试（MCT）的核心组成部分。OSFI（2019）规定了以下全国范围内的地震准备金公式：
\[ \text{全国范围} \mathrm{PML}_{1/500} = \left (\text{加拿大东部} \ \mathrm{PML}_{1/500}^{1.5} + \text{加拿大西部} \ \mathrm{PML}_{1/500}^{1.5}\right )^{\frac {1}{1.5}}, \]
其中，先根据保单持有人的免赔额减少单个损失，然后将损失在投资组合中汇总，之后的PML是在应用任何再保险赔偿之前计算得出的。$\mathrm{PML}_{1/500}$对应于500年一遇事件，代表年最大值的第$99.8^{\text{th}}$百分位。对于精算师来说，最大的挑战不仅在于地震的随机性质，还在于将多种数据源结合成一个连贯、可复制的风险评估框架的实际任务。这个问题影响模型校准、情景生成和监管资本估计，特别是当灾难性事件不频繁但潜在严重时。在时空窗口内的地震风险可以表示为一个复合点过程。设$A\subset \mathbb{R}^2$为研究区域，$T$表示一年的时间范围。总损失为
\[ \begin{equation*} S(A,T)=\sum _{i=1}^{N(A,T)} X_i, \end{equation*} \]
其中$N(A,T)$是在$A\times T$区域内发生的显著地震事件的数量，$X_i$表示第$i$次地震事件产生的总保险损失，通过在应用个别保单免赔额后汇总所有受影响暴露的保单持有人的损失得到。第2节给出了“显著”的精确定义。在本文中，我们使用时空点过程（STPP）对$N(A,T)$进行建模，然后通过基于MMI的损害概率矩阵将模拟的震动转换为损失$X_i$。地震引起的保险风险是联邦监管的P&C保险公司MCT的组成部分之一。虽然商业地震模型（如AIR和HazCan）被广泛用于监管和内部资本评估，但它们是专有的、昂贵的，并且不能公开复制（Ulmi等人，参考Ulmi, Wagner, Wojtarowicz, Bancroft, Hastings, Chow, Rivard, Prieto, Journeay, Struik和Nastev2014）。同时，在地震活动的时空建模方面也取得了实质性进展，特别是在北美。区域地震可能性模型项目（Field，参考Field2007）和地震可预测性研究合作组织（Jordan，参考Jordan2006）推进了基于断层的地震预测，而模型评估统计测试（如Number-test（N-test）、Likelihood-test（L-test）和Likelihood-ratio-test（R-test）有助于评估整体模型拟合度（Field，参考Field2007；Jordan，参考Jordan2006；Schorlemmer等人，参考Schorlemmer, Gerstenberger, Wiemer, Jackson和Rhoades2007）。这些商业平台通常结合了统计和基于物理的组件，包括断层破裂建模、波传播和地面运动预测。相比之下，本文的重点是一个透明、完全可复制的、数据驱动的精算框架，用于损失和资本评估，旨在补充而非替代基于物理的地震危险模型。在地震预测中常用的公开记录的统计评估工具，如N-test和L-test，主要评估全球模型拟合度，并对空间局部化的模型缺陷提供有限的见解。在本文中，我们使用空间解析的残差诊断来揭示可能通过总体统计无法检测到的局部模型缺陷。已经提出了基于残差的方法来评估STPP模型的拟合优度。一种常见的方法涉及基于像素的残差方法，其中空间域被划分为规则的网格，并计算每个像素中观察到的事件数量和预期事件数量之间的差异（见Baddeley等人，参考Baddeley, Turner, M?ller和Hazelton2005；Zhuang，参考Zhuang2006）。这些诊断提供了模型拟合的局部评估，并有助于检测低估或高估的区域。然而，当预期强度接近零时，这些方法可能会不稳定，这在地震应用中很常见，因为许多网格单元可能不包含事件。此外，网格分辨率的选择是任意的，可能与数据的空间结构不完全对齐。为了克服这些限制，Bray等人（参考Bray, Wong, Barr和Schoenberg2014）引入了Voronoi残差方法，其中空间根据观察到的事件模式被划分为不规则的数据驱动的多边形。每个Voronoi单元恰好包含一个观察到的点，并包含最接近该点的区域。这种数据驱动的分区适应了局部事件密度，使其特别适用于建模不均匀的过程。正如Bray等人（参考Bray, Wong, Barr和Schoenberg2014）所指出的，标准残差诊断（如Pearson和偏差残差）可以自然地扩展到Voronoi单元，遵循Baddeley等人（参考Baddeley, Turner, M?ller和Hazelton2005）为基于像素的模型制定的通用框架。在这项研究中，我们在新颖的精算背景下应用和实施了Voronoi残差诊断。具体来说，我们在Voronoi单元上应用Pearson和偏差残差来评估加拿大显著地震点模式的时空模型的拟合优度。此外，我们利用这些诊断与端到端的地震保险应用模拟框架相结合。据我们所知，这是首次在基于模拟的地震保险风险建模和资本评估框架中开源、完全可复制的Voronoi基残差诊断的集成。我们还提供了一个交互式Web应用程序，用于在用户指定的位置模拟显著地震的发生和财务影响。最后，我们回顾了OSFI的MCT公式（方程1），并提出了基于我们时空模型的潜在替代方案。本文的其余部分组织如下。第2节描述了分析中使用的地震数据来源。第3节介绍了模拟方法、暴露建模以及加拿大地震风险的另一种MCT方法。第4节报告了模拟结果，并将我们的结果与其他文献中的模拟地震结果进行比较，总结了我们的方法的结果，并将OSFI的地震风险MCT方法与我们的提议进行了比较。第5节总结了本文。除非另有说明，所有货币价值均以加拿大元表示。

2. 数据
本研究使用了由加拿大地质调查局（GSC）在开放文件报告8285中发布的、来自历史和仪器来源的合并的显著地震数据集（Lamontagne等人，参考Lamontagne, Halchuk, Cassidy和Rogers2018）。该数据集涵盖了1600年至2017年的时期，包括172次主震事件。每条记录都包含了关键的地震属性：发生时间、震中坐标（纬度和经度）、深度、震级以及区域分类，还有可用的影响指标，如滑坡、海啸、结构损坏、死亡人数和最大修正麦加利强度（MMI）。这些属性是为后续的地面震动、损害和保险财务损失建模提供了实证基础。该数据集区分了前仪器时代事件（1600年至19世纪末）和仪器时代事件（20世纪起），前者是从档案文件、感觉报告和古地震证据中重建的，后者是通过具有更高空间和时间精度的地震网络记录的。在这项研究中，我们仅关注该目录中的主震子集，排除了余震。这种选择保持了与监管资本建模实践一致的事件为基础的地震风险解释，其中极端的投资组合损失通常由主震驱动。如果将余震作为独立的损失事件包括在内，就有可能重复计算同一暴露组合的损害，因为余震在条件上依赖于主震，并且通常影响相同的地理区域。在所提出的框架中，首先在空间和时间上模拟地震的发生，然后通过第3.1节中描述的地面震动模块确定每个事件的严重程度。因此，模拟的年度损失分布反映了给定年份影响投资组合的最严重破坏事件，无论它在地震学上被分类为主震还是余震序列的一部分，因为对于地震保险背景下的偿付能力资本估计来说，最严重的破坏事件是关注的数量。图1展示了研究期间加拿大显著地震的空间概览。如果事件满足以下至少一个标准，则被认为是显著的：(i) 矩震级超过6，或 (ii) 报告的震动达到MMI V级或更高。选择这些阈值是为了关注具有触发地震保险覆盖的现实潜力的事件，无论是通过高能量释放还是通过足以造成结构损坏的显著震动。附录A中的表8提供了MMI级别的定义供参考。在精算背景下，震级和强度之间的区别特别相关：震级量化了地震源释放的能量，而强度则因位置而异，更好地反映了保险损失的可能性和严重性。因此，MMI是地震保险中用于暴露建模和地理定价的宝贵指标。图1显示了1600至2017年间加拿大发生的重大地震。圆圈的大小和颜色与矩震级成正比。3. 地震模拟方法本节解释了用于通过基于模拟的框架来估计与地震相关的保险损失的方法，该框架结合了地震数据、暴露数据和保险合同数据。同时，还提出了一个替代公式来计算加拿大财产保险公司地震风险的MCT（Measurable Construction Loss）。3.1 计算暴露和模拟地震损失使用以下五步框架来评估研究区域的地震风险：(i) 收集建筑物清单并计算暴露量，(ii) 模拟地震在空间和时间上的发生情况，(iii) 估计地面震动强度，(iv) 估计地震引起的损失率，(v) 估计地震损失和保险索赔支付。分析是在加拿大市的统计单位——普查分区（CSD）层面进行的，加拿大共有5,162个这样的分区（加拿大统计局，2017年）。步骤(ii)至(v)在算法1和算法2中有概述。3.1.1 收集建筑物清单并计算暴露量居民区暴露量是根据2016年加拿大人口普查（加拿大统计局，2016年）中的建筑物类型数量以及加拿大住房统计计划（加拿大统计局，2020a年）中的平均面积估计值来计算的。对于每种建筑物类型，暴露量的计算公式为：\[ \text{暴露量}=\text{建筑物数量} \times \text{平均面积} \times \text{每平方英尺的替换成本} \]。替换成本来自HazCan（Ulmi等人，参考文献Ulmi, Wagner, Wojtarowicz, Bancroft, Hastings, Chow, Rivard, Prieto, Journeay, Struik和Nastev 2014），并使用建筑施工价格指数（BCPI）（加拿大统计局，2020b年）调整至2021年6月的水平。为了考虑区域性和估计的不确定性，我们在替换成本中引入了±10%的随机误差。居民建筑物的内容暴露量假设为结构替换价值的50%（Ulmi等人，参考文献Ulmi, Wagner, Wojtarowicz, Bancroft, Hastings, Chow, Rivard, Prieto, Journeay, Struik和Nastev 2014；FEMA，2013年）。然后根据建筑材料（例如，木材、钢铁、混凝土）对暴露量进行细分，使用HAZUS中的通用建筑方案映射（FEMA（2013年）中的表5.1）。因此，暴露量是按建筑物类别和建筑材料分别计算的，并且根据应用到空间变化震动强度上的特定类型损害函数来生成财务损失。虽然没有明确地模拟单个保险单，但通过建筑物类型、建筑材料、区域面积差异以及替换成本的随机变化，在CSD层面捕捉到了投资组合的异质性。非居民区暴露量是间接估算的，因为目前没有全面的数据库。我们使用建筑许可证数据（加拿大统计局，2020c年），并应用非居民区与居民区许可证之间的年度比率（例如，商业、机构），来推断非居民区的面积。然后类似的计算方法计算这些建筑的暴露量，建筑材料和内容假设遵循HAZUS框架。图2总结了各省的计算暴露总量，为了保护CatIQ（加拿大的损失和暴露指数提供商）的数据隐私，数值被进行了处理。加拿大西部包括不列颠哥伦比亚省（BC）、阿尔伯塔省（AB）、萨斯喀彻温省（SK）、曼尼托巴省（MB）、西北地区（NT）和育空地区（YT），而加拿大东部包括纽芬兰和拉布拉多省（NL）、新斯科舍省（NS）、爱德华王子岛（PE）、新不伦瑞克省（NB）、魁北克省（QC）和努纳武特省（NU）。图2显示了各省的总暴露量（左侧）以及分为东部和西部加拿大（右侧）。垂直虚线分隔了这两个区域。图3提供了每个CSD区域内每平方公里的总暴露量、建筑物和内容的空间地图。附录B提供了关于暴露量估算的更多细节。3.1.2 模拟地震在空间和时间上的发生情况本节简要总结了本文开发的精算地震损失和资本建模框架中使用的STPP（Stochastic Point Process）和残差诊断工具。STPP是一种随机过程，用于模拟事件在空间和时间上的发生情况。其应用包括疾病和流行病的传播，以及地震、海啸和火山爆发等自然灾害。点过程在时间上的研究非常深入；参见Cox和Isham（参考文献Cox and Isham 1980）；Daley和Vere-Jones（参考文献Vere-Jones 2003）；在空间上的研究则参见Cressie（参考文献Cressie 2015）；Moller和Waagepetersen（参考文献Moller and Waagepetersen 2003）；Diggle（参考文献Diggle 2013）。在地震学中，STPP被广泛用于模拟地震模式，特别是余震序列；参见Ogata（参考文献Ogata 1998）；Zhuang等人（参考文献Zhuang, Ogata和Vere-Jones 2002）。STPP定义在子集 \(S\subseteq \mathbb{R}^2\times \mathbb{R}^{+}\) 上，其中每个实现由时空点 \(\lbrace (\boldsymbol{x}_i,t_i) \, : \, i =1, \ldots , n \rbrace\) 组成，其中 \((\boldsymbol{x}_i,t_i) \in A \times T\) 对于某个已知的空间区域 \(A\) 和时间周期 \(T\)，\(\boldsymbol{x}_i\in A\subset \mathbb{R}^{2}\) 表示空间坐标（经度和纬度），\(t_i\in T \subset \mathbb{R}^{+}\) 表示发生时间。一阶属性由其时空强度函数描述：\[ \lambda (\boldsymbol{x},t)=\lim _{\lvert d\boldsymbol{x}\rvert ,\lvert dt\rvert \to 0}\frac {\mathbb{E}\left [ N(d\boldsymbol{x},dt) \right ]}{\lvert d\boldsymbol{x}\rvert \lvert dt\rvert }\]，其中 \(d\boldsymbol{x}\) 表示围绕位置 \(x\) 的一个小空间区域，\(\lvert d\boldsymbol{x}\rvert\) 表示包含时间点 \(t\) 的一个小时间间隔，\(\lvert dt\rvert\) 表示其长度，\(N(d\boldsymbol{x},dt)\) 表示在 \(d\boldsymbol{x}\times dt\) 中的事件数量。因此，\(\lambda (\boldsymbol{x},t)\) 表示单位面积单位时间内的预期事件数量。如果 \(\lambda (\boldsymbol{x},t)=\lambda\) 是常数，则该过程被认为是均匀的。边际强度可以通过对时间或空间进行积分来定义：\[ \lambda _A(\boldsymbol{x})=\int _T\lambda (\boldsymbol{x},t) dt, \quad \text{和} \quad \lambda _X(t)=\int _A\lambda (\boldsymbol{x},t) d\boldsymbol{x} \]。STPP的模型拟合可以使用Schorlemmer等人（参考文献Schorlemmer, Gerstenberger, Wiemer, Jackson和Rhoades 2007）引入的N检验和L检验来评估。N检验将观察到的事件数量与拟合模型的模拟结果进行比较，而L检验评估观察到的对数似然是否超过大多数模拟的对数似然。这些检验评估了整体拟合情况，但对空间或时间的偏差提供的洞察有限。为了提供更细致的诊断，Baddeley等人（参考文献Baddeley, Turner, M?ller和Hazelton 2005）；Zhuang（参考文献Zhuang 2006）引入了基于将域离散化为规则像素网格的残差分析方法。对于一个拟合的强度 \(\hat \lambda (\boldsymbol{x},t)\)，像素 \(B_i \in A\times T\) 的原始残差为：\[ r^R(B_i) = N(B_i) - \int _{B_i}\hat \lambda (\boldsymbol{x},t) dt d\boldsymbol{x} \]，其中 \(N(B_i)\) 是像素 \(B_i\) 中观察到的事件数量，\(B_i\) 是一个可测量的集合，满足 \(B_i\subset A\times T\)，对于 \(i = 1,\ldots ,m\)，其中 \(m\) 是预设的像素总数。这比较了规则网格上均匀分布的像素内观察到的和估计的点总数。原始残差可以缩放以提高可解释性。例如，Pearson残差将原始残差缩放到接近单位方差，类似于线性模型中的Pearson残差，即：\[ r^P(B_i) = \sum _{(t_j,\boldsymbol{x}_j)\in B_i} \dfrac {1}{\sqrt {\hat \lambda (\boldsymbol{x}_j,t_j)}} - \int _{B_i}\sqrt {\hat \lambda (\boldsymbol{x},t)}dt d\boldsymbol{x} \]，对于所有 \(\hat \lambda (\boldsymbol{x}_i,t_i)\gt 0\)。有关缩放网格残差的其他方法，以及Baddeley等人（参考文献Baddeley, Turner, M?ller和Hazelton 2005）和Baddeley等人（参考文献Baddeley, M?ller和Pakes 2008）对其属性的描述，请参见他们的研究。类似于回归模型的偏差，Wong和Schoenberg（参考文献Wong和Schoenberg 2009）提出使用偏差残差来比较两个拟合模型 \(\hat \lambda _1\) 和 \(\hat \lambda _2\) 的对数似然贡献：\[ r^D(B_i) = \left (\sum _{(t_j,\boldsymbol{x}_j)\in B_i} \log (\hat \lambda _1(\boldsymbol{x}_j,t_j)) - \int _{B_i} \hat \lambda _1(\boldsymbol{x},t) dt d\boldsymbol{x} \right ) \quad - \left (\sum _{(t_j,\boldsymbol{x}_j)\in B_i} \log (\hat \lambda _2(\boldsymbol{x}_j,t_j)) - \int _{B_i} \hat \lambda _2(\boldsymbol{x},t) dt d\boldsymbol{x} \right ) \]。尽管基于像素的残差提供了信息，但当每个像素的预期事件数量较小时，它们会变得不稳定，这在地震建模中是一个常见问题。为了稳定估计，Bray等人（参考文献Bray, Wong, Barr和Schoenberg 2014）提出了Voronoi残差分析方法，该方法使用Voronoi镶嵌来划分空间域。每个Voronoi单元 \(C_i\) 恰好包含一个观察到的事件，并包含所有比 \(\boldsymbol{x}_i\) 更靠近的点，对于所有 \(i=1, 2, \ldots , n\)。图4展示了图1中地震位置的Voronoi镶嵌。有关Voronoi镶嵌及其属性的更多细节，参见Boots等人（参考文献Boots, Okabe和Sugihara 1999）。对于一个拟合的强度 \(\widehat \lambda (\boldsymbol{x},t)\)，原始Voronoi残差计算为：\[ r^R_V(C_i) = 1 - \int _{C_i}\hat \lambda (\boldsymbol{x},t) dt d\boldsymbol{x} \]。每个残差衡量了单元 \(C_i\) 中观察到的计数（一个）与预期计数之间的差异。在此基础上，我们应用Pearson Voronoi残差，其缩放方式与基于像素的残差类似，定义为：\[ r^P_V(C_i) = \sum _{(t_j,\boldsymbol{x}_j)\in C_i} \dfrac {1}{\sqrt {\hat \lambda (\boldsymbol{x}_j,t_j)}} - \int _{C_i}\sqrt {\hat \lambda (\boldsymbol{x},t)}dt d\boldsymbol{x} \]。这些残差在 \(\hat \lambda (\boldsymbol{x}_i,t_i)\gt 0\) 的条件下是明确定义的，并且在模型适当的情况下会产生近似标准正态值。为了在局部空间层面上比较两个竞争模型，我们按照Bray等人（参考文献Bray, Wong, Barr和Schoenberg 2014）和Baddeley等人（参考文献Baddeley, Turner, M?ller和Hazelton 2005）的框架，将偏差残差概念应用于Voronoi多边形。类似于线性模型，得到的残差可以称为偏差Voronoi残差，定义为：\[ r^D_V(C_i) = \left ( \sum _{(t_j,\boldsymbol{x}_j)\in C_i} \log (\hat \lambda _1(\boldsymbol{x}_j,t_j)) - \int _{C_i} \hat \lambda _1(\boldsymbol{x},t) dt d\boldsymbol{x} \right ) \]。具有正（负）值的多边形表明 \(\hat \lambda _1\)（\(\hat \lambda _2\)）提供了局部改进的拟合。将这些残差汇总在所有多边形上可以得出一个对数似然比分数，从而在全局上比较竞争模型的适用性。现在，我们将这些残差诊断应用于视觉和空间上评估候选时空模型的适当性，作为更广泛的地震损失和资本模拟框架中的一个验证步骤。为了模拟未来的事件，我们基于1900年至2017年间加拿大发生的137次重大地震（Lamontagne等人，参考文献Lamontagne, Halchuk, Cassidy和Rogers 2018）构建了一个时空点模式数据集。由于历史数据的不完整性可能会扭曲模型的时间组成部分，因此排除了1900年之前的事件。我们假设空间是连续的，时间是离散的（每年的），并使用基于核的方法（Diggle，参考文献Diggle 1985）来估计时空强度函数 \(\hat {\lambda }(\boldsymbol{x}, t)\)。这种选择的动机是核平滑在空间和STPP建模中的稳健性和可解释性得到了广泛验证。我们还采用了一种实际的方法，假设STPP是可分离的，即 \(\lambda (\boldsymbol{x}, t)\) 可以按照（3）分解为：\[ \lambda (\boldsymbol{x}, t)=\lambda _{\boldsymbol{x}}(\boldsymbol{x})\lambda _{T}(t) \]，对于所有 \((\boldsymbol{x},t) \in A\times T\)。这个假设是合理的，因为我们的兴趣在于罕见但影响严重的事件（例如，\(\mathrm{PML}_{1/500}\)），在这些事件中，确切的时间不如风险暴露的空间分布重要。空间组成部分 \(\lambda _{\boldsymbol{x}}(\boldsymbol{x})\) 是Voronoi残差分析的主要焦点。模型拟合、模拟和残差计算是在R语言中使用的splancs、spatstat、deldir和stpp包（Rowlingson和Diggle，参考文献Rowlingson和Diggle 2022；Baddeley和Turner，参考文献Baddeley和Turner 2005；Turner，参考文献Turner 2021；Gabriel等人，参考文献Gabriel, Diggle, Rowlingson和Rodriguez-Cortes 2022；R Core Team，2022），以及为Voronoi基残差诊断开发的自定义函数来进行的。为了估计强度函数的时间组成部分 \(\lambda _{T}(t)\)，我们应用了具有带宽 \(h_T\) 的高斯核密度估计器，带宽 \(h_T\) 是根据Silverman的经验法则（Silverman，参考文献Silverman 1986）选择的，即：\[ h_T=0.9\alpha n^{-1/5} \]，其中 \(\alpha =\min \lbrace \text{样本标准差}, \text{样本四分位数范围}/1.34\rbrace\)。对于空间强度函数 \(\lambda _{\boldsymbol{x}}(\boldsymbol{x})\)，我们使用四次（双权重）核（Diggle，参考文献Diggle 2013）。带宽选择在平衡偏差和方差方面起着关键作用：较小的带宽能够降低偏差，但会导致较高的方差；而较大的带宽虽然会提高估计的平滑度，但会增加偏差。对于地震等非均匀过程而言，不存在普遍最优的带宽，这促使人们使用标准的数据驱动的带宽选择标准，而不是针对特定问题的自定义估计器。我们使用相同的核函数但不同的带宽来拟合两个空间强度模型。

- 模型 $H_1$：带宽 $h_1$ 旨在最小化估计的均方误差 $\hat \lambda _{\boldsymbol{x}}(\boldsymbol{x})$。
- 模型 $H_2$：带宽 $h_2$ 旨在最大化留一法验证得分（leave-one-out likelihood cross-validation score）：
$$\begin{equation*} LCV(h_2)=\sum _l\log (\lambda _{-l}(\boldsymbol{x}_l))-\int _A \lambda (\boldsymbol{u}) \text{d} \boldsymbol{u}, \end{equation*}
$$
其中 $\lambda _{-l}(\boldsymbol{x}_l)$ 是排除该点后，在 $\boldsymbol{x}_l$ 处的强度估计值，$\lambda (\boldsymbol{u})$ 是在域 $A$ 上的平滑估计值。Baddeley 等人（参考文献 Baddeley, Rubak 和 Turner2015）讨论了在点模式表现出强聚类（如地震发生时典型的情况）的情况下使用似然交叉验证的方法。最小化均方误差被作为全局带宽选择的主要标准，而基于留一法验证的带宽提供了一种替代的、数据驱动的平滑度控制方法。随后使用 Voronoi 残差来进行空间诊断验证和比较不同模型。

此外，我们还拟合了一个均匀的泊松过程，记为模型 $P$，其中强度 $\lambda (\boldsymbol{x})=\lambda$ 对所有 $\boldsymbol{x} \in A$ 都是恒定的，作为模型比较的基准。图 5 展示了模型 $P$（a）、$H_1$（b）和 $H_2$（c）的原始 Voronoi 残差。理想情况下，模型的残差应接近零。由于原始 Voronoi 残差的最大值为 1（因为每个单元格包含一个观测到的事件），较大的负值表示模型估计过高。模型 $P$ 显然低估了加拿大东部和西部高地震风险区域的风险，而在低风险区域则高估了风险。鉴于假设的均匀性，这是可以预料的。$H_1$ 和 $H_2$ 的表现都显著更好，大多数区域的残差都接近零。这些视觉模式通过下面讨论的偏差 Voronoi 残差得到了定量验证。附录 C 中的图 14 和图 15 以及图 16–18 中的 Pearson Voronoi 残差图也支持了这个结论。

图 6 比较了模型 $H_1$ 和 $H_2$ 的偏差 Voronoi 残差：加拿大（a）、西部加拿大（b）、东部加拿大（c）。正值表示模型 $H_1$ 在该区域的局部对数似然值高于模型 $H_2$。为了进一步诊断模型拟合的空间差异，我们计算了偏差 Voronoi 残差（见公式（2）），并在图 6 中绘制了模型 $H_1$ 和 $H_2$ 在每个 Voronoi 单元格中的对数似然值差异。模型 $H_1$ 在所有区域的局部拟合都明显优于 $H_2$，特别是在北部和中部加拿大等低风险区域，模型 $H_2$ 有过度预测地震风险的倾向。在高风险区域，包括 Queen Charlotte 断层（西部）和 St. Lawrence 古裂谷（东部），$H_1$ 的表现也略好。最大的偏差残差出现在 Queen Charlotte 断层区域（图 6b）。模型 $H_1$ 和 $H_2$ 之间的总对数似然比值比为 210.67，表明 $H_1$ 提供了系统上更好的局部拟合。与模型 $P$ 相比，$H_1$ 的得分为 586.56，而 $H_2$ 的得分为 375.89，这表明两种非均匀模型都显著改善了相对于均匀基线的拟合效果。

综上所述，均方误差标准、似然交叉验证和基于 Voronoi 的残差诊断方法使我们选择模型 $H_1$ 作为后续地震损失和资本模拟的首选空间强度模型。最终拟合的空间-时间强度函数是通过结合估计的空间和时间组成部分获得的：
$$\begin{equation*}\hat \lambda (\boldsymbol{x}, t)=\hat \lambda _{\boldsymbol{x}}(\boldsymbol{x}) \cdot \hat \lambda _T(t).\end{equation*}$$
然后使用 R 语言中的 stpp 包通过拟合的模型进行地震模拟，从而生成多年内一系列合成性的显著地震事件。

3.1.3 估计地面震动强度
加拿大地震风险定期由 GSC 评估，并在《国家建筑规范》的修订中反映这些评估结果。地面震动强度是地震风险建模的关键输入，通常通过地面运动指数（如峰值地面加速度 PGA、峰值地面速度或谱加速度）的衰减关系来表示。这一组成部分代表了我们数据驱动模拟框架中的物理信息层，允许物理地面运动过程直接影响精算损失和资本计算。在所提出的框架中，地震震级不是从独立的参数化震级分布中模拟出来的，而是通过物理信息丰富的地面震动组成部分隐式地引入的。对于每个模拟的事件位置，PGA 是从拟合到最近网格位置的离散 GSC 地震危险超越点的广义帕累托分布（GPD）中抽样的。然后使用 Wald 等人（参考文献 Wald, Quitoriano, Heaton 和 Kanamori1999）的经验关系将 PGA 转换为 MMI。地面震动随距离震中的衰减通过 Bakun 等人（参考文献 Bakun, Johnston 和 Hopper2003）以及 Bakun 和 Wentworth（参考文献 Bakun 和 Wentworth1997）的衰减关系进行建模，这些关系提供了强度、震级和震源距离之间的明确函数联系。我们使用 2015 年《加拿大国家建筑规范》（NBCC）中为坚实土壤站点提供的 GSC 地震危险值，覆盖了加拿大及其周边地区的网格（Halchuk 等人，参考文献 Halchuk, Adams 和 Allen2015）。对于每个网格点，提供了八个年度超越概率下的平均 PGA 值：0.02、0.01375、0.0100、0.00445、0.0021、0.0010、0.0005 和 0.000404。在每个网格点处拟合一个广义帕累托分布（GPD），以插值不同回报周期下的 PGA 值。GPD 的分位数函数在附录 D 的公式（9）中给出。对于每个模拟的地震震中，我们确定最近的 PGA 网格点，并使用其拟合的 GPD 抽样一个地面运动值。然后使用以下经验关系估计相应的 MMI 级别：
$$\begin{equation*}\text{MMI} = 3.66 \log _{10}(\text{PGA}) -1.66,\end{equation*}
$$
标准误差为 1.08（单位：cm/s2，Wald 等人，参考文献 Wald, Quitoriano, Heaton 和 Kanamori1999）。为了模拟震动的空间分布，我们使用 Bakun 等人（参考文献 Bakun, Johnston 和 Hopper2003）和 Bakun 和 Wentworth（参考文献 Bakun 和 Wentworth1997）开发的衰减方程来估计矩震级作为 MMI 和震中距离 d（单位：km）的函数：
- 加拿大东部（Bakun 等人，参考文献 Bakun, Johnston 和 Hopper2003）：
$$\begin{equation} \text{M} = \frac {1}{1.68}\left (\text{MMI} - 1.41 + 0.00345\text{d} + 2.08\log _{10}(\text{d})\right ),\end{equation}$$
- 加拿大西部（Bakun 和 Wentworth，参考文献 Bakun 和 Wentworth1997）：
$$\begin{equation} \text{M} = \frac {1}{1.09}\left (\text{MMI} - 5.07 + 3.69\log _{10}(\text{d})\right ),\end{equation}$$
为了与显著地震的定义保持一致（Lamontagne 等人，参考文献 Lamontagne, Halchuk, Cassidy 和 Rogers2018），只有那些在将 PGA 转换为 MMI 并应用上述衰减关系后推断出的震级超过 $M\gt 6$ 的模拟事件才会被保留用于损失分析。由于单个地震可以在不同距离产生多个强度值，我们通过反转上述关系来确定每个 MMI 级别对应的距离 d，从而构建等震线图。图 7 展示了一个模拟地震的等震线图样本。同心阴影区域对应不同的 MMI 阈值，说明了强度随距离震中的衰减情况。我们通过引入基于 Bakun 和 Wentworth（参考文献 Bakun 和 Wentworth1997）；Bakun 等人（参考文献 Bakun, Johnston 和 Hopper2003）；Wald 等人（参考文献 Wald, Quitoriano, Heaton 和 Kanamori1999）报告的标准误差的随机噪声来考虑参数不确定性。对于每个 MMI 级别，我们确定了受影响的 CSD 并计算每个 CSD 的面积落在相应等震线圆内的比例。算法 1 概述了模拟地震和生成相应等震线图的程序。

3.1.4 估计地震引起的损失率
地震可能导致财务和非财务损失。本研究关注的是由于建筑物及其内部的物理损坏而导致的直接财务损失。次生灾害（例如滑坡、海啸、火灾）和非财务损失（例如人员伤亡、业务中断）不在本次分析的范围内。结构和非结构损坏分别进行建模：
- 结构损坏指的是建筑物的承重部件（例如墙壁、屋顶）。
- 非结构损坏分为：
- 对位移敏感的部件（例如隔断、装饰层），
- 对加速度敏感的部件（例如机械系统、电梯），
- 建筑物内部物品（例如家具、电子设备）。
为了量化损失，我们依赖于基于 MMI 的损失概率矩阵（DPM），这些矩阵为每种建筑类型分配了在每个 MMI 级别下处于七种定义的损坏状态（例如无损坏、轻微损坏、中度损坏、完全损坏等）之一的概率。DPM 由 Rojahn 等人（参考文献 Rojahn 和 Sharpe1985）开发，并由 Thibert（参考文献 Thibert2008）在假设地面坚实、几何形状规则以及 1990 年之前的建筑规范的基础上扩展到了 BC 地区。表 1 展示了一个 DPM 的示例，展示了木结构住宅建筑在 VI 至 XII MMI 级别下的概率。每个损坏状态都与一个损坏因子范围相关联，该范围表示预期损失的建筑价值比例。

表 1. 木结构住宅建筑的结构损坏 DPM（Thibert，参考文献 Thibert2008）
*** 表示非常小的概率值，几乎为 0。为了将损失量化为暴露比例，我们计算了每个 CSD、建筑类型和损坏组件的平均损坏因子（MDF）。MDF 是在给定 MMI 级别下预期损坏的百分比，它是基于 DPM 中相关概率的损坏因子的加权总和。为了引入模拟变异性，损坏因子是从它们对应的范围内随机抽样的。分别为每种损坏类型（结构（S）、对加速度敏感的（AS）非结构、对位移敏感的（DS）非结构以及建筑物内部物品（BldgC）计算了不同的 MDF。这些 MDF 在第 3.1(v) 节中用于估计保险损失。

3.1.5 估计地震损失和保险理赔支付
对于每个模拟事件，通过将平均损坏因子（MDF）应用于分配给结构（S）、对位移敏感的（DS）、对加速度敏感的（AS）和内部物品（BldgC）组件的暴露值来计算 CSD $j$ 和建筑类别 $k$ 的组件损失。遵循 Onur 等人（参考文献 Onur, Ventura 和 Finn2005）的方法，暴露值被划分为 25% 分配给 S，37.5% 分配给 DS，37.5% 分配给 AS；内部物品的暴露值单独处理。表示为：
$$\begin{equation*} E_{j,k}=\text{building exposure}_{j,k}, \quad C_{j,k}=\text{contents exposure}_{j,k},\end{equation*}$$
在 MMI 级别 $m$ 下的组件损失分别为：
$$\begin{align*} L_{j,k,S}^{(m)} &= \mathrm{MDF}_{S,k}^{(m)} \times 0.25\,E_{j,k},\\[3pt] L_{j,k,DS}^{(m)} &= \mathrm{MDF}_{DS,k}^{(m)} \times 0.375\,E_{j,k},\\[3pt] L_{j,k,AS}^{(m)} &= \mathrm{MDF}_{AS,k}^{(m)} \times 0.375\,E_{j,k},\\[3pt] L_{j,k,BldgC}^{(m)}&= \mathrm{MDF}_{BldgC,k}^{(m)} \times C_{j,k}, \end{align*}$$
总损失为：
$$\begin{equation*} \mathrm{Loss}_{j,k}^{(m)} = L_{j,k,S}^{(m)} + L_{j,k,DS}^{(m)} + L_{j,k,AS}^{(m)} + L_{j,k,BldgC}^{(m)}. \end{equation*$$
当一个 CSD 涵盖多个 MMI 级别时，每个 $L_{j,k,*}^{(m)}$ 都会根据在级别 $m$ 的面积比例进行加权。表 2 显示了用于保险理赔计算的市场渗透率 $\pi _j$、免赔额 $d_j$ 和保险限额 $u_j$。对于每个 CSD：
$$\begin{align*} \mathrm{Claim}_{j,k} = \begin{cases} 0, &\mathrm{Loss}_{j,k}\le d_j,\\[3pt] \pi _j\,(\mathrm{Loss}_{j,k}-d_j), &\ d_j\lt \mathrm{Loss}_{j,k}\le u_j,\\[3pt] \pi _j\,(u_j-d_j), &\ \mathrm{Loss}_{j,k}\gt u_j, \end{cases} \end{align*}$$
其中：
$$\begin{equation*} d_j = \delta _j \sum _{k} E_{j,k}, \quad u_j = \lambda _j \sum _{k} E_{j,k}, \end{equation*$$
$\delta _j$ 和 $\lambda _j$ 分别表示 CSD $j$ 的免赔额和限额百分比。表 2 总结了每个 CSD 的总损失和理赔金额。

附录 E 详细介绍了算法 1 和 2，这些算法用于模拟地震事件和计算损失及理赔。第 3.1(ii)–3.1(v) 节共同构成了一个端到端的、完全可复制的精算流程，该流程在一个模拟框架内 link 了地震发生的时空特性、物理信息丰富的地面震动、按建筑类型划分的暴露值以及保险损失和理赔建模。

3.2 地震风险的最小资本测试
OSFI 当前的监管资本框架在广泛的区域级别（加拿大东部和西部）汇总了地震风险。相比之下，所提出的框架在更细致的空间分辨率（省份和县辖区（CSD）上计算保险损失，并引入了一种资本聚合方法，该方法明确考虑了各地区之间的模拟空间依赖结构。加拿大金融监管局（OSFI）在方程（1）中的MCT聚合是基于加拿大东部和西部每500年发生一次的最大可能损失（PML）来表示的。在精算应用中，PML被定义为年度最大损失分布的高分位数。设$X_1,\ldots ,X_n$为独立同分布（i.i.d.）的年度损失，其分布函数为$F$，并定义$M_n=\max \{X_1,\ldots ,X_n\}$。虽然实践中出现了诸如$(1+\theta )\,\mathbb{E}[M_n]$或$\mathbb{E}[M_n]+\theta \sqrt {\mathrm{Var}(M_n)}$这样的经验法则近似值（参见Wilkinson, Reference Wilkinson1982; Kremer, Reference Kremer1990, Reference Kremer1994），但我们采用了分位数定义：对于小的$\epsilon \gt 0$，有$\mathbb{P}(M_n \le \mathrm{PML}_\epsilon )=1-\epsilon$，这样$\mathrm{PML}_\epsilon$就是超过该阈值的损失，其发生的概率非常低。在超过阈值的峰值框架下，高于高阈值$u$的超越值用GPD $G_{\xi ,\sigma }$来近似，而超越次数$N$近似服从率为$\lambda$的泊松分布，详见附录D。然后，$N$次超越值的最大值近似服从广义极值（GEV）分布$H_{\xi ,\mu ,\psi }$，其中$\mu =\frac {\sigma }{\xi }(\lambda ^\xi -1)$和$\psi =\sigma \lambda ^\xi$（参见Cebrián等人，Reference Cebrián, Denuit and Lambert2003）。解方程（6）得到$\mathrm{PML}_\epsilon = u + \frac {\sigma }{\xi }\left [\left ( -\frac {\lambda }{\ln (1-\epsilon )}\right )^{\xi } -1 \right ]$。受到欧盟委员会2015年提出的偿付能力II框架的启发，我们提出通过反映保险地震损失跨区域依赖性的方差-协方差结构，来聚合所有加拿大省份和地区的PML：$\text{Country-wide} \ \mathrm{PML}_{1/500} = \sqrt {\sum _{r,s}\text{CorrEQ}_{r,s}\times \mathrm{PML}_{1/500,r} \times \mathrm{PML}_{1/500,s}}$，其中$\text{CorrEQ}_{r,s}$表示地区$r$和$s$之间模拟的年度保险损失的相关性，$\mathrm{PML}_{1/500,r}$和$\mathrm{PML}_{1/500,s}$是扣除保单持有人自赔后并应用于任何再保险回收之前的相应地区PML。$\sum _{r,s}$表示对所有地区的求和。对于每个地区$r$，我们通过（a）模拟年度最大值的经验$(1-1/500)$分位数，以及（b）将$\hat {\mu },\hat {\sigma },\hat {\xi },\hat {\lambda }$代入方程（7）得到的EVT估计来计算$\mathrm{PML}_{1/500,r}$。相关矩阵$\text{CorrEQ}$是根据模拟的年度保险损失估计的，包括没有地震损失的年份的零值。单一事件的同步影响自然会在各省之间引起正相关性。我们检查了皮尔逊相关系数、肯德尔$\tau$和相关性（对于零膨胀的连续变量，参见Pimentel等人，Reference Pimentel, Niewiadomska-Bugaj and Wang2015）。

4. 结果
本节展示了三组发现。首先，我们将我们的方法产生的暴露估计与CatIQ的保险暴露指数进行了对比，见第4.1节。其次，我们通过与Onur等人（Reference Onur, Ventura and Finn2005）的比较，对外部验证了模拟的财务损失和保险理赔支付，见第4.2节。第三，我们详细报告了地震模拟的输出，包括分布总结和PML估计，见第4.3节。

4.1 暴露值的比较
CatIQ通过一个订阅平台汇总了加拿大市场的保险公司报告的暴露和损失数据。我们通过将表2中的市场渗透率应用于第3.1节构建的替代成本暴露，将我们的保险暴露估计与CatIQ 2020年的暴露数据进行对比。所有金额均以2021年6月的加元表示。为了尊重CatIQ的保密性，图8中的$y$轴值被屏蔽了。总体而言，我们的公开数据方法产生的住宅暴露量与CatIQ的保险公司报告的值相当，这对于这个开源风险评估框架来说是足够的。对于非住宅暴露，差异较大，这反映了（i）缺乏全面的非住宅建筑国家清单，（ii）我们依赖建筑许可比率将住宅建筑数量转换为非住宅建筑数量，以及（iii）将映射假设应用于HAZUS占用类别；详见第3.1节。图8展示了基于2021年6月加元的计算保险暴露与CatIQ保险公司报告的暴露的对比。

4.2 与以往地震损失研究的比较
我们通过将我们的直接损失估计与Onur等人（Reference Onur, Ventura and Finn2005）的比较来外部验证我们的估计，他们评估了维多利亚市和温哥华市在MMI VIII均匀震动情景下的全市影响。他们的研究报告称，2001年维多利亚市的总货币损失为4.3亿加元，温哥华市为35亿加元。通过BCPI（加拿大统计局，政府，2020b）调整到2021年6月后，这些数字分别为8.35亿加元和68亿加元。在相同的MMI VIII假设下，我们的方法得出的维多利亚市直接损失为9.41亿加元，温哥华市为57.3亿加元。与Onur等人（Reference Onur, Ventura and Finn2005）使用BCPI调整后的数值相比，维多利亚市大约高出12.7%，温哥华市则低15.7%。差异来源于几个方面：（i）Onur等人（Reference Onur, Ventura and Finn2005）没有包括内容损失，而我们的框架明确考虑了内容损失；（ii）建筑成本输入不同：Onur等人（Reference Onur, Ventura and Finn2005）使用的是当地承包商的估计，而我们的方法是基于HAZUS的替代成本并使用BCPI进行了调整。综合考虑这些因素，可以合理解释两者之间的差异方向和幅度，同时在数量级上保持一致。因此，我们将这种比较视为对损失模块的外部合理性检查。

4.3 模拟结果
我们报告了我们的地震模型的模拟输出，并将通过方程（8）提出的聚合方法得到的加拿大全境的PML与OSFI公式（方程（1）得到的结果进行了比较。我们使用第3节中的程序模拟了10万年的地震活动，选择这个时间范围是为了减少蒙特卡洛变异性。

表3总结了至少发生一次重大事件的年份比例以及发生造成损害的地震事件的年份比例。

图9展示了200年模拟事件的一个代表性样本。其空间分布和震级与图1中的历史模式非常相似。图9显示了200年的模拟地震。圆圈的大小和颜色与矩震级成正比。使用第3.1节中的暴露框架和算法1和2中的模拟程序，我们计算了每个模拟事件在县辖区（CSD）层面的直接财务损失和保险理赔支付。图10显示了在CSD发生地震的情况下，损失和理赔的条件期望值，而图11和图12显示了省级的总结。活跃地震带外的CSD的条件期望值为零。在省级层面，条件期望损失随着暴露量的增加而增加，见图2。相比之下，条件期望理赔不仅取决于暴露量，还取决于市场渗透率和合同条款，如免赔额和限额，因此它们的跨省级排名可能与基于暴露量的排名不同。在加拿大东部，较高的非住宅建筑渗透率导致商业/工业类的条件理赔相对较大。最后，东西部之间的差异与衰减特性一致：对于震级为6的事件，在加拿大东部，MMI VI可以延伸到大约200公里，而在加拿大西部则仅约为33公里（参见美国地质调查局，2018年）。

图10基于10万年的模拟年份，显示了每个CSD的平均财务损失（a）和保险理赔（b）。垂直虚线分隔了东部和西部省份。

表4显示了加拿大各省和地区之间模拟财务损失的皮尔逊相关系数。表10、11和12（附录F）中还呈现了肯德尔$\tau$和相关性分析的结果。正如预期的那样，相关性在区域集群内部最强，例如QC-ON和BC与邻近地区之间，这反映了大型事件的同步影响。

表5报告了超过阈值$u$的超越值的GPD和泊松参数估计，其中$u$被选为0.95分位数，当超越值稀少时为0.90，并提供了标准误差。诊断图（未显示）支持GPD的适用性。

表6和表7展示了使用模拟年度最大值的经验分位数（“Simulated”）和方程（7）中的EVT公式（“Estimated”）计算出的各省和地区的PML $_{1/x}$值（以十亿美元为单位）。此外，还根据OSFI的方程（1）和方程（8）以及皮尔逊相关性和肯德尔$\tau$报告了加拿大全境的PML值。模拟和EVT估计的结果紧密一致，证实了模型的适应性。我们观察到，基于相关性的聚合方法得出的值与OSFI的方法相当，而使用肯德尔$\tau$的方法更为保守。

表6和表7基于10万年的模拟地震，展示了加拿大各省和地区之间的模拟财务损失的皮尔逊相关性。图12基于10万年的模拟年份，显示了每个省在发生地震条件下的平均保险理赔。垂直虚线分隔了东部和西部省份。

表5报告了超过阈值$u$的GPD和泊松参数估计，其中$u$被选为0.95分位数，对于超越值稀少的省份则为0.90，并提供了标准误差。诊断图（未显示）支持GPD的适用性。这笔金额代表了相当大的政府支出，例如，它是2020-2024年新冠支持及恢复措施预算的九倍。其他省份也出现了类似的差距，特别是在加拿大东部，那里的住宅保险覆盖率仍然很低。有几个限制因素值得注意。由于缺乏全国性的全面清单，非住宅建筑的暴露风险是通过建筑许可证比例估算的，这给商业和工业损失带来了不确定性。强度模型依赖于可分离性以及完整性不一的历史数据，尽管这适合我们的概率重大损失（PML）研究，但潜在的非平稳性和未明确考虑余震触发因素需要进一步研究。在目前的框架下，损失是基于 annual （年度）数据进行模拟的，因此损失分布反映了给定年份内影响投资组合的最严重事件，无论该事件是被归类为主震还是余震，这对于偿付能力资本的估算来说是相关的。结果还受到政策参数的影响，包括免赔额和限额以及保险覆盖率，因此随着市场条件的变化，输入数据应该进行更新。未来的工作应该纳入更详细的暴露清单，如市政评估记录，扩展危害模型以包括次要风险和场地放大效应，并模拟时间上的非平稳性，并在聚合公式中考虑尾部依赖性。尽管有这些注意事项，所提出的可复制方法仍为加拿大地震风险的定价、资本评估和压力测试提供了一个实用的精算工具箱。

**补充材料**：本文的补充材料可以在以下链接找到：https://doi.org/0.1017/S1748499526100293

**致谢**：这项研究部分得到了精算师协会James C. Hickman学者计划的支持。第二作者感谢加拿大自然科学与工程研究委员会RGPIN-2016-04114项目的支持以及HEC蒙特利尔基金会的支持。第三作者感谢加拿大自然科学与工程研究委员会RGPIN-2015-05447项目的支持。

**利益冲突**：作者声明没有利益冲突。

**附录A. 修改后的Mercalli强度定义**
表8. 修改后的Mercalli强度定义（Wood和Neumann，参考文献Wood和Neumann1931）

**附录B. 建筑清单的收集和暴露量的计算**
我们定义如下：
- **总建筑面积** = 建筑类型单位数量 × 该类型单位的平均建筑面积
- **建筑暴露量** = 总建筑面积 × 每平方米的平均重建成本
- **建筑内容暴露量** = 建筑暴露量 × 内部价值的百分比

接下来，我们描述用于建筑清单和重建成本的数据来源和估算方法。

**住宅建筑**：
按住宅类别划分的计数数据来自加拿大统计局（2016年），并使用Ulmi等人（参考文献Ulmi, Wagner, Wojtarowicz, Bancroft, Hastings, Chow, Rivard, Prieto, Journeay, Struik和Nastev2014）的表C.4转换为HAZUS占用代码。按类别和CSD划分的平均居住面积数据来自加拿大住房统计计划（CHSP，加拿大统计局，2020a年），数据涵盖不列颠哥伦比亚省、安大略省和新斯科舍省。对于缺失的CSD级别的建筑面积，我们使用同一省份中家庭收入中位数最接近的CSD的数据进行估算。如果两者都缺失，则使用该省份的加权平均值，权重为每个CSD中的房屋数量。对于没有CHSP建筑面积的省份，我们使用家庭收入值最接近的CSD的建筑面积进行估算。具体来说：
- 安大略省的建筑面积用于魁北克省、曼尼托巴省和努纳武特省。
- 新不伦瑞克省的建筑面积用于新斯科舍省和佩内布斯科省。
- 不列颠哥伦比亚省的建筑面积用于萨斯喀彻温省、阿尔伯塔省、育空省和纽芬兰与拉布拉多省。

对于没有建筑面积信息的住宅类型，我们假设“其他单栋房屋”的建筑面积与“独立单栋房屋”相同，“复式公寓”的建筑面积与“半独立房屋”相同。

**建筑重建成本**：
建筑重建成本数据来自Ulmi等人（参考文献Ulmi, Wagner, Wojtarowicz, Bancroft, Hastings, Chow, Rivard, Prieto, Journeay, Struik和Nastev2014），最初的数据来自RSMeans，然后使用建筑建造价格指数（BCPI，加拿大统计局，2020b年）将其调整为2021年6月的水平。BCPI是针对选定的普查大都市区（CMA）报告的，因此我们将每个CMA的比率应用到相应的省份。对于没有BCPI的省份/地区，我们采用了观察到的最低比率，即阿尔伯塔省的埃德蒙顿比率。住宅类别的内部价值设定为建筑重建成本的50%（Ulmi等人，参考文献Ulmi, Wagner, Wojtarowicz, Bancroft, Hastings, Chow, Rivard, Prieto, Journeay, Struik和Nastev2014；FEMA，2013）。

**非住宅建筑**：
由于缺乏全国性清单，我们根据加拿大统计局的建筑许可证数据（加拿大统计局，2020c年）对住宅建筑暴露量进行比例调整。我们计算了2011-2019年间机构/政府、商业和工业许可证与住宅许可证的平均年比率，并将这些比率应用于各地区非住宅建筑的总暴露量。表9显示了许可证标签与HAZUS占用代码之间的对应关系。对于非住宅建筑，我们同样使用FEMA（2013）的表5.1将占用类型转换为建筑类型。非住宅类别的内部价值遵循FEMA（2013）中的百分比。

**附录C. 适合的STPP模型的残差分析**
本附录展示了三个强度模型的Voronoi残差诊断结果：均匀泊松模型$P$和两个经过核平滑处理的非均匀模型$H_1$和$H_2$。我们展示了加拿大全国及各地区的原始残差图和Pearson残差图。所有图表都使用相同的颜色刻度以便于视觉比较。

**附录D. 极值理论**
极值理论（EVT）是统计学的一个分支，用于描述分布的尾部行为。主要使用两种方法：（i）块最大值，用于模拟连续块内独立同分布（i.i.d.）观测值的最大值；（ii）超过阈值的峰值，用于模拟超过高阈值$u$的情况。基于渐近性的广义帕累托分布（GPD）用于拟合从地震危险曲线获得的峰值地面加速度（PGA）的极值，以实现极端地面运动模拟。

**附录E. 算法**
- **算法1：地震损失和索赔模拟**
- **算法2：计算损失和索赔**

**附录F. 损失和索赔的相关性**
表10. 基于100,000年模拟地震数据的加拿大各省和地区之间保险索赔的Pearson相关系数
表11.基于100,000年模拟地震数据的加拿大各省和地区间财务损失的肯德尔tau值
表12. 基于100,000年模拟地震数据的加拿大各省和地区间保险理赔的肯德尔tau值

脚注：
1. 如Denuit和Mesfioui（参考文献Denuit和Mesfioui2017）所讨论的，对于零膨胀连续变量，肯德尔tau的值存在一定的范围限制，相关系数不会介于$[-1,1]$之间，而是处于一个更小的范围内。在我们的模拟中，这种限制导致全国范围内的PML（Point-Mass Estimates）结果始终低于皮尔逊相关系数或标准肯德尔tau值，因此我们并未采用这种处理方式。
2. 进行了第二次共计100,000年的模拟后，汇总统计数据几乎没有发生变化。