基于GA-BP神经网络进行呼吸运动预测的研究

黄志业 陈武凡 周凌宏 徐子海 陈超敏*（南方医科大学生物医学工程学院，广州 5055）

2（解放军303医院放射治疗中心，南宁 530021）

引言

在胸腹部肿瘤的放射治疗过程中，由于心跳、呼吸运动等的影响，使得肿瘤组织总是随着它们的运动而产生位置的偏移，这是提高胸腹部肿瘤放射治疗疗效的重要障碍。为了解决呼吸运动的影响，在传统的放射治疗和IMRT中通常采用肿瘤靶区射野扩边［1］、呼吸保持［2］、呼吸门控［3－4］等技术。这些技术的产生，对胸腹部肿瘤的放射治疗有着积极的意义。但是，人们很快发现，这些技术并不能很好地帮助人们达到期望的目的。扩边增加了肿瘤周围健康组织的照射体积，对周围的健康组织造成额外的损伤；呼吸保持方法要求病人在治疗过程中保持某一呼吸状态，对于那些有肺功能缺损的病人来说耐受性很差；门控技术要求放射线和病人呼吸周期同步，只在一个很小的呼吸期内打开射线进行照射，因此增加了治疗的时间。所以，人们正在积极地探索新的解决方法。目前，实时地跟踪和预测肿瘤运动轨迹，是研究补偿呼吸运动影响的重要方向，因为它可以克服扩边、呼吸保持和门控等技术的限制。赛博刀（cyberknife）系统就是实时跟踪技术在进行肿瘤放射治疗中的成功应用。

实时跟踪的最大挑战是治疗系统存在的延时，它必须在射线重新对准肿瘤的时间延时内预测出肿瘤的运动轨迹，并通过控制环把肿瘤的坐标发送到位置调整系统，重新调整射线和肿瘤之间的位置［5］。由于治疗系统的延时总是存在的，所以要实时地跟踪肿瘤的位置，需要提前预测出肿瘤在将来一段时间内的运动轨迹，以作为系统延时的补偿。本研究采用遗传算法（genetic algorithm，GA）和 BP（back propagation）神经网络结合，形成GA-BP神经网络算法，并把呼吸运动信号作为时间序列信号进行呼吸运动的预测。遗传算法（GA）是模拟自然演化过程搜索全局最优解的方法，它的群体性搜索策略使其不易陷入局部最优；同时，遗传算法基本上不用搜索空间的知识或其他辅助信息，而仅用适应度函数值来评估个体，在此基础上进行遗传操作。适应度函数不仅不受连续可微的约束，而且其定义域可以任意设定，使得遗传算法的应用范围大大扩展，方便与其他算法相结合［6－7］。BP（Back Propagation）神经网络算法本质上为梯度下降法，所以容易陷入局部最优，而且BP网络上每个节点的权值和阈值都会影响网络的最后输出，存在着收敛速度慢、学习过程中常常发生振荡等缺点。

鉴于遗传算法的优点和BP网络存在的缺陷，有学者就提出了将两者结合使用的设想，希望利用遗传算法不易陷入局部最优的特性，补偿BP神经网络梯度下降易于陷入局部最优的缺陷［8-10］。本研究基于这样的思想，对GA-BP神经网络算法在呼吸运动预测上做了初步的探讨。

1 GA-BP网络预测原理

利用GA-BP网络预测涉及两个重要的过程：首先应用遗传算法（GA），优化BP网络的权值和阈值的初始值；然后使用BP网络，根据获得的优化初始值进行网络训练。通过利用遗传算法能够搜索全局最优的优点，结合BP网络能进行任意维数输入输出映射的特性，达到提高网络预测精度的目的。

1.1 遗传算法优化BP神经网络的权值和阈值

在遗传算法中，首先要进行基因表述，即对遗传算法中的基因进行编码，因为基因包含了待求问题的所有参数。在遗传算法优化BP神经网络的权值和阈值中，遗传算法的基因包括了 BP神经网络各层间的权值和阈值。以3层BP网络的权值和阈值为例，可以用式子来表示遗传算法的一个基因，即式中，w和W分别表示神经网络输入层到隐层和隐层到输出层的权重，b和B分别表示神经网络隐层和输出层的阈值。

其次，设计合理的适应度函数，这是进行遗传算法选择运算操作的重要依据。进化论中的适应度，是表示某一个体对环境的适应能力，也表示该个体繁殖后代的能力。遗传算法的适应度函数也叫评价函数，是用来判断群体中个体优劣程度的指标，可根据所求问题的目标函数来进行评估。适应度函数的选取应具有规范性、合理性和通用性，可采用公式来加以定义，有

式中，ei表示第i个个体样本输入得到的实际输出和期望输出差的平方和。Tkj、okj分别表示第k个样本在第j个输出神经元的期望输出和实际输出，l和p分别表示训练样本数目和输出神经元数目，fi表示第i个个体的适应度值。

最后，实现遗传算法的选择、交叉和变异3个基本遗传操作。在简单遗传算法（SGA）中，选择的策略采用轮盘赌的方法，虽然每一个个体都可以获得复制的机会传到下一代，但并不能很好地体现优秀个体的竞争力和遗传算法优胜劣汰的基本准则。采用有条件的最佳保留策略，可以更好地体现优秀个体的竞争力，同时也能有效地防止算法的“早熟收敛”。

交叉操作是遗传算法中产生新个体的主要方式。在遗传算法中，交叉操作决定着算法的全局搜索性能。在算法搜索的前期，由于大部分个体分散在问题域中，应该强化算法的全局搜索求最优解，即在前期给予较大的概率Pc，从群体中选择个体进行交叉，形成新的个体，有利于算法的全局搜索，但是在算法的后期，假设个体的平均适应度达到了某个值，即认为群体中的大部分个体都聚集在全局最优解的附近，这时候如果仍然采用大概率进行交叉操作，则很可能破坏优秀的基因结构，不利于算法的稳定性。所以，根据具体的情况，自适应地改变交叉概率Pc的大小，即

式中，Pc表示交叉概率；fmax和favg分别表示群体个体的最大适应度值和平均适应度值。

变异操作的作用是保持群体中基因的多样性，偶然（Pm的概率很小）的变异在遗传算法中起着辅助和微调的作用，遗传算法中的变异操作决定着算法的局部搜索性能。基于以上的假设，在算法搜索的前期，由于大部分个体分散在问题域中，应该强化算法的全局搜索来求最优解，所以只需要很小的Pm对个体进行局部的微调；在算法的后期，假设个体的平均适应度达到了某个值，群体中的大部分个体都聚集在全局最优解的附近，这时候可以加大局部搜索的概率，提高局部搜索的效率。在局部进行个体的调整，使得在最优解附近的个体进一步靠近。所以，根据具体的情况，自适应地改变变异概率Pm的大小，有式中，Pm表示变异概率，fmax和favg分别表示群体个体的最大适应度值和平均适应度值。

上述遗传算法根据特定的网络结构编码，选择相应的适应度函数，自适应地改变交叉概率（Pc）和变异概率（Pm）大小，经过反复迭代地操作，最后输出优化的BP网络权值和阈值，并将此输出作为下一步BP网络训练的初始条件。

1.2 BP神经网络结构和预测

用遗传算法得到的输出作为BP神经网络的初始权重和阈值，因为遗传算法全局搜索的优点使得到的输出可以保证是全局最优解或接近最优解，而不同于直接使用BP神经网络，初始的权值和阈值是随机获得的。优化的初始权值和阈值可以帮助克服BP神经网络梯度下降容易陷入局部最优的缺点，使得BP网络依据优化初始权值和阈值在最优解附近进行搜索，以此提高网络的全局性能和网络的预测精度。

本研究采用3层BP神经网络，其拓扑结构如图1所示。该结构有3个输入神经元、1个输出神经元，由于BP神经网络对隐层神经元数目的确定尚无理论上的指导，所以按经验取6。图1所示的BP神经网络拓扑结构，首先把呼吸运动位移数据D分解成每三个连续数据作为一组神经网络输入序列S，然后利用S的每一列呼吸运动数据通过如图1所示的网络，在对应预期呼吸运动位移输出T的指导下进行网络权值和阈值的学习，最后形成BP神经网络输入输出的映射关系为

式中，λ表示预测超前的步数，采用30Hz的数据采样频率，在试验中分别用到100、300和500 ms的延时序列，所以在这些不同的预测延时序列中，λ分别取3、9和15。

图1 BP神经网络拓扑结构Fig.1 Topology of BP neural network

在图1中，S为式（6）表示的呼吸运动位移输入序列，圆圈和里面的符号表示各层网络节点和序号，w表示输入层神经元到隐层神经元的权重，W表示隐层神经元到输出神经元的权重，b和B分别表示隐层和输出层神经元的阈值，［DtDt－1Dt－2］表示一组 BP 神经网络输入，Dt+λ表示输入［DtDt－1Dt－2］序列得到的预期运动位移T中的输出值，λ表示预测超前的步数。

BP神经网络预测的基本思想是利用各个神经元的权值和阈值记忆序列的发展模式。在合适的样本下通过训练得到的权值和阈值在下次有相似的序列输入，可以得到对应预期的输出。在神经网络训练和预测中，样本的选取对最后的结果有非常重要的影响，一般要求训练的样本必须大于一个周期，并且能够体现总体数据的基本特征。训练样本数据选择越合理，训练后的神经网络权值和阈值包含越准确的信息，预测的精度也越高。

2 实验结果与分析

实验是在Windows XP操作平台、VC6.0和Matlab 7.6软件系统下进行的，用VC6.0对胸外标记物随呼吸产生运动位移的视频进行检测、跟踪，并记录标记物的运动坐标数据；利用Matlab 7.6对获得的运动坐标数据进行预处理，并使用了Matlab 7.6提供的神经网络工具箱（NNET Toolbox）进行GA-BP、BP神经网络的训练和仿真预测。图2（a）是所获得的其中一个样本数据归一化后的波形，由于在视频采集的过程中存在着一定的噪声，所以在实验的过程中还应该对实验数据进行平滑滤波，以获得更佳的实验效果，如图2（b）所示。

图2 样本数据滤波前后的波形。（a）滤波前的波形；（b）滤波后的波形Fig.2 Sample data waveform before and after filtering.（a）waveform before filtering；（b）waveform filtered

2.1 实验的基本步骤

步骤1：在受试者胸左下部贴一标记物，平躺在正常呼吸的情况下；设置摄像机采集频率为30Hz，拍摄一段标记物随呼吸产生运动位移的视频。然后，根据所得到的视频，在 VC6.0下利用帧间差分法的运动物体检测跟踪方法［11－12］，跟踪并记录呼吸过程中胸外标记物的运动坐标数据D。

步骤2：在Matlab7.6下，对获取的标记物运动坐标数据D进行预处理，包括归一化、平滑滤波。

步骤3：取已经过预处理的某一段数据（大于一个周期），并将其分解成如式（6）左边S所示的序列，作为GA-BP和BP网络（神经网络的拓扑结构如图1所示）的训练输入样本。在实验中，用100、300、500ms延时的预期输出如式（6）右边T所示，在采样频率为30Hz下 分别为3、9和15。对训练输入样本S进行训练，并获得相应的训练网络。

步骤4：利用所获得的数据和训练网络，分别进行GA-BP网络和纯BP网络的仿真预测。

步骤5：对比GA-BP网络和纯BP网络的预测结果。

2.2 实验结果及分析

实验对9个样本数据进行了测试，为了使比较有统一的标准，采用预测结果的均方误差（MSE）作为比较预测结果的准则。实验分别对GA-BP网络和BP网络在不同延时条件下的9个样本数据进行了预测，结果如表1所示。从表中可以看出，在相同的延时条件下，用GA-BP网络预测的结果比BP网络预测的结果更佳。而在同一种方法中，较小的延时具有更小的预测误差，这是因为预测的点与最相邻的点相关性最高，随着点之间的距离加大，它们的相关性也会逐渐减小［13］，致使长时间预测精度下降。

表1 GA-BP网络和BP网络的预测结果比较Tab.1 Comparison of GA-BP network's and BP network's predict error

图3给出了两个样本（样本1和样本5）在不同延时下两种预测算法的预测曲线。其中，（a）和（b）是在100ms下两个样本分别用两种预测算法得到的曲线，（c）和（d）是在300 ms下两个样本分别用两种预测算法得到的曲线，（e）和（f）是在500ms下两个样本分别用两种预测算法得到的曲线。

图4给出了两种算法在不同延时下样本平均预测误差的曲线。在图中可以更清楚地看出，GA-BP网络在不同的延时条件下对预测精度都有一定的提高，而且随着预测延时的变大，预测误差也逐渐变大。

图3 不同延时下的预测曲线。（a）样本1在100 ms时；（b）样本5在100 ms时；（c）样本1在300 ms时；（d）样本5在300 ms时；（e）样本 1在 500 ms时；（f）样本 5在 500 ms时Fig.3 Prediction curve under different delay condition.（a）sample 1 with 100 ms；（b）sample 5 with 100 ms；（c）sample 1 with 300 ms；（d）sample 5 with 300 ms；（e）sample 1 with 500 ms；（f）sample 5 with 500 ms

图4 不同延时条件下的样本平均预测误差。（a）100 ms；（b）300 100 ms；（c）500 100 msFig.4 Average prediction error underdifferent delay condition.（a）100 ms；（b）300 ms；（c）500 ms

3 结论

本研究对GA-BP网络结合的原理和优势进行了探讨。利用遗传算法（GA）良好的全局搜索能力，优化BP网络权值和阈值的初始值，从而克服BP网络容易陷入局部最优、收敛速度慢和学习过程中容易发生振荡的不足，使得BP网络在权值和阈值的训练学习中能更好地达到全局最优解。通过以上GA-BP网络预测和BP网络预测的实验结果分析，在相同的延时条件下，GA-BP神经网络预测呼吸运动相比于纯BP神经网络精度有所提高，特别是在较小延时条件下，GA-BP网络的预测结果优势更加明显。

上述预测是在获取外标记物运动坐标的基础上进行的，外标记物运动是随胸腹部在呼吸时运动而发生，与内部肿瘤的实际运动必然存在着差异。内部肿瘤的运动非常复杂而且差异性很大，其运动规律还需要进一步研究和探索。本研究对外标记物随呼吸过程的运动预测取得了相对满意的结果，如果在以后的工作中可以建立起外标记物运动和肿瘤运动之间的精确关系，那么上述方法在胸腹部肿瘤放射治疗中将会有更加实际的用途。

［1］Allen AM，Siracuse KM，Hayman JA，et al.Evaluation of the influence of breathing on the movement and modeling of lung tumors［ J］.International Journal of Radiation Oncology，Biology，Physics，2004，58（4）：1251 －1257.

［2］Treea A，Brocka J，McNaira H，et al.Improvement in tumor control probability with active breathing controland dose escalation：A modeling study［J］.Radiotherapy and Oncology，2009，91（3）：325－329.

［3］Minohara S，Kanai T，Endo M，et al.Respiratory gated irradiation system for heavy-ion radiotherapy［J］.International Journal of Radiation Oncology，2000，47（4）：1097－1103.

［4］Emery R，Rodriguez L，Barsa J，et al.Clinical experience using respiratory gated radiation therapy：comparison of freebreathing and breath-hold techniques［J］.International Journal of Radiation Oncology，2004，60（2）：419－426.

［5］Tchoupo G，Docef A.Nonlinear set membership time series prediction ofbreathing ［A］.In：Luciano Sbaiz， Olivier Cuisenaire，eds.Proceedings of 16th European Signal Processing Conference（EUSIPCO 2008）［C］.Lausanne：EURASIP，2008.25－29.

［6］雷英杰，张善文，李续武，等.MATLAB遗传算法工具箱及应用［M］，西安：西安科技大学出版社，2009.30－65.

［7］Man KF，Tang KS，Wong SK.Genetic algorithms：concepts and designs［M］.London：Springer，2001.

［8］Yang Yang，Li Kaiyang.Neural network based on GA-BP algorithm and its application in the protein secondary structure prediction ［J］.Chinese Journal of Biomedical Engineering，2006，15（1）：1－9.

［9］Belew RK，McInerney J，Schraudolph NN.Evolving networks：Using the genetic algorithm with connectionist learning［R］.CS90 －174，1990.

［10］Kitano H.Empirical studies on the speed of convergence of neural network training using genetic algorithms［R］.AAAI－90－118，1990.

［11］吴大鹏，程卫平，于盛林.基于帧间差分法和运动估计的Camshift目标跟踪算法［J］.光电工程，2010，37（1）：55－60.

［12］Kim Changick，Hwang Jenq-Neng.Fast and Automatic Video Object Segmentation and Tracking for Content-Based Applications［J］.IEEE Transactions on Circuits and Systems for Video Technology，2002，12（2）：122－129.

［13］Wu Huanmei， Betty Salzberg， Gregory C Sharp， etal.Subsequence matching on structured time series data［A］.In：Widom J，Chirkova R，eds.Proceedings of the ACM SIGMOD International Conference on Management of Data［C］.New York：ACM，2005，682－693.