射血分数保留型心力衰竭的集总参数与有限元建模

摘要

本研究基于集总参数法和有限元分析，提出了两种射血分数保留型心力衰竭的计算模型。这些模型用于评估压力超负荷和心室顺应性降低所引起的心室及关联血管系统的血流动力学变化。

摘要

在心血管疾病计算建模领域，科研工作主要集中在射血分数降低型心力衰竭（HFrEF），而广泛忽视了射血分数保留型心力衰竭（HFpEF），而后者近年来已成为全球最主要的心力衰竭类型。鉴于目前HFpEF的in silico模型 representation 极为匮乏，本文提出了两种不同的计算模型，用于模拟由左心室压力超负荷引起的HFpEF血流动力学变化。首先，采用数值求解器开发了一个面向对象的集总参数模型。该模型基于零维（0D）Windkessel类网络，依赖于构成元件的几何和力学特性，具有计算成本低的优势。其次，利用有限元分析（FEA）软件包实现了多维仿真。FEA模型结合了三维（3D）电-机械心脏响应、结构形变以及基于流体腔室的血流动力学的多物理场模型，并采用简化的集总参数模型来定义不同流体腔室之间的流量交换特征。通过这两种方法，均成功模拟了压力超负荷引起的左心室及近端血管系统的急性和慢性血流动力学变化。具体而言，压力超负荷通过减小主动脉瓣口面积进行建模，而慢性重构则通过降低左心室壁的顺应性来模拟。与HFpEF的科学和临床文献一致，两个模型的结果均显示：（i）左心室与主动脉之间的跨主动脉压力梯度急性升高，同时每搏输出量减少；（ii）左心室舒张末期容积慢性下降，提示存在舒张功能障碍。最后，FEA模型表明，在整个心动周期中，HFpEF心肌组织的应力显著高于健康心肌组织。

引言

心力衰竭是全球主要的死亡原因之一，当心脏无法充分泵血或充盈以满足机体的代谢需求时即可发生。临床上使用射血分数（即每次收缩时从左心室射出的血液占储存血量的相对比例）将心力衰竭分为两类：（i）射血分数降低的心力衰竭（HFrEF）和（ii）射血分数保留的心力衰竭（HFpEF），分别对应射血分数低于或高于45%的情况^1,2,3。HFpEF的症状常由左心室压力超负荷引发，后者可由多种疾病导致，包括主动脉瓣狭窄、高血压和左心室流出道梗阻^3,4,5,6,7。压力超负荷引发一系列分子和细胞水平的异常，导致左心室壁增厚（向心性重构），最终引起心室壁僵硬度增加或顺应性下降^8,9,10。这些生物力学改变显著影响心血管血流动力学，表现为舒张末期压力-容积关系升高以及舒张末期容积减少¹¹。

心血管系统的计算建模推动了人们对生理和病理状态下血压与血流的理解，并促进了诊断和治疗策略的发展¹²。计算机 模型可分为低维或高维模型，前者利用解析方法评估整体血流动力学特性，计算需求较低；后者则在二维或三维域内提供更广泛的多尺度、多物理场心血管力学与血流动力学描述¹³。集总参数Windkessel模型是低维描述中最常见的形式。该模型基于电路类比（欧姆定律），通过电阻、电容和电感元件的组合来模拟心血管系统的整体血流动力学行为¹⁴。本课题组最近提出了一种在液压域中的替代性Windkessel模型，相较于传统的电路类比模型，该方法能更直观地模拟大血管、心腔及瓣膜几何结构和力学特性的变化。该仿真基于一种面向对象的数值求解器（见材料表）开发，能够捕捉正常血流动力学、心肺耦合的生理效应、单心室生理状态下呼吸驱动的血流，以及主动脉缩窄引起的血流动力学改变。该模型在集总参数模型的基础上进一步拓展，提供了一种物理意义上更直观的方法，可用于模拟包括心力衰竭在内的多种病理状态¹⁵。

高维模型基于有限元分析（FEA）来计算时空血流动力学和流固相互作用。这些模型能够提供局部血流场的详细且精确的描述；然而，由于其计算效率较低，不适用于整个心血管系统的模拟研究^16,17。本研究采用了一款软件包（见材料表），作为具有解剖学精度的成人四腔心有限元分析平台，该平台集成了电-机械响应、结构形变以及基于流体腔室的血流动力学。改进后的人类心脏模型还包含一个简化的集总参数模型，用于定义不同流体腔室之间的流量交换，以及对心肌组织的完整力学表征^18,19。

已有多种集总参数模型和有限元分析（FEA）模型被提出，用于表征心力衰竭的血流动力学异常并评估治疗策略，特别是在射血分数降低型心力衰竭（HFrEF）中应用机械循环辅助装置的背景下^{20,21,22,23,24}。因此，一系列复杂程度不同的0D集总参数模型已成功通过优化二元件或三元件的电模拟Windkessel系统，再现了人体心脏在生理状态及HFrEF状态下的血流动力学特征^20,21,23,24。这些模型大多采用时变弹性公式构建单心室或双心室模型，以模拟心脏的收缩功能，并利用非线性的舒张末期压力-容积关系来描述心室充盈过程^25,26,27。一些综合性模型能够表征复杂的心血管网络，并模拟心房与心室的泵血功能，已被用作器械测试的平台。然而，尽管关于HFrEF领域的文献已相当丰富，针对射血分数保留型心力衰竭（HFpEF）的计算机模型却极少被提出^{20,22,28,29,30,31}。

由 Burkhoff 等人³² 和 Granegger 等人²⁸ 近期开发的一种射血分数保留型心力衰竭（HFpEF）血流动力学的低维模型，能够捕捉四腔心的压力-容积（PV）环，全面再现 HFpEF 各种表型的血流动力学特征。此外，他们利用其计算机模拟平台评估了针对 HFpEF 的机械循环辅助装置的可行性，开创了将计算方法应用于 HFpEF 生理研究及器械开发的先河。然而，这些模型仍无法捕捉疾病进展过程中观察到的血流和压力的动态变化。Kadry 等人³⁰ 的一项最新研究通过在低维模型上调节心肌的主动舒张功能和左心室的被动刚度，实现了对不同舒张功能障碍表型的模拟。该研究基于心肌的主动与被动特性，提供了对舒张功能障碍的全面血流动力学分析。类似地，高维模型的研究主要集中在射血分数降低型心力衰竭（HFrEF）^{19,33,34,35,36,37}。Bakir 等人³³ 提出了一种完全耦合的心脏流体-电-机械有限元分析（FEA）模型，用于预测 HFrEF 的血流动力学特征以及左心室辅助装置（LVAD）的治疗效果。该双心室（或两腔）模型采用耦合的 Windkessel 电路来模拟健康心脏、HFrEF 以及接受 LVAD 支持的 HFrEF 的血流动力学状态^33,37。

同样，Sack 等人³⁵构建了一个双心室模型以研究右心室功能障碍。该双心室几何结构来源于患者磁共振成像（MRI）数据，并通过图像分割构建模型的有限元网格，用于分析心室辅助装置（VAD）支持下衰竭右心室的血流动力学³⁵。已有研究开发了四腔室有限元分析（FEA）心脏模型，以提高对心脏电机械行为建模的准确性^19,34。与双心室描述相比，基于MRI的四腔室心脏模型能更准确地反映心血管解剖结构¹⁸。本研究采用的心脏模型是四腔室有限元分析模型的一个成熟范例。与集总参数模型和双心室有限元模型不同，该模型能够捕捉疾病进展过程中发生的血流动力学变化^34,37。例如，Genet 等人³⁴利用同一平台实现了对射血分数降低型心力衰竭（HFrEF）和射血分数保留型心力衰竭（HFpEF）中观察到的心肌重塑过程的数值生长建模。然而，这些模型仅评估心肌肥厚对结构力学的影响，未能全面描述相关的血流动力学变化。

为了解决射血分数保留型心力衰竭（HFpEF）在计算机模拟中的缺失问题 本研究中采用的模型，即该研究团队先前开发的集总参数模型¹⁵ 并对有限元分析模型进行调整以模拟射血分数保留型心力衰竭（HFpEF）的血流动力学特征。为此，首先将展示各模型在基线状态下模拟心血管血流动力学的能力。随后评估由狭窄引起的左心室压力超负荷，以及由于心脏重构导致的左心室顺应性降低——这是HFpEF的典型特征——所产生的影响。

方案

1. 0D 集总参数模型

1. 模拟设置
   注意：在数值求解器环境中（参见材料表），构建如图1所示的计算域。该域由四腔心脏、上半身、腹部、下半身和胸腔 compartments，以及近端血管结构组成，包括主动脉、肺动脉和上、下腔静脉。本模拟中使用的基本元件属于默认液压元件库的一部分。详细信息请参见补充文件。
   1. 浏览液压元件库，找到所需元件：液压管道、恒容液压腔室、线性阻力、离心泵、单向阀、可变截面孔口以及自定义液压流体。
      1. 将液压管道元件拖入工作区。
         ​注意：这些元件用于模拟血管和心腔中的摩擦损失、管壁顺应性以及流体可压缩性。通过该模块，采用达西–威斯巴赫公式计算压力损失；而由于管壁顺应性引起的直径变化则取决于顺应性比例常数、管腔内压力和时间常数。最后，流体的可压缩性由介质的体积模量定义。
      2. 插入恒容液压腔室元件，以定义管壁顺应性和流体可压缩性。
         ​注意：该模块不考虑因摩擦引起的压力损失。
      3. 添加线性阻力元件，以定义流动阻力。
         注意：该阻力与血管几何特性无关，类似于电路类比中的Windkessel模型所使用的电阻元件。还应插入其他模块，如离心泵、单向阀、可变截面孔口和自定义液压流体元件，以向系统提供所需的压力输入，模拟心脏瓣膜对血流的影响，并定义血液的力学特性。通过这些元件，可全面表征心血管系统在生理和病理状态下的行为。离心泵的输入信号见图S1A。
      4. 通过自定义的可变顺应性腔室元件，对每个心腔的收缩性进行建模。
         注意：该元件接受顺应性作为随时间变化的用户自定义输入信号，其基础为时变弹性模型（图S1B-D）。
   2. 为每个元件提供相应参数，参数值如表S1所示，亦见于Rosalia等人的研究¹⁵。
   3. 对于需要时间变化的用户自定义输入信号的模块（左心室泵、可变顺应性元件和可变截面孔口模块），分别插入物理信号（PS）重复序列模块。
      注意：本模拟所使用的输入信号见图S1。
   4. 选择默认的ODE 23t 隐式求解器，并运行模拟100秒，以达到稳态。

2. 有限元分析模型

1. 模拟设置
   注意：该有限元分析（FEA）模型采用顺序耦合的电-机械分析方法。在此模型中，首先进行电学分析，然后将所得的电势作为后续力学分析的激励源。因此，模拟设置包含两个工作域：电学域（ELEC）和力学域（MECH），这两个域已在有限元分析软件中预定义（材料表）¹⁸。因此，下文仅描述分析流程。该有限元模型在电学和力学材料建模中使用了以下用户子程序：HETVAL, VUANISOHYPER 和 UAMP¹⁸。
   1. 进入 ELEC 域，使用 Standard 模块中的预定义温度过程进行电学分析。
      1. 使用一个名为 BEAT 的单步分析。将心动周期持续时间设为 500 ms，并对代表窦房结（SA 节点）的节点集施加电势脉冲（节点集：R_Atrium-1.SA_NODE）。
      2. 查看默认电波形，其在 200 ms 内从 -80 mV 平滑上升至 20 mV，采用平滑步幅定义，如模型指南所述¹⁸。在电学分析中使用默认的材料常数以调整房室延迟（AV delay）。
      3. 启动 Job 模块，并创建一个名为 heart-elec 的作业。
   2. 完成电学分析设置后，进入 MECH 域，执行基于流体腔的力学分析。
      注意：力学模拟在电学分析之后进行，电学分析所得的电势将作为力学分析的激励源。力学分析包含多个步骤。
      1. 使用三个主要步骤：PRE-LOAD, BEAT1 和 RECOVERY1。在 PRE-LOAD 步骤中，检查心脏预应力状态的边界条件。使用 0.3 s 作为步长时间，线性增加流体腔室内的压力。
         注意：预定义的流体腔压力值见 表 S3。心脏的预应力状态已在正常心脏模拟设置中定义，初始节点条件包含在外部模拟文件中，详见 表 S5。每当修改边界条件时，均需通过逆向力学模拟重新计算零应力状态，具体步骤见 3.2.2–3.2.4。
      2. 在 BEAT1 步骤中，使用 0.5 s 作为步长时间以模拟收缩过程。
      3. 在 RECOVERY1 步骤中，选择 0.5 s 用于心脏舒张和心室充盈，对应心率为 60 bpm。
      4. 启用后续步骤 BEATX 和 RECOVERYX，以模拟多个心动周期，直至达到稳态。
         注意：三个心动周期足以达到稳态。在 24 核处理器（3.2 GHz × 24）上，一个周期的模拟耗时约 8 小时。
      5. 启动 Job 模块，并创建一个名为 heart-mech 的作业，启用 双精度 选项。
2. 检查简化的集中参数 Windkessel 模型
   注意：该有限元模型的力学域包含一个血流模型，该模型基于简化的集中参数电路，由基于表面的流体腔和流体交换组合而成，如 图 2¹⁸ 所示。
   1. 使用上述注释中提到的 Windkessel 表示法运行模拟。检查血流模型的表示方式，以调整阻性元件和容性元件的数值，分别对应流动阻力和结构顺应性。
   2. 检查四个心腔的三维有限元表示，并确保其几何位置准确。
   3. 检查心脏装配体，并切换到 Interaction 模块，调整四个心腔各自的顺应性和收缩性参数。
      注意：Interaction 模块中的默认值已配置为模拟理想化的健康人的心动周期¹⁸。
   4. 在 Interaction 模块中检查以下静水压流体腔：CAV-AORTA、CAV-LA、CAV-LV、CAV-PULMONARY_TRUNK、CAV-RA、CAV-RV、CAV-SVC、CAV-ARTERIAL-COMP、CAV-PULMONARY-COMP 和 CAV-VENOUS-COMP（表 S3）。
   5. 将顺应性腔（CAV-ARTERIAL-COMP、CAV-PULMONARY-COMP 和 CAV-VENOUS-COMP）作为立方体积使用，因其分别代表动脉、静脉和肺循环的顺应性。
   6. 将三个顺应性立方体积连接至接地弹簧，并检查其刚度值，以模拟动脉、静脉和肺循环中的压力-容积响应。
   7. 检查以下静水压流体腔之间的流体交换定义：动脉-静脉、静脉-右心房、右心房-右心室、右心室-肺系统、肺系统-左心房、左心房-左心室、左心室-主动脉（表 S4）。
   8. 调整黏性阻力系数，以修改每个流体交换连接中的血流模型（详见 补充文件 中关于黏性阻力效应的说明）。
3. 多物理场模拟
   1. 在工作目录中定位 CAE 数据库 文件。
      注意：本方案中的有限元模型以数据库形式提供，文件名为 LH-Human-Model-Beta-V2_1.cae。
   2. 将输入文件、对象文件和库文件插入工作目录以运行模拟。完整输入文件和库文件列表见 表 S5。
   3. 启动有限元分析模型模拟软件（参见 材料表）。
      注意：如需了解与后续版本的兼容性，请咨询软件供应商¹⁸。
   4. 按照第 2.2 和 2.3 节所述，在 ELEC 和 MECH 域中检查部件、装配体和边界条件。
   5. 首先运行名为 heart-elec 的电学模拟作业，具体操作见第 2.1.1.3 节。目视检查电势结果，确认 heart-elec 模拟按预期运行。然后确保结果文件 heart-elec.odb 存在于工作目录中。
   6. 切换至 MECH 域，进入第二阶段模拟。检查力学模拟中所用材料常数的数值，以建模所需的被动和主动心脏响应。
   7. 确保用于力学分析的材料库文件使用 HYBRID- 字符串命名。如需修改心腔的材料响应，可调整相应的混合材料文件，或在 CAE 模块的 Materials 部分定义新的材料行为以完全替换原有材料响应。
      注意：有关内置本构关系的详细信息，请参阅用户指南¹⁸。
   8. 在 PRE-LOAD 步骤中，设置静水压腔的压力，以获得所需的生理行为。使用内置的平滑幅值选项，从零逐步升至目标压力水平，具体操作见步骤 2.1.2.1。
   9. 禁用步骤 2.1.2.1 中定义的压力边界条件，以在循环系统内保持总体血容量恒定的条件下运行血流模型。运行名为 heart-mech 的模拟作业，具体操作见第 2.1.2.5 节。

3. 主动脉瓣狭窄

注意：主动脉瓣狭窄常常是射血分数保留型心力衰竭（HFpEF）的驱动因素，因为它会导致左心室壁的压力超负荷，最终引起向心性重构及顺应性丧失。主动脉瓣狭窄中观察到的血流动力学变化常逐步发展为主动脉瓣狭窄中所见的HFpEF样改变。

1. 集总参数模型
   1. 修改位于左心室腔室内的主动脉瓣处的 PS 重复序列元件中的输入信号。模拟与基线相比孔口面积减少 70% 的情况（表 S6）。
      注：输入值将表示每个心动周期中狭窄瓣膜的孔口面积。可通过将主动脉瓣 PS 元件的 start output values 向量乘以一个小数来轻松调整孔口面积值，该小数表示相对于原始值的最终孔口面积。在本研究中，使用 0.3 的系数以实现 70% 的狭窄。
2. 有限元分析（FEA）模型
   1. 修改 LINK-LV-ARTERIAL 参数的流体交换定义。
      注：该参数具有一个针对左心室与主动脉之间血流调校的黏性阻力系数。可通过修改有效交换面积来调节血流，从而建立适当的主动脉瓣狭窄模型（表 S7）。
   2. 找到工具箱文件夹，并将该文件夹内的文件复制到主工作目录中。
   3. 通过执行工具箱文件¹⁸进行逆向力学模拟。为此，在流体腔中将左心室和左心房的吸力压力调整为 6 mmHg，以调节其在主动脉瓣狭窄模型中的初始容积状态。运行 inversePreliminary.py 函数。
      ​注：每当边界条件发生改变时，均需通过逆向力学模拟重新计算无应力状态。
   4. 完成逆向力学模拟后，运行后处理函数：calcNodeCoords.py 和 straight_mv_chordae.py。对其他流动参数使用默认值，并按照第 2.1.2.5 节所述运行新的力学模拟。

4. 射血分数保留型心力衰竭（HFpEF）血流动力学

注意：为了模拟慢性重塑的影响，对左心的心脏力学特性进行了修改。

1. 集总参数模型
   1. 修改左心室顺应性元件的左心室舒张期顺应性，以模拟压力超负荷引起的室壁僵硬，使用表 S8中的舒张末期顺应性值。
      注：假设顺应性从收缩末期到舒张末期呈线性下降。
   2. 将左心室泵的泄漏阻力增加至 18 × 10⁶ Pa s m^-3（表 S8），以反映射血分数保留型心力衰竭（HFpEF）中观察到的升高的左心室压力。
2. 有限元分析（FEA）模型
   1. 编辑左心室几何结构的主动材料属性。增加刚度分量，以调节构成模型中影响纤维方向和层面方向应力分量的主动组织响应。
      1. 修改mech-mat-LV_ACTIVE文件中左心室的材料响应。
         ​注：可调节左心室腔的刚度大小，以实现适当的舒张期顺应性效应。
      2. 在各向异性超弹性本构关系中，增加刚度参数a和b，以模拟射血分数保留型心力衰竭（HFpEF）生理状态下的刚度增加。
      3. 在预加载（PRE-LOAD）步骤中，将左心室和左心房的流体空腔压力设置为 20 mmHg。
      4. 执行逆向力学模拟，以获得左心室和左心房的容积状态。从heart-mech-inverse.odb文件¹⁸中导出节点坐标。
      5. 执行后处理函数：calcNodeCoords.py 和straight_mv_chordae.py，操作方法如步骤 3.2.4 所述。在工作目录中定位生成的新节点输入文件，并按照第 2.1.2.5 节所述执行新的力学模拟。

结果

基线模拟的结果如图3所示。该图展示了左心室和主动脉的压力与容积波形（图3A），以及左心室的压力-容积环（图3B）。两种计算机模型显示出相似的主动脉和左心室血流动力学特征，均处于生理范围之内。在心室排空和充盈阶段，两种平台预测结果之间存在细微差异，其中有限元分析（FEA）模型相较于集总参数模型能够更好地捕捉非线性效应。在生理条件下，此类非线性效应主要源于心肌组织的超弹性响应，因此多域高阶计算模型能更准确地再现这些现象¹⁸。

对于主动脉瓣狭窄，获取了心室和主动脉的血流动力学参数，因为这种情况常导致左心室压力超负荷，并最终发展为射血分数保留的心力衰竭（HFpEF）。在主动脉瓣口面积减少70%的情况下，两种模型的压力和容积波形如图4所示。狭窄导致主动脉瓣跨瓣压差升高。在本研究考虑的70%狭窄程度下，集总参数...

讨论

本研究提出的集总参数模型和有限元分析（FEA）平台在生理条件下重现了心血管血流动力学，涵盖了狭窄引起的压力超负荷急性期以及慢性射血分数保留型心力衰竭（HFpEF）阶段。通过捕捉压力超负荷在HFpEF急性期和慢性期发展过程中的作用，这些模型的结果与HFpEF的临床文献一致，包括由于主动脉瓣狭窄导致的跨主动脉压力梯度出现、左心室压力升高，以及因心室壁僵硬引起的舒张末期容积减少⁴¹。此外，该FEA模型能够在整个心动周期中捕捉到HFpEF心肌应力的升高。为确保这些模拟的正确设置，必须严格遵循上述方案部分中列出的步骤。对于集总参数模型，必须准确重建液压元件网络，如图1所示，并将规定的数值作为输入参数输入（表S1和表S2）。此外，求解器模块必须被定义，并连接至网络中的任意节点。

FEA 模型的运行需要与求解器一起打包的所有模拟文件，这些文件列在表 S5中

材料

本文使用的材料清单

姓名	公司	目录编号	评论
Abaqus 软件	Dassault Systèmes Simulia Corp.		使用版本：2018；有限元分析仿真软件
HETVAL	Dassault Systèmes Simulia Corp.		使用版本：2018
液压（等温）库	MathWorks		使用版本：2020a
Living Heart 人体模型	Dassault Systèmes Simulia Corp.		使用版本：V2_1，解剖结构精确的成人四腔心有限元分析平台
MATLAB	MathWorks		使用版本：2020a，面向对象的数值求解器
SIMSCAPE FLUIDS	MathWorks
UAMP	Dassault Systèmes Simulia Corp.		使用版本：2018
VUANISOHYPER	Dassault Systèmes Simulia Corp.		使用版本：2018