基于交替方向乘子法(ADMM)直接重建心脏11C-acetate PET动力学参数图

何升级，刘 帅，张 辉，霍 力，尚 斐*

(1.北京理工大学生命学院生物医学工程系，北京 100081；2.清华大学医学院生物医学工程系，北京 100084； 3.中国医学科学院 北京协和医学院 北京协和医院核医学科，北京 100730)

11C-acetate是临床非侵入式定量评估心肌代谢和血流灌注的重要示踪剂之一[1-2]。目前多采用单组织房室模型(即单室模型)计算心脏11C-acetate PET动态成像动力学参数[3-4]，一般分为两步：①基于PET投影数据重建PET动态图像；②采用动力学模型逐像素拟合时间-活度曲线，得到动力学参数图。为提高动力学参数图像质量，有学者[5]将图像重建与动力学模型相结合而提出直接重建算法，即根据投影数据直接重建动力学参数图，以降低噪声对参数拟合的影响。现有相关研究[6-8]多采用图模型(如Patlak和Logan模型)重建动力学参数图。本研究观察基于交替方向乘子法(alternating direction method of multipliers, ADMM)直接重建心脏11C-acetate PET动力学参数图的应用价值。

1 资料与方法

1.1 数据来源 PET/CT数据来自1例有30年饮酒史、无其他心脏疾病，于2017年7月在北京协和医院接受PET/CT检查的51岁男性早期酒精性心肌病患者。本研究获医院伦理委员会批准(伦理审批号：20160012)，检查前患者签署知情同意书。

1.2 仪器与方法 采用PoleStar m660 PET/CT扫描仪。嘱患者仰卧，先行低剂量肺部CT扫描，管电压120 kV，管电流140 mA；之后经静脉注射11C-acetate 740 MBq，行40 min动态PET心脏扫描；基于CT图像进行衰减校正后，采用有序子集最大期望值法(ordered subsets expectation maximization, OSEM)及飞行时间技术重建PET图像，图像大小192×192×117，分辨率3.147 mm×3.147 mm×1.866 mm，共得到53帧重建图像，包括15帧10 s、15帧30 s、16帧60 s及7帧120 s图像。

1.3 构建模型

1.3.1 PET动态图像重建模型 将第n帧PET图像记为xn∈N×1(N表示图像的像素数目，表示实数域)，将测得的第n帧投影数据记为yn∈M×1，其期望为M×1(M表示响应线数目)；则与xn的关系为：

(1)

式中，E为期望符号，P∈M×N为系统矩阵，rn∈M×1表示第n帧随机和散射事件的期望；则动态图像重建似然方程L为：

(2)

式中，y为动态投影数据，y=[y1,y2,…,yT]；x为动态图像，x=[x1,x2,…,xT]；T为图像重建帧数；i表示投影数据射线的索引。

1.3.2 心脏11C-acetate单室模型 见图1。

图1 心脏11C-acetate单室模型

其中，CP(t)和CT(t)分别表示血液和心肌中的11C-acetate浓度随时间t的变化值；K1和k2表示11C-acetate在血液和心肌之间交换的速率常数，是临床量化评估心肌的重要指标，K1与心肌血流量相关，k2与心肌耗氧量相关[1,4]。该模型用微分方程表示为：

(3)

对公式(3)求解，可得：

CT(t)=K1CP(t)⊗e-k2t

(4)

式中，⊗代表卷积运算；由于PET成像存在容积效应，需要对心肌区域进行溢出及恢复校正[9]，校正公式为：

CS(t)=(1-vl)×CT(t)+vl×CP(t)

(5)

式中，CS(t)为PET图像中实际测得的心肌区域示踪剂浓度，vl∈[0,1]为血液容积的比例因子。

图2 直接重建算法流程图

对公式(4)和(5)进行拟合，可得到参数K1、k2和vl，参数集合为θ=[K1,k2,vl]。

1.3.3 直接重建流程 将重建PET动态图像和单室模型参数求解两个独立步骤合二为一，以实现直接重建参数图像。动力学模型中的参数θ与动态图像x关系为x=K(θ)，其中K(θ)表示上述心脏11C-acetate单室模型；结合公式(1)可知，参数θ与投影数据的关系为：

(6)

(7)

对θ则可通过如下公式进行优化并求解：

(8)

(9)

再通过增广拉格朗日法将公式(9)转化为无约束优化问题，对应优化目标函数为：

(10)

最后通过ADMM框架进行迭代求解：

(11)

(12)

μn+1=μn+xn+1-K(θn+1)

(13)

式中，ρ为惩罚因子，μ为对偶变量，n表示迭代次数。求解流程见图2。

1.4 动力学参数重建

1.4.1 参考标准 对本例PET图像所示心脏区域进行裁剪，裁剪后图像大小为64×64×48；再由1名具有10年工作经验的放射科主治医师勾画左心室心外膜和心内膜，获得左心室血池和心肌区域，并基于左心室血池获取CP(t)，根据公式(2)和(3)计算得到K1、k2和vl参数图；最后根据美国心脏协会提出的17节段模型[10]将左心室心肌区域参数图映射生成参数图，并以之作为标准，间接重建采用OSEM方法重建动态图像。

1.4.2 直接重建 采用MATLAB R2017a软件，操作系统为64位Windows 10，CPU为英特尔i7-9700F，16 GB内存。基于单室模型由参数图获得53帧理想动态图像，并对其进行120个角度(0～180°平均划分)投影后，加入泊松噪声；采用开源NiftyRec工具箱[11]实现基于ADMM的直接重建算法，总迭代次数15，惩罚项ρ为2.0×10-7，子迭代次数2，初始μ0为0；并通过滤波反投影法重建得到初始图像x0，通过11C-acetate单室模型对x0进行逐像素曲线拟合，得到初始参数图θ0。见图3。

1.5 图像评估 通过偏差(bias)评估图像重建结果，其定义式为：

(14)

式中，Iest为通过直接或间接重建方法计算得到的值，Igs为参考标准。

1.6 统计学分析 采用线性回归分析参考标准与直接、间接重建图像的相关性，并计算Pearson相关系数；以0<│r│<0.3为弱相关，0.3≤│r│<0.5为低度相关，0.5≤│r│<0.8为中度相关，│r│≥0.8为高度相关。P<0.05为差异有统计学意义。

2 结果

间接重建耗法耗时23 min，直接重建法耗时140 min。

图3 心脏11C-acetate PET动力学参数图直接重建流程图

图4 参考标准、间接重建及直接重建参数图及其偏差

参考标准与间接重建得到的K1、k2和vl值的平均偏差分别为2.97%、1.61%和10.19%，与通过直接重建得到的K1、k2和vl值的平均偏差分别为1.16%、1.06%和4.30%。见图4。

参考标准与通过间接重建得到的K1(r=0.94,P<0.001)、k2(r=0.94,P<0.001)和vl值(r=0.91,P<0.001)及与通过直接重建得到的K1(r=0.99,P<0.001)、k2(r=0.97,P<0.001)和vl值(r=0.99,P<0.001)均高度相关，且参考标准与直接重建法的相关系数略大。见图5。

3 讨论

心脏11C-acetate PET动态成像有助于评估心肌代谢和血流灌注，目前多采用图模型(如Patlak和Logan模型)重建其动力学参数图[6-8]。GONG等[7]结合18F-FDG PET数据及MRI所示结构信息，以基于ADMM的非监督深度学习法直接重建，得到脑部Patlak动力学参数图。GALLEZOT等[8]对传统Logan模型进行改良，采用最大期望值法(expectation maximization, EM)直接重建，得到脑部11C-PBR28动力学参数图。上述研究结果均表明，通过图模型直接重建法得到的参数图质量较高。

图5 参考标准与间接重建K1(A)、k2(B)和vl(C)参数图及与直接重建K1(D)、k2(E)和vl(F)参数图的相关性分析图

SHI等[12]发现，计算11C-acetate动力学参数时，基于OSEM重建图像的稳定性和鲁棒性更好。本研究采用OSEM间接重建心脏11C-acetate PET动力学参数图，并基于ADMM实现直接重建，结果表明直接及间接重建得到的K1、k2和vl均与参考标准高度相关；直接重建法与参考标准的相关系数略大于间接重建，且直接重建得到的K1、k2和vl与参考标准间的偏差更小。分析原因，主要在于间接重建过程中动态图像易受噪声干扰而在拟合时出现奇异值，降低计算动力学参数的准确性，同时由于参数拟合时忽略了动态图像像素间的相关性，导致所获参数图常具有较大噪声；而在直接重建迭代过程中，每一次迭代均包含图像重建和动力学参数拟合两个步骤，且优化过程中的θ与x相互约束，加之重建图像时考虑了空间相关性，使得参数图的拟合误差较小，抗噪性能则有所提高。

本研究间接重建法耗时23 min，直接重建法耗时140 min，主要原因在于采用单室模型进行直接重建，其动力学参数与PET数据的关系为非线性，使算法较为复杂，增加了时间成本；但非线性单室模型是以示踪剂在人体内的代谢过程为基础而建立的动力学数学模型，与人体内生化过程紧密相关，有助于得到更具生理意义指标，以定量评估靶器官功能。

目前已有学者尝试将正则化方法，如全变分正则化[13]、核方法[14]等，加入直接重建框架中，以提高重建图像质量，且已知采用深度神经网络重建PET图像可提高成像质量[15]。后续研究可将深度神经网络嵌入直接重建PET参数图中，以进一步提高其动力学参数图的图像质量。

综上所述，相比间接重建方法，基于ADMM直接重建心脏11C-acetate PET动力学参数图的准确性更高。