预测模型，用方程表示为：
[0047] /)(Α')。(灰'-l)u'(A' -1)
[004引其中，參从)为Ξ维列向量，表示屯、脏表面目标点P在k时刻的Ξ维空间坐标预测值； Q化-1)是3 X 8维的模型设计矩阵，可表示为：
[0049] Q(k-l) = [p化-4) P化-3) P化-2) P化-1) hi化-1) h2化-1) h3化-1) 1]
[0050] 其中，P化-4)，p化-3)，p化-2)，p化-1)是屯、脏表面目标点P在k之前的4个历史测量 值山化-1)山化-l)，h3化-1)山化-1)是4个辅助点在k-1时刻的ii量值;W化-1)是k-1时刻 的模型参数，是由8个权值系数组成的列向量。
[0化1 ] S1.2、初始化GLM模型:令k= 1，对k= 1时刻之前的屯、脏表面目标点P和辅助点测量 值、模型参数及其方差矩阵进行初始化：
[0化2] p(-3)=p(-2)=p(-l)=p(0)=0
[0化3] hi(0) =h2(0) =h3(0) =0
[0化4] w(0)=0
[0055] V(0) = 1000018X8
[0化6] 其中，l8x康示8X8维的单位矩阵。
[0057] S2、判断k时刻屯、脏表面目标点P是否测量成功，若测量系统提供测量值，则测量成 功，获得该时刻屯、脏表面目标点P的测量值P化），然后执行步骤S3;反之测量失败，则执行步 骤S4。
[0化引 S3、更新GLM模型参数
[0059] 基于迭代最小二乘滤波(RLS)原理，利用当前获取的屯、脏表面目标点P的测量值P 化），更新模型参数W化)及其方差矩阵Wk)，具体包括如下步骤：
[0060] S3.1、计算模型当前的预测误差和计算增益矩阵：
[0061 ] 预测误差：6 =八^')-知又')
[006^ 增益矩阵:K = Wk-l)QT(k-l)[Al3x3+Q(k-l)Wk-l)QT(k-l)]
[0063] 其中遗忘因子λ = 〇. 98,13X3表示3X3维的单位矩阵；
[0064] S3.2、利用预测误差和增益矩阵更新模型参数及其方差矩阵 [00化]模型参数:W化）=w化-1)+Κ · e
[0066]方差矩阵:V化）=λ-?ν化-1 )-λ-?Κ9 化-1 )Wk-l)
[0067]待模型参数更新完毕，执行步骤S5。
[006引S4、基于当前GLM模型获得屯、脏表面目标点P的预测值，具体如下：
[0069] /WA) = y(A -i)W' -I)
[0070] 用预测值代替测量值，即令抑&)=参(^，然后执行步骤S5。
[0071] S5、更新GLM模型设计矩阵
[0072] 利用当前时刻的屯、脏表面目标点P的测量值P化)和4个辅助点测量值hi化），h2化）， h3化），h4化)更新模型设计矩阵，得到
[0073] Q(k) = [p化-3) P化-2) P化-1) p(k) hi化）h2化）h3化）1]
[0074] 然后执行步骤S6。
[00对 S6、将当前时刻值k加1，即：k = k+l，再返回步骤S2,进入下一时刻k+1的处理流程。
[0076] 图2是本发明实施例中屯、脏表面目标点和辅助点的位置示意图。
[0077] 在本实施例中，测量系统采用基于立体内窥镜的视觉测量系统，对屯、脏表面1个目 标点和周围3个辅助点进行测量，3个辅助点选取了特征较明显的像素点，相较目标点更容 易测量。实际应用中，对于辅助点也可W采用完全相同的预测方法(此时，原目标点和其它2 个辅助点则变为辅助点）。运种处理可W确保当辅助点在某一时刻无测量时，由其预测值替 代，对目标点的GLM预测模型进行有效更新。
[0078] 图3是本发明实施例中屯、脏表面目标点测量失败时的运动预测结果示意图。
[0079] 在本实施例中，目标点在第908帖~987帖，即k = 908~k = 987期间，因手术器械的 遮挡，无法获得测量值，而3个辅助点未受遮挡，仍可正常测量。图3中显示了测量系统获得 的目标点和辅助点的运动曲线，并在目标点运动曲线中显示了其在第908帖~987帖受遮挡 期间，利用本发明方法得到的预测值的曲线。在本实施例中，根据图1中的处理流程：
[0080] 在第908帖，测量系统无法提供目标点测量值时，经步骤S2判断后，执行步骤S4,基 于当前GLM模型进行预测，获得目标点预测值，具体方程为：
[0081 ] 扔撕轉二试907)1<9的）
[0082] 其中，Q(907)和w(907)是在第907帖正常更新后的模型设计矩阵和模型参数。接下 来，在执行步骤S5时，由于第908帖无测量值，故令=风9()8)，利用目标点的预测值和 辅助点的测量值111(908)，}12(908)，}13(908)，}14(908)更新设计矩阵，得到：
[0083] Q(908) = [p(905) p(906) p(907) p(908) hi(908) h2(908) hsOOS) 1]
[0084] 而模型参数未更新，仍为前一帖的值，即：w(908)=w(907)，其方差矩阵仍有：V (908) =V(907)。
[0085] 在第909帖，与第908帖一样，在步骤4,目标点的经GLM模型方程预测得 i>(9〇9# ρ (9齡）（，利用预测值(令p(9〇9) =對9〇9))和辅助点测量值hi(909)山(909)，h3 (909) ，h4(909)对模型参数矩阵更新，得到Q(909)，而模型参数仍未更新，有w(909)=w (908)和 V(909)=V(908)。
[0086] 直到，第987帖，测量系统恢复，重新获得目标点测量值p(987)，此时，在步骤S3中， 恢复对模型参数矩阵的正常更新。
[0087] 尽管上面对本发明说明性的【具体实施方式】进行了描述，W便于本技术领域的技术 人员理解本发明，但应该清楚，本发明不限于【具体实施方式】的范围，对本技术领域的普通技 术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，运些 变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。
【主权项】
1. 一种基于时空相关性的心脏表面目标点运动预测方法，其特征在于，包括以下步骤： (1) 、建立并初始化线性预测模型GLM (1.1) 、建立GLM模型：利用心脏表面目标点p的N个历史测量值和Μ个辅助点建立L=M+N +1阶的GLM模型，用方程表示为：其中，/^)为三维列向量，表示心脏表面目标点P在k时刻的三维空间坐标预测值;Q(k-1)是3 X L维的模型设计矩阵，可表示为：它是由心脏表面目标点P在k时刻之前的N个历史测量值:p(k-N)，p(k-N+l)，. . .，p(k-1)和Μ个辅助点点在k-1时刻的测量值111仏-1)，112仏-1)，...，1?仏-1)，以及一个元素全为1 的三维列向量组成;w (k-1)是k-1时刻的模型参数，是由L个权值系数组成的列向量； (1.2) 、初始化GLM模型:令k = 1，将k = 1之前的N个历史测量值和0时刻的Μ个辅助点测 量值都初始化为〇向量，即： ρ(卜 Ν)=ρ(2_Ν) =…= ρ(0)=0 hi(0) =h2(0) = ??? = hM(0) = 0 并将〇时刻的模型参数初始化为〇向量，即w(0)=0,将方差矩阵初始化为V(0) = 1000 OIlxl，其中Ilxl表示L X L维的单位矩阵； (2) 、判断k时刻心脏表面目标点p是否测量成功，若测量系统提供测量值，则测量成功， 获得该时刻心脏表面目标点P的测量值P(k)，然后执行步骤(3);反之测量失败，则执行步骤 (4)； (3) 、更新GLM模型参数 基于迭代最小二乘滤波(RLS)原理，利用当前获取的心脏表面目标点p的测量值p(k)， 更新模型参数w(k)及其方差矩阵V(k)，待更新完毕后，跳入步骤(5); (4) 、基于当前GLM模型获得心脏表面目标点p的预测值，具体如下：用预测值代替测量值，即令P(k)=p(k)，然后执行步骤(5); (5) 、更新GLM模型设计矩阵 利用当前时刻的心脏表面目标点P的测量值P (k)和Μ个辅助点的测量值111仏-1)，112仏-1)，. . .，hM(k-l)更新模型设计矩阵，得到(6) 、当前时刻值k加1，即：k = k+l，再返回步骤(2)，进入下一时刻k+Ι的处理流程。2. 根据权利要求1所述的基于时空相关性的心脏表面目标点运动预测方法，其特征在 于，所述的步骤(2)中，更新模型参数w(k)及其方差矩阵V(k)的具体方法为： (2.1 )、计算模型当前的预测误差和计算增益矩阵： 预测误差：=- ##) 增益矩阵:K = V(k~l)QT(k~l) [Al3x3+Q(k~l)V(k-l)QT(k-l)] 其中，λ为遗忘因子，I3X3代表3 X 3维的单位矩阵； (2.2)、利用预测误差和增益矩阵更新模型参数及其方差矩阵 模型参数:w(k)=w(k_l)+K · e 方差矩阵:ν(?〇=λ-^(k-D-A-kQU-DvU-i)。
【专利摘要】本发明公开了一种基于时空相关性的心脏表面目标点运动预测方法，利用心脏表面目标点运动轨迹的时间自相关性，以及与周围辅助点之间的空间相关性，建立预测目标点的线性模型GLM，再基于当前GLM模型下获得预测值，从而解决了长时间跨度预测时，由于时间自相关性减弱导致的预测误差急剧增大的问题。
【IPC分类】G06T7/20, G06T7/00
【公开号】CN105631864
【申请号】CN201510979738
【发明人】杨波, 郑文锋, 刘珊, 谢建军
【申请人】电子科技大学
【公开日】2016年6月1日
【申请日】2015年12月23日