2020, Vol. 41
2. 沈阳东软智能医疗科技研究院有限公司, 辽宁 沈阳 110167
2. Neusoft Research of Intelligent Healthcare Technology, Co., Ltd., Shenyang 110167, China
中心主动脉压力波形(central aortic pressure wave, CAPW)是评估人体心血管功能、高血压治疗用药效果、动脉硬化的发展，以及冠心病的病变程度的重要指标[1], 在临床诊断中有重要的应用价值.目前主动脉压的测量主要是有创测量.有创测量因其对病人造成创伤、费用昂贵等原因难以推广, 故人们提出了许多无创测量方法, 其中主要的是传递函数法[2].
传递函数法一般测量上肢动脉压力波形, 即利用上肢动脉压力波形通过函数转换或模型推导得到中心主动脉压力波形.在众多传递函数方法中, 由澳大利亚的AtCor Medical公司研发的SphygmoCor系统能够比较准确地利用桡动脉压力波形重建CAPW, 是一种性能非常好的传递函数法[3].传递函数法需要通过先验知识建立传递函数, 但不同个体的心血管系统并不同, 所以该方法在实际应用中有一定局限, 精度也需要进一步提高.
近年来, 盲系统辨识技术由于其能够动态、实时、无创地重建中心动脉压力波形, 被广泛应用于生物医学领域[4].本文选择了基于子空间方法的主动脉压力波形重建[5-11], 同步测量了多种外周动脉压力波形(peripheral arterial pressure waves, PAPW), 并按照不同组合进行主动脉压力波形重建; 以AtCor设备获得的主动脉压力波形为标准[3], 使用动态时间规整算法进行比较, 探究不同的外周动脉压力波形组合对重建中心动脉压力波形的影响, 以期在临床应用中, 能够找出一种合理的方法准确地进行主动脉压力波形重建.
1 盲系统辨识方法 1.1 盲系统辨识基本理论盲系统辨识(blind system identification, BSI)是一种诞生于30多年前的信号处理技术.盲系统辨识技术就是仅通过输出信号来获取信道的脉冲响应的过程.近年来, 盲系统辨识技术已经广泛应用于许多不同的场景, 如移动通信、地质勘探、图像重建及医学应用.至今已经有了众多方法能够实现盲系统辨识, 包括子空间法[5](subspace, SS)、交叉相关法(cross-relation, CR)、两步极大似然法(two-step maximum likelihood, TSML)等.在生物医学领域, 主动脉压力波形的重建是一个重要的研究方向.主动脉有着丰富的心血管系统的病理生理学信息[2], 因此无创测量主动脉压力波形尤为重要.
可辨识系统需要满足以下三个必要条件[10]:
① 所有通道间没有公共零极点, 即各通道特性要足够地不同;
② 输入信号要足够复杂以激励各个通道, 其模式个数不小于L+2(L为系统传递函数阶数), 例如:输入序列不能为零、常数或单个正弦信号;
③ 输出信号足够长, 每一通道输出的信号长度N>2L+1.
本文根据盲系统辨识的必要条件, 探究不同外周动脉组合对主动脉压力波形重建的影响.
1.2 基于子空间的盲系统辨识算法对于一个单输入多输出(single input multiple output, SIMO)系统, 将每一个输出定义为一个长度为N的向量, 系统中M个信道的输出就是一个M×N的矩阵.设y(n)为M个信道上的观测样本, 输入序列定义为长度是N的向量s(n), v(n)是M个子信道上的噪声.输入信号s(n)中各元素独立同分布且均值为零, 信道噪声v(n)为高斯白噪声.同样地, M个信道的系统传输矩阵就可以表示为一个M×N的矩阵HM.那么一个SIMO系统就可以表示为
|
(1) |
|
图 1 SIMO系统模型 Fig.1 SIMO system model |
子空间法是针对输出矩阵y(n)的自相关协方差矩阵Ry和输入信号s(n)的自相关协方差矩阵Rs，利用信号子空间和噪声子空间的正交关系求解.对式(1)两边分别求自相关协方差矩阵, 得到式(2):
|
(2) |
式中Rv为v(n)的自相关协方差矩阵.
输入序列s(n)是独立同分布的随机变量序列, Rs满秩.对Ry进行奇异值分解, 可以得到信号子空间和噪声子空间的矩阵.系统传输矩阵和噪声子空间UN中的任意向量正交, 因此能得到式(3):
|
(3) |
实际中只能得到相关矩阵Ry特征向量的估计, 因此式(3)只能近似求解.根据最小二乘法, 定义目标函数:
|
(4) |
式中h是信道系数向量.进一步对目标函数求解, 可得
|
(5) |
对Q进行特征值分解, 最小特征值对应的特异化特征向量就是信道系数向量h的估计.根据信道系数向量, 就能够计算得到输入信号.
2 研究方案 2.1 数据采集本文采用的主动脉压力波形是由SphygmoCor系统(AtCor Medical公司, 澳大利亚)采集了测试对象桡动脉压力波形后重建的, 该系统重建主动脉压力波形使用了传递函数法, 此方法被视为传递函数法重建主动脉波形的“金标准”.AtCor设备采集一路桡动脉压力波形进行中心动脉压力波形的重建.经过验证, AtCor使用左右两侧桡动脉进行重建的结果无差异, 故本文统一采用右侧桡动脉脉搏波的重建结果作为实验对照组.
盲系统辨识所用的桡动脉脉搏波和足背动脉脉搏波是实时同步数据[12], 测试对象采取平卧位, 并使用了相同的压力传感器.该传感器的返回值是压力数据, 因此首先要对数据重新标定.本文的压力标定由水银血压计测得肱动脉压力, 根据临床测量的平卧位外周动脉压力关系, 对桡动脉压力波形及足背动脉压力波形进行校准[13-14].
在采集过程中, 首先使用AtCor设备获取测试对象的主动脉压力波形的估计, 作为对照.然后同步测量不同组合的两路外周动脉压力波形:左侧桡动脉和左侧足背动脉、左侧桡动脉和右侧足背动脉、右侧桡动脉和左侧足背动脉、右侧桡动脉和右侧足背动脉, 以及左右两侧的桡动脉、足背动脉.
2.2 重建主动脉压力波形本文采用基于子空间的盲系统辨识算法.该方法是基于Hankel矩阵的秩来定阶, 并结合子空间法用两路外周脉搏波重建主动脉波形[7].根据盲系统辨识的必要条件, 选择桡动脉和足背动脉作为SIMO系统的两个输出信道, 中心动脉脉搏波作为输入信号.脉搏波可以近似地视为独立同分布的随机信号.采集的桡动脉和足背动脉的脉搏波都足够长, 满足可辨识条件③.盲系统辨识的信道参数向量选择2阶, 满足可辨识条件②.如图 2所示, 根据两路外周压力脉搏波波形可以估计出主动脉压力波形.
|
图 2 实测的两路PAPW与估计的CAPW Fig.2 Measured PAPW and estimated CAPW (a)—实测左右桡动脉脉搏波；(b)—实测左右足背动脉脉搏波；(c)—实测桡动脉和足背动脉脉搏波. |
根据估计结果可以看出, 由桡动脉和足背动脉的压力波形进行盲系统辨识能够很好地重建出主动脉压力波形；而使用左右两侧桡动脉和足背动脉的压力波形的重建结果是两路信号的平均, 不能实现对中心动脉压力波形的重建.两侧的桡动脉或足背动脉虽然属于不同的信道, 但生理结构上非常相似, 即信道的特性非常相似, 不满足可辨识条件①.
2.3 用动态时间规整进行重建结果对比AtCor设备给出的主动脉波形是几个周期的平均, 因此首先将盲系统辨识得到的主动脉压力波形进行平均.图 3所示为提取的连续几个周期的脉搏波, 以各周期舒张期压力为对照, 提取各周期的脉搏波并求平均, 再进行比较.
|
图 3 连续提取的若干周期的脉搏波 Fig.3 Some cycles of pulse waves extracted continuously |
由于无法得知主动脉的实际压力, 因此将所有波形均进行归一化处理, 然后用动态时间规整(dynamic time warping, DTW)算法进行对比, 结果如图 4所示.
|
图 4 重建结果与实际测量结果的对比 Fig.4 Reconstructed results compared with measured results |
动态时间规整算法可以找到两个时间序列的最佳配准路径[15], 经常被用来计算两个时间序列的相似程度.该算法通过对一个时间序列进行拉伸或压缩的非线性调整来逼近另一个时间序列, 而这个调整的程度被用来衡量两个时间序列的相似程度.由于脉搏波不是随机信号, 各周期的点数并不相同, 导致动态时间规整结果产生误差, 所以本文对比的误差是动态时间规整数值根据点数所求的平均, 即
|
(4) |
式中:e为平均误差；X和Y是两个任意长度的时间序列; wp为规整路径中的第p个元素; l为序列Y的长度.规整路径中的第一个元素是两序列起始点的路径, 最后一个元素是两序列末尾点的路径.规整路径按照时间顺序选择两序列上点的最短欧式距离.在实验中, X序列是“金标准”, Y序列是盲系统辨识的结果.为了对比不同组合的盲系统辨识结果, 以Y的长度对规整路径的和取平均值.
如图 5所示, 采集30个不同年龄个体的桡动脉与足背动脉脉搏波, 图中lwlf, lwrf, rwlf, rwrf分别表示左侧桡动脉和左侧足背动脉、左侧桡动脉和右侧足背动脉、右侧桡动脉和左侧足背动脉、右侧桡动脉和右侧足背动脉.根据重建结果, 按照相同方法求其动态时间规整的数值.可以看出4种组合下盲系统辨识的结果相似, 动态时间规整的数值十分接近, 因此可以认为4种组合均能很好地重建出主动脉压力波形.对4组结果进行单因素方差分析, 其F值约为0.356, 显著性检验P值达到0.785, 可以认为不同的外周动脉组合对重建中心动脉压力波形无显著性影响.
|
图 5 30个测试者CAPW重建结果对比 Fig.5 Comparison of reconstruction of CAPW in different individuals |
从图 5看出, lwlf的第14组结果、rwlf的第27组结果、lwrf的第29和30组结果因为噪声而与对照组有较大差异, 其他数据都表明使用左侧桡动脉和右侧足背动脉、右侧桡动脉和右侧足背动脉重建得到的中心动脉压力波形十分稳定.如表 1所示, 根据实验结果的均值与方差得出, 右侧桡动脉和右侧足背动脉的组合能够相对准确地得到重建结果.
|
表 1 四个组合下重建结果的方差和均值 Table 1 Variance and mean of different combinations |
为了验证这一结果的稳定性, 本文在同一个体上采集了6组外周动脉波形, 每组数据均间隔足够的时间才进行采集, 其中rwrf因采集噪声的影响产生异常点.从图 6可以看出, 4种组合方式均可重建出中心动脉压力波形, 左右两侧的桡动脉和右侧足背动脉组合使用能够获得相对较好的结果.在同一测试者的6组实验数据中, 右侧桡动脉和右侧足背动脉的重建效果要明显优于其他组合.30组实验中右侧桡动脉和右侧足背动脉的组合整体优于其他组合, 但并不是总呈现这一趋势.这一现象仅从信道参数和特性上难以验证.
|
图 6 同一测试者的CAPW重建结果对比 Fig.6 Comparison of reconstruction of CAPW in the same individual |
1) 左右桡动脉组合或者左右足背动脉组合无法重建出主动脉压力波形; 不同侧、相同位置的脉搏波组合也不能重建出主动脉压力波形.
2) 使用桡动脉和足背动脉的组合能够很好地重建主动脉压力波形，左右两侧不同的组合对辨识结果无显著性影响.
3) 由表 1、图 5和图 6中能够看出, 右侧桡动脉与右侧足背动脉的重建结果更加稳定.
4) 根据盲系统辨识的必要条件, 输入信号, 即采集的桡动脉脉搏波和足背动脉脉搏波, 已经足够复杂, 并且输入长度也足够长, 因此满足辨识条件②和③; 然而, 若使用两路桡动脉脉搏波或两路足背动脉脉搏波进行辨识, 其生理结构上十分相似, 作为盲系统辨识的两路通道会出现相同的零极点, 不满足盲系统辨识条件①.若使用桡动脉脉搏波和足背动脉脉搏波进行盲系统辨识, 满足了必要条件①, ②和③, 能够准确地辨识出中心动脉脉搏波.
| [1] |
Hashimoto J. Central hemodynamics and target organ damage in hypertension[J]. Tohoku Journal of Experimental Medicine, 2014, 233(1): 1-8. DOI:10.1620/tjem.233.1 |
| [2] |
O'Rourke M F. Influence of ventricular ejection on the relationship between central aortic and brachial pressure pulse in man[J]. Cardiovascular Research, 1970, 4(3): 291-300. DOI:10.1093/cvr/4.3.291 |
| [3] |
Yao Y, Xu L, Sun Y, et al. Validation of an adaptive transfer function method to estimate the aortic pressure waveform[J]. IEEE Journal of Biomedical and Health Informatics, 2016, 21(6): 1599-1606. |
| [4] |
Hahn J O, Reisner A T, Asada H H. Blind identification of two-channel IIR systems with application to central cardiovascular monitoring[J]. Journal of Dynamic Systems, Measurement, and Control, 2009, 131(5): 051009-051023. DOI:10.1115/1.3155011 |
| [5] |
Mayyala Q, Abed-Meraim K, Zerguine A. Structure-based subspace method for multichannel blind system identification[J]. IEEE Signal Processing Letters, 2017, 24(8): 1183-1187. DOI:10.1109/LSP.2017.2715418 |
| [6] |
Qiu W, Hua Y.Performance comparison of three methods for blind channel identification [C]//IEEE International Conference on Acoustics, Speech, and Signal Processing.Atlanta, 1996: 2423-2426.
|
| [7] |
刘文彦, 徐礼胜, 李宗鹏, 等. 基于子空间法的主动脉压力脉搏波重建[J]. 东北大学学报(自然科学版), 2019, 40(9): 1240-1245. (Liu Wen-yan, Xu Li-sheng, Li Zong-peng, et al. A subspace algorithm for reconstructing aortic pressure wave[J]. Journal of Northeastern University(Natural Science), 2019, 40(9): 1240-1245.) |
| [8] |
Guţǎ M, Yamamoto N. System identification for passive linear quantum systems[J]. IEEE Transactions on Automatic Control, 2016, 61(4): 921-936. DOI:10.1109/TAC.2015.2448491 |
| [10] |
Mei T. Blind multichannel identification based on Kalman filter and eigenvalue decomposition[J]. International Journal of Speech Technology, 2018, 22(1): 1-11. |
| [10] |
Jia L, Lou J, Yang Z. Blind adaptive identification of 2-channel systems using bias-compensated RLS algorithm[J]. International Journal of Adaptive Control & Signal Processing, 2018, 32(1): 301-315. |
| [11] |
Zhang Y, Asada H. Blind system identification of noncoprime multichannel systems and its application to noninvasive cardiovascular monitoring[J]. ASME Journal of Dynamic Systems:Measurement and Control, 2004, 126(4): 834-847. DOI:10.1115/1.1852460 |
| [12] |
Hahn J O, Reisner A T, Asada H H. Estimation of pulse transit time using two diametric blood pressure waveform measurements[J]. Medical Engineering and Physics, 2010, 32(7): 753-759. DOI:10.1016/j.medengphy.2010.04.017 |
| [13] |
Auricchio F, Conti M, Ferrara A, et al. A clinically applicable stochastic approach for noninvasive estimation of aortic stiffness using computed tomography data[J]. IEEE Transactions on Biomedical Engineering, 2015, 62(1): 176-187. DOI:10.1109/TBME.2014.2343673 |
| [14] |
Bazaral M G, Welch M, Golding L A R, et al. Comparison of brachial and radial arterial pressure monitoring in patients undergoing coronary artery bypass surgery[J]. Anesthesiology, 1990, 73(1): 38-45. DOI:10.1097/00000542-199007000-00007 |
| [15] |
Salvador S, Chan P. Toward accurate dynamic time warping in linear time and space[J]. Intelligent Data Analysis, 2007, 11(5): 561-580. DOI:10.3233/IDA-2007-11508 |