《International Journal of Microwave and Wireless Technologies》：An effective parametric approach for patient-specific head model generation in quantitative microwave imaging
医学微波成像（MMWI）正作为一种前景广阔的非侵入性诊断模式兴起，为传统成像方法提供了安全、便携且经济的替代方案。尽管潜力巨大，MMWI由于其底层逆散射问题的病态和非线性特性而面临重大挑战。目前主要采用两种策略：定性成像（提供快速但近似的重建）和定量成像（旨在恢复组织的完整介电常数分布）。两种方法均显著受益于精确的患者特异性背景模型，这些模型能提高可靠性并降低计算复杂性。本工作基于在EuCAP 2025会议上报告的初步成果，推进了用于定量脑成像的患者特异性头部模型的开发。该研究采用了贴合解剖结构的基函数，并将模型扩展至包含脑脊液（CSF），这是一种具有显著不同介电特性、散射强烈的物质。所提出的框架提高了重建精度和收敛效率，从而为推动定量MMWI走向实际临床应用做出了贡献。
医学微波成像（MMWI）作为一种非侵入性诊断技术，因其低成本和使用非电离辐射而受到越来越多的关注。传统的医学成像模式，包括磁共振成像（MRI）、计算机断层扫描（CT）和正电子发射断层扫描（PET），仍然是临床诊断的金标准，提供高分辨率图像。然而，它们的广泛应用受限于高昂的运营成本、庞大的设备，以及在CT和PET情况下可能暴露于有害电离辐射。相比之下，MMWI提供了一种便携、安全且经济的替代方案，减轻了MRI和CT的许多缺点，尽管空间分辨率较低。在其各种潜在应用中，脑部成像（特别是中风或肿瘤的检测和随访）已成为MMWI发展最有前景的领域之一，突显了其在快速床边诊断评估方面的潜力。
实现MMWI需要求解一个病态且高度非线性的逆散射问题，这带来了巨大的计算挑战。为此，文献中提出了两种主要策略。第一种策略涉及将逆问题公式化为线性化版本，以重建被检查区域的定性实时图像。这种方法主要对检测小异常物有效，前提是有先前关于健康组织的患者特异性信息。然而，在缺乏此类信息的情况下，定性MMWI重建往往不可靠。
第二种策略，称为定量成像，旨在重建整个感兴趣区域内的复介电常数分布，包括中风或肿瘤的分布。由于生物组织表现出不同的介电常数值，这种方法产生反映组织分布的类扫描图像。定量成像涉及求解完整的非线性逆散射问题，通常通过迭代技术，如非精确牛顿方法或失真玻恩迭代法（DBIM）。然而，由于其迭代性质，这些方法通常计算成本过高。更关键的是，它们容易产生假解，即满足测量数据但偏离真实介电常数分布的估计。过高的计算成本和假重建的风险都可以通过在优化过程中采用患者特异性背景模型作为初始估计来缓解。
基于上述讨论，拥有合适的患者特异性模型对于提高定性和定量逆散射方法的性能都至关重要。在定性成像中，模型必须足够精确，以实现对细微组织变化（如肿瘤生长或疾病进展）的快速检测。在定量成像中，目标是捕捉更精细的组织变化，模型作为可靠的初始化，引导优化过程朝向真实解，加速收敛，并显著减少所需迭代次数，使这些算法更适用于实时诊断应用。
患者特异性背景模型的概念最初在EuCAP 2025会议上提出，并发表在会议论文集中，报告了初步结果。在该工作中，使用基于玻恩近似、并以多项式基函数表示对比度的线性逆散射模型重建了初始头部图像。这些图像随后被分割并转换为四区域体模，并采用适当公式化的DBIM方法估计每个区域的介电常数。从MMWI图像生成体模的技术在后续研究中进一步完善和扩展，并进行了可能受益于患者特异性背景模型的中风检测分析。
在本工作中，研究人员引入了相比前述工作的几项关键改进。主要进步在于使用贴合解剖结构的基函数（而非多项式基函数）来建模未知的对比度。具体而言，采用低阶球谐函数的组合来捕捉相对于外部头部几何形状的角度偏差，而径向脉冲函数则用于表示沿径向方向快速变化的介电常数。这种公式化能够更精确地描绘组织边界，提高重建精度，并减少所需基函数的数量，从而缓解了问题的病态性。
此外，与之前的工作相比，考虑了一个更真实且复杂的体模，该体模包含了一个对应于脑脊液（CSF）的额外组织层。最后，通过在散射参数中引入不同信噪比（SNR）水平下的噪声，评估了所提出方法在实际条件下的鲁棒性。
因此，本工作的主要成果是一个改进的四区域体模以及相应的组织介电常数估计值。所提出的体模可以作为健康头部的快速、患者特异性近似，并可作为临床工作流程中的初始建模步骤进行集成，即使仅有有限的先验信息（如外部头部形状），也能实现快速中风检测。这种能力在时间紧迫的场景中尤其有价值，因为那时传统成像模式可能无法轻易获得。
论文其余部分组织如下：“介电常数分布的参数化建模”部分介绍了用于表示介电常数分布的所选参数化模型。“MMWI场景”部分详述了所考虑的MMWI场景。在“失真玻恩迭代法”部分，概述了所应用的DBIM公式。“数值算例”部分展示了数值示例的结果。最后，“结论”部分提供了结论性评述。
介电常数分布的参数化建模
为了缓解逆问题的病态性，使用合适的基函数来表示未知的介电常数分布是有利的，从而降低解空间的维度。在头部上部区域，头部组织近似贴合头部的整体形状。因此，采用遵循这种外部几何形状但同时允许变化的基函数是合理的，这些变化可以有效地用低阶球谐函数建模。头部表面上的每个点可以表示为(1)，其中是l阶m次的球谐函数，L是球谐函数的最大阶数，alm是展开系数。
类似地，介电常数分布，在球坐标系中表示为三维（3D）函数ε(r, θ, φ)，可以用由径向函数（建模径向依赖性）和球谐函数（捕捉角度依赖性）乘积构成的正交基函数级数来近似。该展开形式如式(2)所示，其中cnlm是对应的展开系数，是归一化径向坐标，定义为(3)（即距离头部表面的归一化径向距离），pn(?r)是径向脉冲或指示基函数，n = 1, ..., N，定义如式(4)。
在数值实现中，空间被均匀划分为K个voxel，每个voxel具有恒定的介电常数εk，k = 1, ..., K。以矩阵形式，这可以写为(5)，其中ε是未知介电常数向量(6)，c是未知系数向量(7)，Q是系数的总数Q = N(L+1)2。在(7)中，使用了简化符号cj = cnlm，其中j是对应于三元组(n, l, m)的单一索引，j = 1, ..., Q。此外，H是变换矩阵，其元素计算如(8)，其中(rk, θk, φk)是第k个voxel的球坐标，?rk = [Rs(θk, φk) - rk] / Rs(θk, φk)。
MMWI场景
图1(a)展示了一个MMWI场景，由与头部紧密接触的天线阵列组成。在所考虑的示例中，头部模型由六种均匀组织构成。假设的初始背景模型如图1(b)所示，代表一个具有等于脑组织介电常数、并包含一个对应于皮肤的薄层的均匀头部。其底层物理模型是散射方程，在玻恩近似下，由(9a)给出，其中sij和s?ij分别是在实际系统和假设背景下，第i个和第j个天线之间的传输系数，ε是未知的头部介电常数，εb是已知背景模型的介电常数，Δε是介电常数差，r'是头部内部一点的位置向量，E_i^{i(j)}是在i(j)个天线发射时在背景模型中计算的电场向量，a_{i(j)}是i(j)端口入射功率波的复振幅，M是阵列大小，v′是头部已知的体积，ω是角频率。
在离散情况下，(9)变为(10)，其中t_k是第k个voxel中心的位置向量，εk和ε_{k,b}分别是实际头部和背景模型中第k个voxel的恒定介电常数，Δεk是它们的差值，Δv是每个voxel的体积，K是占用体积被划分成的voxel总数。以矩阵形式，(10)变为(11)，其中Δs如(12a)，L如(12b)，Δε如(12c)。
研究人员使用所采用的基函数(5)表示介电常数差向量(12c)为(13)。将(13)代入(11)，得到联系测量值和多项式系数的最终方程组(14)，其中A = LH。为求解(14)，采用了Tikhonov正则化。相应的代价函数为(15)，其中λ是正则化系数。(15)的解为(16)，估计的介电常数向量为(17)，其中(18)是背景介电常数向量。基于近似解(17)，通过假设所有强度低于所选阈值的voxel对应于低介电常数的颅骨，提取了关于颅骨的信息。一旦识别出颅骨，其内部被假定为大脑。然后通过添加一个贴合外部头部形状的薄表面层来完成模型。最终，获得了一个由四种均匀组织构成的参数化头部模型。
失真玻恩迭代法
一旦获得参数化头部模型，研究人员使用DBIM计算其各区域的复介电常数。
DBIM在MMWI中已被广泛应用，有显著应用如下：一个经过校准的、基于频率分割的DBIM二维（2D）微波层析成像算法用于头部成像，并在头部体模和健康志愿者身上得到验证。一种3D自适应聚类DBIM方法用于中风检测，并在假设成像域中组织数量及其近似电磁特性已知的情况下，使用头部扫描仪系统和头部体模进行了评估。一种使用DBIM和两步迭代收缩阈值（DBIM-TwIST）的中风检测算法，在2D头部体模上进行了实验验证，该体模被天线环包围。一种具有感兴趣区域限制方案的、分辨率增强的DBIM用于微波乳腺成像，在2D数值模型上产生了有希望的结果。扩展玻恩迭代法结合局部非线性近似，在中风和乳腺肿瘤的真实2D数值模型中实现了快速收敛和稳健的异常检测。
与之前的方法不同，研究人员在应用DBIM时利用了先前估计的组织边界，将未知数减少到组织的数量，无需特定的患者特异性先验信息。
DBIM公式
DBIM通过迭代更新背景模型中的介电常数估计值ε_{k,b}，以及与ε_{k,b}估计相对应的Δsij和E_i^{i(j)}来求解(10)。为了启动DBIM算法，必须初始化介电常数ε_{1,b}^{(1)}, ε_{2,b}^{(1)}, ..., ε_{K,b}^{(1)}。
在第l次迭代中，线性方程组(11)为(19)，其中Δs^{(l)}、L^{(l)}和Δε^{(l)}是第l次迭代中的矩阵值。系统(19)使用截断奇异值分解求解。正则化解为(20)，其中u_n和v_n是矩阵L^{(l)}的奇异向量，σ_n是对应的奇异值，n_{\max}是截断索引，根据条件10 log_{10}(σ_{n_{\max}} / σ_1) < -20 dB计算。下一次迭代中的介电常数估计计算如(21)。
当向量Δε^{(l)}的每个成员的幅值小于某个预定义阈值（在本工作中所有实验均设置为0.01）或达到最大迭代次数时，DBIM算法停止。
具有已知域边界的DBIM公式
如果感兴趣的区域可以划分为具有已知边界和均匀介电常数的子区域，则未知参数的数量减少到这些子区域的数量。相应的测量模型可以表示为(22)，其中W是大小为K × S的变换矩阵，S是子区域的数量。W的第(k, j)个元素在第k个voxel属于第j个区域时为1，否则为0。因此，未知向量为(23)，其中Δε_j^{(l)}指第j个区域在第l次迭代中的对比度。
数值算例
作为一个数值示例，考虑了图1(a)所示的真实头部体模。该体模基于NEVA女性解剖模型，并在WIPL-D Pro仿真环境中为MMWI应用进行了调整。体模由六个均匀组织层组成，如图1(a)所示，其各自的介电常数值列于表1中。研究人员假设体模的外表面是已知的，并由图2中的三角网格给出。实际上，仅在头盔下方的头部上部需要几何信息。测量设置包括一个配备21个紧密排列的梯形贴片天线的头盔，每个天线工作在1 GHz，并针对存在头部模型时的性能进行了优化。关于头部体模和测量系统配置的更多细节见参考文献。
论文中展示的所有结果均使用台式计算机生成：Intel Core i7-9700 CPU @3 GHz，NVIDIA GeForce GTX 750 Ti GPU，32 GB RAM，Windows 10操作系统。
参数化模型创建
第一步涉及在定义的感兴趣区域内提取头部表面上的点，其中颅骨和上覆组织近似贴合。然后使用这些点构建一个封闭表面。由于图2所示的三角网格相对粗糙，研究人员首先插值了额外的表面点，随后提取了靠近头盔的点，如图3(a)所示。为确保几何均匀性，然后将每个切片中的边界点映射到极坐标中的均匀角度网格上，为每一层提供一致数量的点。最后，为了完成封闭表面，通过缩放最底层并向下平移来生成额外的点，如图3(b)所示。由此产生的外推表面点作为球面近似的匹配点，结果如图3(c)所示。在计算球谐函数近似时，假设xy平面与头部周长最大的层重合。展开系数用于计算式(2)中的归一化半径。
对于微波成像，使用在封闭表面内且在z_min = -0.02 m上方的均匀3D网格。这导致总共有34,947个边长为4 mm的voxel，如图2(d)所示。较低的区域由于场穿透不足而被省略，这限制了它们用于可靠成像的效用。
关于介电常数展开(2)，沿径向方向，设置N = 25个脉冲用于径向变化，L = 3用于球谐函数展开的最大阶数。脉冲数量的选择是空间分辨率和数值稳定性之间的折衷。球谐函数的阶数保持较低以保持组织的贴合性，同时仍允许一定程度的灵活性。
作为假设的背景模型，使用了图1(b)中的两组织头部体模。外部组织模拟皮肤和脂肪，其介电常数设置为28.50–j8.96，作为皮肤和脂肪介电常数值的平均值。采用的厚度为3.5 mm。内部组织的介电常数设置为大脑的介电常数值。
选定横截面（从最低层到最高层共30个中的部分）的结果如图4所示。在每一对图像中，右侧图像表示重建的对比度（介电常数差的实部），左侧图像显示真实对比度，两者都相对于所采用的背景模型。
通过敏感性分析研究了脉冲数量和球谐函数最大阶数的影响。不同参数选择的重建结果如图5所示。具体来说，将脉冲数量增加到N = 30并将球谐函数最大阶数增加到L = 4仅导致重建的介电常数分布发生轻微变化。这些结果表明，所提出的方法对参数的适度变化具有鲁棒性，而所选值N = 25和L = 3在重建精度和计算效率之间取得了良好的平衡。
值得注意的是，所提出公式化中的未知系数总数由N(L+1)2给出，对于N = 25和L = 3，产生400个未知数。作为比较，在研究人员之前基于多项式基函数的工作中，未知数为(Q+1)3 = 133 = 2197，其中Q表示最大多项式阶数。未知系数数量的大幅减少改善了逆问题的条件数，并有助于重建的稳定性，而不会损害精度。
此外，评估了所提出方法对测量噪声的鲁棒性。为此，在原始六区域模型的散射参数中引入了加性白高斯噪声（AWGN）。
在评估可接受的噪声水平时，研究人员考虑初始两层背景模型的信号与建模误差功率比约为11 dB，表明起始点存在显著的建模失配。此外，在中风诊断等应用中，与病理区域相关的信号通常比测量信号低50-70 dB。因此，噪声水平超过此范围的MWI系统将不适用于可靠的中风检测。
尽管如此，在分析中，研究人员考虑了低至10 dB的SNR值以评估所提出方法的鲁棒性。图6展示了一个远离表面的选定横截面切片在SNR为20 dB和10 dB时的重建头部图像，以及使用无噪声数据获得的相应重建。
获得的结果表明，重建的MWI图像中明显的质量下降开始出现在约10 dB SNR时（此时噪声水平变得与初始粗糙两层模型相关的建模误差相当）。
在下一步中，将连续图像离散化以提取颅骨边界。为了确定最佳阈值，使用了直方图。图7显示了一个顶部横截面（z = z 27）的直方图。为了更清晰，研究人员展示了估计的介电常数实部而非对比度。在最靠近天线的层（即顶部层），直方图始终显示三个明显的区域：一个中心在46左右，对应于软组织；另一个在50以上，对应于脑组织；第三个在43以下，与颅骨相关。基于这些观察，选择了阈值TH = 44。
在使用从直方图导出的阈值对图像进行二值化后，识别出颅骨的外部和内部边界，然后用于生成颅骨模型。一旦生成颅骨模型，通过将颅骨插入图1(b)的背景模型中，获得近似的头部模型。该模型包含四个不同的区域：（1）颅骨，（2）大脑，（3）皮肤与脂肪，以及（4）占据剩余体积的软组织。在对应于有效场穿透深度的深度处，大脑和颅骨区域被压平，因为发现来自该深度以下层的散射场可忽略不计。获得的体模如图8所示。
DBIM结果
为了完成患者特异性模型，需要计算四个区域的介电常数。为此，采用了“具有已知域边界的DBIM公式”小节中描述的DBIM流程。使用与上一小节相同的搜索空间，包含头部内34,947个体素。一旦执行DBIM并获得介电常数值，就使用WIPL-D Pro计算患者特异性模型的散射矩阵。为了评估所获得的模型对原始六区域模型的近似准确性，计算了散射参数的相对均方根误差（RRMSE），记为δ：(24)，其中skj和s_{kj}^{ref}分别是使用从DBIM获得的模型和具有六种组织的原始模型（图1(a)）计算的散射参数。
作为DBIM的第一次测试，研究人员假设来自图1(a)的六区域体模的边界是精确已知的。所有组织的介电常数初始值设置为参考文献中的平均值（34.64–j12.68）。在这种情况下，DBIM在五次迭代后收敛到表1中的实际值。
此外，为了评估使用四种组织模型所能获得的最大精度，从具有六种组织的原始模型创建了两个参数化模型，假设组织间边界已知。具体来说，从颅骨中移除了面骨，并使用球谐函数对真实颅骨边界进行了参数化，直到两个不同的切割平面。因此，获得了参数化模型1（图9(a)）和参数化模型2（图9(b)）。
对于这两个参数化模型和图7中的患者特异性模型，DBIM在10次迭代内收敛（患者特异性模型7次，每个参数化模型6次）。每个WIPL-D模型大约有33,000个未知数，因此每次DBIM迭代耗时约50分钟。获得的介电常数值和RRMSE列于表2。
当使用精确的颅骨边界进行参数化时，RRMSE值分别为3.54%和3.55%。另一方面，对于从重建的介电常数创建的患者特异性模型，获得了4.27%的RRMSE。作为比较，通过均匀头部体模近似原始模型所能获得的最小RRMSE为12.27%。因此，使用相同先验信息的患者特异性四组织模型实现了近三倍更低的误差值。
对于用于生成患者特异性模型的初始两层模型，RRMSE为26.3%。等效地，信号与建模误差功率比从初始两层模型的大约11 dB改善到患者特异性模型的约27.4 dB。
此外，通过将AWGN加入到原始六区域模型的模拟散射参数中，评估了DBIM对测量噪声的鲁棒性。考虑了对应于SNR值为20、15和10 dB的噪声水平。
为了将噪声对DBIM算法的影响与其对构建四区域体模的影响隔离开来，采用了在无噪声条件下获得的组织几何结构。相应的介电常数估计值和RRMSE值报告于表3中。在所有情况下，DBIM算法在10次迭代内收敛。
如表3所示，重建精度随着SNR的降低而逐渐下降，反映了噪声影响的增长。对于SNR为20和15 dB的水平，估计的介电常数值仍然接近无噪声情况下获得的值，RRMSE仅适度增加。即使在SNR为10 dB时，DBIM算法也能可靠地收敛，尽管内部组织（如颅骨和大脑）的介电常数偏差变得更加明显。这些发现与图6中呈现的MWI结果一致，其中相对于总信号的10 dB SNR似乎代表了一个阈值，超过该阈值会发生更显著的退化。考虑到中风相关信号比总信号弱50 dB以上，可以得出结论，所研究的噪声水平不太可能显著影响患者特异性背景模型的主要应用。
结论
本工作提出了一种构建患者特异性头部模型以增强定量MMWI的先进方法。该框架整合了通过插值表面点的球谐函数近似获得的头部几何形状的解剖贴合参数化。这种方法确保了感兴趣区域几何上一致的表示，特别是在与微波测量相关的上颅部分，并为直接从测量数据构建个体化模型提供了实用途径。
所提出的方法在基于NEVA女性解剖模型的现实六组织头部体模基础上进行了数值评估。仿真在WIPL-D Pro环境中进行，使用工作于1 GHz的21天线头盔阵列。首先应用线性化重建来估计介电常数对比度的实部，然后进行基于直方图的分割以描绘颅骨边界并定义四个不同的组织区域。随后，使用DBIM定量估计这些组织的介电常数值。由此产生的患者特异性四组织模型获得了4.27%的RRMSE，与使用精确颅骨边界获得的3.54%接近，并显著优于均匀头部近似（12.27% RRMSE）。
这些结果证明，所提出的建模策略显著提高了DBIM的准确性，同时保持了计算效率。能够直接从微波数据生成可靠的个体化背景模型，代表了增强定量逆散射算法初始化和收敛的重要一步，从而推动MMWI走向实际临床应用。
