Title
题目
Accurate and efficient cardiac digital twin from surface ECGs: Insights intoidentifiability of ventricular conduction system
从体表心电图获取准确高效的心脏数字孪生体：关于心室传导系统可识别性的见解
01
文献速递介绍
心脏数字孪生（CDT）技术旨在构建与受试者心脏1:1的计算机虚拟复制品（Niederer等人，2021年）。基于先进的生物物理模拟，心脏数字孪生通过来自观测数据的持续或周期性实时更新（Laubenbacher等人，2024年），精准追踪患者心脏的生理状态。基于所用建模技术的机械特性，心脏数字孪生范式的一个关键要素是内在假设，即虚拟心脏与实体心脏紧密相连，以至于任何刺激或扰动都会产生相同的应急反应，无论这种情况发生在现实空间还是虚拟空间。符合这一假设的心脏数字孪生为个性化医疗提供了巨大的变革潜力，因为它们提供了安全有效的方式来评估患者的心脏功能，以支持诊断、分层和优化治疗方案（Corral-Acero等人，2020年；Cluitmans等人，2024年）。 大规模构建高保真心脏数字孪生这一前景广阔的愿景面临着多项技术挑战。必须实施复杂的建模工作流程，该流程包括两个不同的阶段——称为“解剖孪生”和“功能孪生”，分别致力于从医学图像构建解剖结构精确的心脏模型，以及随后校准模型的大量参数，以同化模型与患者心脏之间的电生理功能。在自动化解剖孪生方面已取得重大进展，能够大规模构建患者心脏解剖结构（Crozier等人，2016年）。功能孪生方法的发展较为滞后，且明显更具挑战性。绝大多数使用标为心脏数字孪生的模型进行的模拟研究都回避了功能孪生。相反，所有患者都统一使用相同的平均参数，并且避免与电描记图（EGMs）或心电图等临床可观测信号进行同类比较，因为这些比较会揭示虚拟模型与物理现实之间的巨大差异（Sung等人，2022年；Bishop和Plank，2025年）。近期开发增强型功能孪生技术的努力有望打破当前模型校准的限制。一项关键进展是引入了更高效的生物物理模型评估方法，如今这些方法能够以实时性能实现器官尺度的心脏电生理（EP）模拟（Pezzuto等人，2017年；Neic等人，2017年）。这些方法已用于概念验证研究，通过采样（Pezzuto等人，2021年；Gillette等人，2021b；Camps等人，2024年；Álvarez Barrientos等人，2025年）和优化方法（Grandits等人，2021a，2024年），将患者心脏模型校准到非侵入性临床心电图数据。因此，功能孪生的核心重大挑战仍然存在——从有限的临床数据（理想情况下是无创获取的数据）中拟合模型的高维参数空间，以稳健、独特且可扩展的方式准确地复现临床记录，达到高度保真的同类比较效果。 一个主要目标是从患者的体表心电图推断心室（心脏的主要泵血腔室）的电活动序列，并在电波形传播的生理精确模型及相关心电图（特别是QRS波群）中复现该序列（Pezzuto等人，2021年；Grandits等人，2021a；Gillette等人，2022年；Li等人，2022年；Grandits等人，2024年；Camps等人，2024年）。心室激活序列及其在体表观察到的QRS波群中的反映是由称为希氏束-浦肯野系统（HPS）的心室传导系统驱动的，该系统启动心室的激活。希氏束-浦肯野系统包括通过房室结接收心房信号的希氏束、分为附着于浦肯野网络的束支的左束支和右束支、贯穿心内膜和心内膜下组织的快速传导的浦肯野纤维网络（Spach等人，1963年；Myerburg等人，1972年；Vigmond和Stuyvers，2016年）。浦肯野网络在末端连接点（称为浦肯野-心肌连接（PMJs））与心室心肌耦合，这些连接点作为最早心室激活的起始部位，因此在确定心室激活模式和体表心电图方面起着关键作用。然而，希氏束-浦肯野系统的结构以及浦肯野-心肌连接的分布具有高度的患者特异性，通常无法通过无创方法获取（Pullan等人，2010年），即使采用侵入性先进标测方法（Palamara等人，2014年），也只能以有限的准确性观测到。由于其作为控制心室激活的关键结构的重要性，以及作为传导障碍治疗目标的重要性（Sharma等人，2018年），揭示其结构和特性的方法是心脏电生理学中最基本的挑战之一。然而，从非常有限、稀疏且嘈杂的数据（如心电图）中推断诸如心室激活序列之类的极高维对象，是一个不适定问题。因此，目前尚不清楚仅通过标准12导联心电图能够在多大程度上且以何种精度推断心室激活序列（Grandits等人，2022年）。因此，必须考虑不同的心室激活图可能在12导联心电图中产生相同QRS波群的可能性（Colli Franzone等人，2014年；Li等人，2025年）。 基于Eikonal和反应-Eikonal模型的快速电生理模拟器为解决推断心室激活序列（理想情况下是唯一的）的可识别性挑战提供了一种方法（Colli Franzone和Guerri，1993年；Pezzuto等人，2017年；Neic等人，2017年）。这些模拟器通过比较模拟心电图与记录的心电图，允许快速探索参数空间以识别候选浦肯野-心肌连接集（Gillette等人，2021b；Camps等人，2024年）。基于采样的方法已成功在局灶性活动（Pezzuto等人，2022年；Meisenzahl等人，2024年）、束支传导阻滞（Pezzuto等人，2021年；Álvarez Barrientos等人，2025年）以及正常窦性心律期间（Gillette等人，2021b；Camps等人，2024年；Cardone-Noott等人，2016年；Barber等人，2021年）从心电图数据中识别心室激活模式。然而，由于需要大量样本，这些方法仅适用于少数参数，或与快速仿真器结合使用时才可行（Pezzuto等人，2022年；Salvador等人，2023年；Camps等人，2024年）。在这里，希氏束-浦肯野系统通常通过少量浦肯野-心肌连接结合快速心内膜层来近似（Pezzuto等人，2021年）。但总体而言，随着参数数量的增加，基于采样的方法和仿真器的扩展性不佳（Larson等人，2019年），因此它们在更一般的情况下推断心室传导系统的能力有限。 在这项工作中，我们利用Geodesic-BP，这是一种从心电图数据中识别浦肯野-心肌连接的快速、基于梯度的方法（Grandits等人，2024年）。（Geodesic-BP的概述见图1。）Geodesic-BP能够在单个高端图形处理器（GPU）上在30分钟内实现对心电图的高精度拟合。请注意，虽然该方法也可以在中央处理器（CPU）上计算，但运行时间会增加到6小时以上。由于是基于梯度的，它随着参数数量的增加而高效扩展，并保证收敛到局部最优解。这使其成为研究可识别性问题的理想平台。一般来说，Geodesic-BP会收敛到损失函数的局部最小值，该损失函数用于衡量记录的心电图与模拟的心电图之间的不匹配。多个最小值的存在表明缺乏可识别性，而这些最小值在参数空间中的分布为潜在的正则化策略提供了见解。 我们通过使用不同的初始参数条件并改变心电图导联的数量（包括多达高密度体表电位图），对单个病例的心电图进行数百次前所未有的高精度拟合来解决这个问题。然后，我们量化识别过程中的不确定性。我们的结果通过使用不同的电生理模拟器在计算机中生成的高保真基准（GT）数据进行定量验证（Gillette等人，2022年）。最后，我们提出了基于生理的浦肯野-心肌连接约束，这显著提高了从12导联心电图识别心室传导系统的可识别性。未使用来自基准模型的确切约束，以阐明希氏束-浦肯野系统的生理假设在实现心电图精确匹配的能力方面可能发挥的作用。 这是第一项证明从标准临床12导联心电图以前所未有的精度推断希氏束-浦肯野系统可行性的研究。我们的方法产生了高度准确的心电图匹配，所识别的心室激活时间的不确定性分布较窄，远低于最先进的侵入性标测技术的测量不确定性。我们相信，我们新颖、稳健且可扩展的优化和正向电生理方法是创建具有可证明保真度的可信心脏数字孪生的关键技术，因此也是将心脏数字孪生转化为临床应用的关键技术。
Abatract
摘要
Digital twins for cardiac electrophysiology are an enabling technology for precision cardiology. Current forwardmodels are advanced enough to simulate the cardiac electric activity under different pathophysiologicalconditions and accurately replicate clinical signals like torso electrocardiograms (ECGs). In this work, weaddress the challenge of matching subject-specific QRS complexes using anatomically accurate, physiologicallygrounded cardiac digital twins. By fitting the initial conditions of a cardiac propagation model, our noninvasive method predicts activation patterns during sinus rhythm. For the first time, we demonstrate thatdistinct activation maps can generate identical surface ECGs. To address this non-uniqueness, we introduce aphysiological prior based on the distribution of Purkinje-muscle junctions. Additionally, we develop a digitaltwin ensemble for probabilistic inference of cardiac activation. Our approach marks a significant advancementin the calibration of cardiac digital twins and enhances their credibility for clinical application.
心脏电生理数字孪生技术是精准心脏病学的支撑性技术。当前的正向模型已足够先进，能够模拟不同病理生理条件下的心脏电活动，并精确复现躯干心电图（ECGs）等临床信号。在本研究中，我们致力于解决一个挑战，即利用解剖结构准确、生理基础扎实的心脏数字孪生体来匹配受试者特异性QRS波群。通过拟合心脏传导模型的初始条件，我们的非侵入性方法可预测窦性心律期间的激动模式。我们首次证实，不同的激动图可生成完全相同的体表心电图。为解决这一非唯一性问题，我们引入了一种基于浦肯野-心肌连接分布的生理先验。此外，我们还开发了一种数字孪生集合，用于心脏激动的概率推断。我们的方法标志着在心脏数字孪生体校准方面取得了重大进展，并提升了其在临床应用中的可信度。
Method
方法
2.1. Electrophysiology model
The forward model for the ECG is based on the eikonal equationwith the lead field method (Neic et al., 2017; Pezzuto et al., 2017;Gillette et al., 2021b). Considering a domain 𝛺 ⊂ R3 representingthe active ventricular myocardium,
2.1. 电生理模型心电图的正向模型基于带有导联场法的程函方程（Neic等人，2017年；Pezzuto等人，2017年；Gillette等人，2021b）。考虑一个代表活跃心室心肌的域\(\Omega \subset \mathbb{R}3\)，
Results
结果
3.1. Identifiability in the general case
From a theoretical standpoint, there is no guarantee that the HPSgoverning ventricular activation can be uniquely identified from thesurface ECG. Optimization algorithms such as Geodesic-BP might converge towards different local minima of the loss function, dependingon the initial positions and timings of the PMJs. Here, we empiricallyquantify the variability of optimized solutions in the unrestricted case,without imposing any physiological constraints on the permissible positions of PMJs, and by randomly sampling the initial guess for the PMJsin the optimization algorithm (see Section 2.4).
Optimization results using Geodesic-BP for fitting the ECG areshown in Fig. 5 (left panel) where the distribution obtained from20 optimized ECGs (𝜇 , 𝜎 ) is overlaid on the GT ECG we aim torecover. The tight envelope around the measured ECG demonstratesthat Geodesic-BP is able to fit the ECG with very high fidelity withoutscaling. The maximum absolute error between optimized and GT ECGwas <2.8 × 10−2 mV (relative error 5.098 930 131 663 843 %) accordingto (17) (average: 2.07 × 10−2 mV/4.066 672 407 563 584 %). Similarly, thePearson correlation coefficient of the ECGs of all samples w.r.t. theGT was >0.994. Convergence is also rapidly achieved with dist fallingbelow 1×10−1 mV within ≤100 iterations for most samples (see also Fig.11).
However, the low errors observed in the ECG did not necessarilycorrespond to small errors in the ventricular activation map 𝜏. Thisis illustrated in Fig. 5 (right panel) where the two cases with themost extreme difference (dist𝜏 ) in the activation map are shown. Thevariability in reconstructed position of PMJs was significant betweenthese two samples, as readily evidenced by fundamentally differentinitial activation sites. For instance, a site of earliest activation in the LVmid-anterior endocardium in one run (top) is a site of latest activationin another run (bottom). Quantitatively, the absolute error in LAT was23.12 ms (average between the samples: 18.32 ms).While the 12-lead ECG was very similar across all samples, a morepronounced variability was observed in the BSPM, as shown in Fig.Again, we show the two BSPMs that were most different accordingto (18) (maximum absolute difference: 0.11 mV), together with theGT BSPM. It is worth noting that these two samples with maximumdistance in BSPM shown in Fig. 6 are not the same samples that yieldeda maximum distance in the activation maps shown in Fig. 5.Overall, the dipole-shaped torso potential was well-captured byall samples, indicating that the overall activation pattern was similaracross all samples. The least correlated areas were found in the sternal region of the BSPM, where reconstructed BSPMs exhibited morecomplex patterns than the GT. Selecting one point in this region showsdifferent extracellular potentials between the two samples 𝜙𝑖 , 𝜙𝑗 andthe GT 𝜙GT (bottom panel of Fig. 6). The correlation between theextracellular potentials at this unseen point was comparatively low,with correlation coefficients against the GT potentials of only 0.69 and0.64, respectively.
3.1. 一般情况下的可识别性 从理论角度而言，无法保证仅通过体表心电图就能唯一确定控制心室激活的希氏束-浦肯野系统（HPS）。像测地线反向传播（Geodesic-BP）这类优化算法，可能会收敛到损失函数的不同局部最小值，这取决于浦肯野-心肌连接（PMJs）的初始位置和时间。在此，我们通过在优化算法中对PMJs的初始猜测值进行随机采样（详见2.4节），在不对PMJs的允许位置施加任何生理约束的无限制情况下，从经验上量化优化解的变异性。 图5（左图）展示了使用测地线反向传播拟合心电图的优化结果，其中20次优化得到的心电图分布（(\varPsi\mu)，(\varPsi\sigma)）叠加在我们旨在恢复的基准（GT）心电图上。围绕实测心电图的紧密包络表明，测地线反向传播能够在不缩放的情况下以极高的保真度拟合心电图。根据式（17），优化后的心电图与基准心电图之间的最大绝对误差小于2.8×10⁻² mV（相对误差5.098930131663843%），平均值为2.07×10⁻² mV/4.066672407563584%。同样，所有样本的心电图与基准心电图的皮尔逊相关系数均大于0.994。收敛速度也很快，对于大多数样本，dist(\varPsi)在≤100次迭代内就降至1×10⁻¹ mV以下（另见图11）。 然而，心电图中观察到的低误差并不一定对应心室激活图(\tau)中的小误差。图5（右图）对此进行了说明，其中展示了激活图中差异最大（dist(\tau)）的两个案例。这两个样本之间，PMJs重建位置的变异性显著，从根本不同的初始激活位点就可明显看出。例如，在一次运行中，左室前中的心内膜是最早激活的位点（上），而在另一次运行中，该位点却是最晚激活的（下）。从数量上看，局部激活时间（LAT）的绝对误差为23.12 ms（样本间平均值为18.32 ms）。 尽管所有样本的12导联心电图都非常相似，但在体表电位图（BSPM）中观察到了更明显的变异性，如图6所示。我们再次展示了根据式（18）得出的差异最大的两个体表电位图（最大绝对差异：0.11 mV）以及基准体表电位图。值得注意的是，图6中所示的这两个体表电位图差异最大的样本，与图5中所示的激活图差异最大的样本并非同一组。 总体而言，所有样本都很好地捕捉到了偶极子形的躯干电位，这表明所有样本的整体激活模式相似。相关性最低的区域出现在体表电位图的胸骨区域，该区域重建的体表电位图比基准图呈现出更复杂的模式。选择该区域的一个点，可以看出两个样本(\phi_i)、(\phi_j)与基准(\phi_{GT})之间的细胞外电位存在差异（图6的下图）。在这个未观测点上，细胞外电位之间的相关性相对较低，与基准电位的相关系数分别仅为0.69和0.64。
Figure
图
Fig. 1. Geodesic-BP: fast definition of a CDT from the surface ECG. Geodesic-BP (Grandits et al., 2024) optimizes the distribution of PMJs until the mismatch between recordedand simulated ECG is minimized. Physiological constraints for the PMJs are automatically imposed during the optimization process.
图1. 测地线反向传播（Geodesic-BP）：通过体表心电图快速构建心脏数字孪生体（CDT）。测地线反向传播（Grandits等人，2024年）对浦肯野-心肌连接（PMJs）的分布进行优化，直至记录的心电图与模拟的心电图之间的不匹配度降至最低。在优化过程中，会自动施加针对浦肯野-心肌连接的生理约束。
Fig. 2. Endocardial surface restrictions limits the domain of possible locations for PMJs. The PMJs are restricted to lie inside the subdomain 𝑆 which is spanned by a band,equidistant to the endocardium 𝑆?
图2. 心内膜表面限制限定了浦肯野-心肌连接（PMJs）可能位置的范围。浦肯野-心肌连接被限定在子域(S)内，该子域由一条与心内膜(S^?)等距的带状区域构成。
Fig. 3. The in silico GT solution generated from the physiologically-detailed HPS. The LAT map 𝜏 is shown together with the HPS in blue (top). The bottom row shows thecorresponding, calculated 12-lead ECG.
图3. 基于生理细节丰富的希氏束-浦肯野系统（HPS）生成的计算机模拟基准（GT）解决方案。上图展示了局部激活时间（LAT）图(\tau)以及蓝色的希氏束-浦肯野系统（HPS）。下图为相应计算得到的12导联心电图。
Fig. 4. Electrode locations of the 12-lead ECG (left) compared to BSPM configurations comprising 32, 64 and 128 electrodes from an anterior (top) and posterior (bottom) view.For each BSPM electrode the uni-polar lead field was computed using the WCT as reference. Electrodes for computing the WCT are shown in red
图4. 左侧为12导联心电图的电极位置，右侧为32导、64导和128导体表电位图（BSPM）的电极配置，分别从正面（上排）和背面（下排）视角展示。对于每个体表电位图电极，以威尔逊中心端（WCT）为参考计算单极导联场。用于计算威尔逊中心端的电极以红色标注。
Fig. 5. Identifiability problems of inverse ECG methods. Left panel: Using Geodesic-BP, the identification of optimal PMJs locations and timings to fit a given ECG (GT) can beachieved with high fidelity. The ECG distribution obtained from 20 optimization runs with different initializations forms a tight envelope around the ECG to be reconstructed(left, 𝜇 ± 𝜎 ). Right panel: The two samples with the largest difference in LAT are shown that reveal significant variability in the solution space which is indicative of limitedidentifiability in the general unrestricted case. We highlight a few of the major differences in LAT (red dashed ellipses). This prompts for quantification of this variability and fora-priori constraints to reduce the non-uniqueness of the solution space.
图5. 心电图逆问题方法的可识别性问题。左图：利用测地线反向传播（Geodesic-BP），能够以高保真度确定拟合特定心电图（基准GT）的最优浦肯野-心肌连接（PMJs）位置和时间。20次不同初始值优化运行得到的心电图分布，在待重建心电图周围形成紧密包络（左图，(\varPsi\mu \pm \varPsi\sigma)）。右图：展示了局部激活时间（LAT）差异最大的两个样本，揭示了解空间存在显著变异性，这表明在一般无约束情况下可识别性有限。我们用红色虚线椭圆突出了局部激活时间的一些主要差异。这促使我们对这种变异性进行量化，并采用先验约束来减少解空间的非唯一性。
Fig. 6. Shown are the two BSPMs of maximum distance among the 20 samples next to the sought-after GT solution at a single instant in time (𝑡 = 60 ms, top), after normalizationwith respect to the mean unipolar potential (see Section 2.6). Differences in potential fields are witnessed in space across the body surface, as well as over time, when comparedat a single unseen electrode location (bottom panel). At all ECG electrode locations used for optimization no differences are apparent, all potential traces are predicted withapproximately zero loss (see Fig. 5, left panel)
图6. 展示了20个样本中差异最大的两个体表电位图（BSPM），以及某一时刻（t = 60 ms，上排）期望得到的基准（GT）解决方案，已通过单极电位均值进行归一化处理（详见2.6节）。对比可见，在体表空间以及时间维度上，电位场存在差异，下图展示了单个未参与优化的电极位置处的差异情况。在所有用于优化的心电图电极位置上，未观察到明显差异，所有电位轨迹的预测误差几乎为零（见图5左图）。
Fig. 7. Preview of the distribution of the 20 samples, reconstructed from 12-lead ECG and BSPM with 128 electrodes as observations, both with and without physiologicalconstraints. We visualize the sample most closely representing the mean 𝜏𝜇 (top 2 rows) along with the deviation between the samples 𝜏𝜎 (bottom 2 rows, see Section 2.9 for theexact definitions). The average deviation over the entire domain, *̄𝜏𝜎 , is indicated at the top right corner of each 𝜏𝜎 model.
图7. 展示了20个样本的分布情况，这些样本是通过12导联心电图和128导体表电位图（BSPM）作为观测数据重建得到的，包含有生理约束和无生理约束两种情况。我们可视化了最接近均值(\tau\mu)的样本（上方2行），以及样本之间的偏差(\tau\sigma)（下方2行，具体定义见2.9节）。每个(\tau\sigma)模型的右上角标注了整个域上的平均偏差(\bar{\tau}\sigma)。
Fig. 8. Comparison between pseudo-bidomain solutions on the torso in terms of errors in potential reconstruction dist𝜙 (18) (left panel) and in the ventricular activation sequencedist𝜏(19) (right panel). Each pair of boxplots shows the errors for unrestricted (blue) and restricted (brown) cases, respectively. High
图8. 躯干伪双域模型解决方案在电位重建误差dist\(\phi\)（式18）（左图）和心室激活序列误差dist\(\tau\)（式19）（右图）方面的比较。每组箱线图分别展示了无约束（蓝色）和有约束（棕色）情况下的误差。高
Fig. 9. Comparison between LATs as reconstructed with an increasing number of torso electrodes, both for ignoring (left panel) and enforcing (right panel) physiological restrictions.For each boxplot, the average dist𝜏 with respect to the GT of the 20 samples is shown as a number at the bottom line
图9. 展示了在忽略（左图）和施加（右图）生理限制两种情况下，随着躯干电极数量增加，局部激活时间（LATs）的重建结果对比。对于每个箱线图，底部标注的数值为20个样本相对于基准（GT）的平均dist\(\tau\)。
Fig. 10. Hyperparameter study of the number of PMJs. The 12-lead ECG was fitted using a variable number of initial PMJs and the error in reconstructing the ECG and theactivation map 𝜏 was measured against the GT model in terms of dist and dist𝜏 . A zoom-in on the results for 𝑁 ≥ 500 is shown in the inset. In the top panel optimization resultsare shown for a number of PMJs of 𝑁 ∈ {10, 50, 300, 5000}, comprising the reconstructed LAT maps and the reconstructed Einthoven II ECG (solid) next to the GT ECG (dashed)
图10. 浦肯野-心肌连接（PMJs）数量的超参数研究。使用不同数量的初始浦肯野-心肌连接拟合12导联心电图，并根据dist(\varPsi)和dist(\tau)来衡量心电图和激活图(\tau)相对于基准（GT）模型的重建误差。插图展示了N ≥ 500时结果的放大图。上图展示了浦肯野-心肌连接数量为N ∈ {10, 50, 300, 5000}时的优化结果，包括重建的局部激活时间（LAT）图以及重建的标准Ⅱ导联心电图（实线）和基准心电图（虚线）。
Fig. 11. Convergence history of Geodesic-BP: We show the evolution of the solution in terms of the ECG on the left 𝑦-axis (logarithmic), over the iterations (first row) andapproximate computational time (second row) on the 𝑥-axis. The top 4 plots show the computed lead II ECG in the current iteration (solid) next to the GT ECG (dashed).Additionally, we show the averaged movement of the initial sites 𝒙𝑖 w.r.t. the last iteration in red (right 𝑦-axis, also in logarithmic scale).
图11. 测地线反向传播（Geodesic-BP）的收敛过程：左侧纵轴（对数刻度）展示了心电图方面的解的演变情况，横轴分别为迭代次数（第一行）和近似计算时间（第二行）。上方4幅图展示了当前迭代中计算得到的Ⅱ导联心电图（实线）与基准（GT）心电图（虚线）。此外，红色曲线（右侧纵轴，同样为对数刻度）展示了初始位点(\boldsymbol{x}_i)相对于最后一次迭代的平均移动距离。
