帕金森病深部脑刺激手术中的图像配准问题研究

帕金森病深部脑刺激手术中的图像配准问题研究

ID:75352309

大小:5.59 MB

页数:92页

时间:2024-02-04

上传者:笑似︶ㄣ無奈
帕金森病深部脑刺激手术中的图像配准问题研究_第1页
帕金森病深部脑刺激手术中的图像配准问题研究_第2页
帕金森病深部脑刺激手术中的图像配准问题研究_第3页
帕金森病深部脑刺激手术中的图像配准问题研究_第4页
帕金森病深部脑刺激手术中的图像配准问题研究_第5页
帕金森病深部脑刺激手术中的图像配准问题研究_第6页
帕金森病深部脑刺激手术中的图像配准问题研究_第7页
帕金森病深部脑刺激手术中的图像配准问题研究_第8页
帕金森病深部脑刺激手术中的图像配准问题研究_第9页
帕金森病深部脑刺激手术中的图像配准问题研究_第10页
资源描述:

《帕金森病深部脑刺激手术中的图像配准问题研究》由会员上传分享,免费在线阅读,更多相关内容在学术论文-天天文库

学校代码:10286分类号:TP391密级:公开UDC:004.9学号:153839帕金森病深部脑刺激手术中的图像配准问题研究研究生姓名:倪杨阳导师姓名:罗守华副教授申请学位类别工学硕士学位授予单位东南大学一级学科名称生物医学工程论文答辩日期2018年6月6日二级学科名称生物医学工程学位授予日期2018年6月日答辩委员会主席沈启松评阅人沈启松周平2018年6月6日 硕士学位论文帕金森病深部脑刺激手术中的图像配准问题研究专业名称:生物医学工程研究生姓名:倪杨阳导师姓名:罗守华副教授本论文受国家重点研发计划(2017YFA0104302)和东南大学苏州纳米技术协同创新中心的支持。 RESEARCHONIMAGEREGISTRATIONFORDEEPBRAINSTIMULATIONSURGERYOFPARKINSON'SDISEASEAThesisSubmittedtoSoutheastUniversityFortheAcademicDegreeofMasterofEngineeringBYNIYang-yangSupervisedbyProf.LUOShou-huaSchoolofBiologicalScience&MedicalEngineeringSoutheastUniversityMay2018 东南大学学位论文独创性声明本人声明所呈交的学位论文是我个人在导师指导下进行的研究工作及取得的研究成果。尽我所知,除了文中特别加以标注和致谢的地方外,论文中不包含其他人已经发表或撰写过的研究成果,也不包含为获得东南大学或其它教育机构的学位或证书而使用过的材料。与我一同工作的同志对本研究所做的任何贡献均已在论文中作了明确的说明并表示了谢意。研究生签名:日期:东南大学学位论文使用授权声明东南大学、中国科学技术信息研究所、国家图书馆、《中国学术期刊(光盘版)》电子杂志社有限公司、万方数据电子出版社、北京万方数据股份有限公司有权保留本人所送交学位论文的复印件和电子文档,可以采用影印、缩印或其他复制手段保存论文。本人电子文档的内容和纸质论文的内容相一致。除在保密期内的保密论文外,允许论文被查阅和借阅,可以公布(包括以电子信息形式刊登)论文的全部内容或中、英文摘要等部分内容。论文的公布(包括以电子信息形式刊登)授权东南大学研究生院办理。研究生签名:导师签名:日期: 摘要摘要帕金森病是一种基底神经节功能障碍引起的神经病症。深部脑刺激(DBS)手术是帕金森后期患者治疗的有效手段。该手术通过对患者脑部植入电极,刺激丘脑底核靶点以实现治疗的目的。DBS手术需要对术前术后的靶点和电极位置进行精确定位,但在现有条件下,对患者脑部的成像存在如下问题:首先,丘脑底核作为靶点体积较小,只有在特定的T2WI或SWI加权的MR图像中可见,且靶点在图像中没有明确边界,难以定位;其次,DBS手术会引起患者脑组织的形变,造成靶点位移,加大术后定位难度;此外,植入电极为金属材质,限制了术后MR成像时间,使得图像分辨率较低,图像中只有电极信息而无靶点信息。这些因素导致实际的靶点和电极位置信息分布在不同模态、不同分辨率的图像上,很难对它们在脑中的位置及相互关系进行直观准确地呈现,给DBS手术计划制定及术后治疗带来困难。为解决上述问题,本文基于ICBM-152脑图谱,对手术前后图像进行了配准方法的研究。实际临床数据在手术过程中存在着模态的不同、空间位移的变化,与标准的人脑图谱之间也存在着坐标空间、尺寸、个体生物学局部特性的较大差异,因此很难使用已有的算法对双方进行配准。针对以上问题,本文所设计实现的脑图谱与临床图像的配准方法包括如下几个阶段,首先,针对临床数据与图谱坐标空间不匹配的问题,使用基于K-L距离的脑中膈面识别算法,获得脑中膈面的空间位置信息,将其坐标空间与图谱坐标空间进行标准化;其次,由于图谱在整体尺寸上与临床图像不同,采用FCM聚类算法提取脑皮质标识点,实现大脑分块映射及对应缩放比例的计算,并使用分段仿射配准方法,对图谱与临床图像进行粗配准;最后,为解决粗配准图谱与图像的局部失配问题,通过使用基于B样条函数的非刚性配准方法,并进行多尺度B样条自由形变优化,确定不同网格尺度下的优化参数,对临床数据与脑图谱进行精配准。上述研究结果虽然较好地实现了ICBM-152图谱与MR图像的配准,但实际DBS手术中,由于条件限制,术前T1WI图像为全脑数据,无靶点信息;术前T2WI图像以及术后图像分别包含靶点及电极信息但缺少标识点信息。因此,配准算法所需的信息分散在不同模态、不同图像质量的数据中。为保证对图像信息的有效利用,本文针对DBS手术图像建立了一个完整的数据处理流程管线,应用于帕金森病DBS计划环节时,实现了术前靶点定位;应用于术后效果评估的环节时,实现了术后的靶点和电极定位。上I 摘要述算法利用实际数据进行了测试,取得了较好的效果。关键词:帕金森病,深部脑刺激手术,丘脑底核,脑图谱,配准II AbstractAbstractParkinson'sdiseaseisaneurologicaldisordercausedbybasalgangliadysfunction.Thedeepbrainstimulation(DBS)surgeryisaneffectivetreatmentforpatientswithadvancedParkinson'sdisease,therapeuticpurposesarereachedbystimulatingthesubthalamicnucleusinthepatient'sbrainthroughtheimplantedelectrode.TheDBSsurgeryrequiresprecisepositioningatthepreoperativeandpostoperativesubthalamicnucleusandtheelectrode.However,undertheexistingsurgicalconditions,itistechnicallychallengingfortheclinicalimagingofpatient’sbrain.Atfirst,thesubthalamicnucleusistinyandonlyvisibleinT2WIorSWIMRimages,thefuzzyboundariesoftheimagesmakesthemdifficulttobelocated.Secondly,theDBSsurgerymayleadtothedeformationofbraintissue,thedisplacementofnucleuswillincreasethedifficultyofpostoperativepositioning.Finally,theimagingtimeandresolutionofthepostoperativeMRimagesarelimitedbytheimplantedmetalelectrode.Therefore,undermultipleimagingmodalitiesduringtheoperation,theactualpositionsofthenucleusandelectrodearedifficulttodiscriminatevisuallyandaccurately,theriskinpreoperativeplanningandpostoperativetreatmentofDBSsurgeryareincreased.Tosolvetheproblemsabove,theregistrationmethodsofpreoperativeandpostoperativeimages,whicharebasedontheICBM-152brainatlas,arestudiedinthisthesis.Sincethemodalitiesandspatialdisplacementoftheclinicalimagingarechangingatvariousstagesofsurgery,thedifferencescomparedwiththestandardhumanbrainatlasexistinthecoordinatespace,sizeandlocalcharacteristicsofindividualbiology,whichmakeithardtoregisterthemusingexistingalgorithms.Inthispaper,weproposeacomprehensiveregistrationalgorithmfrombrainatlastoclinicalimages.Tobemorespecific,firstly,aimingattheproblemofthecoordinatespacemismatchbetweentheclinicaldataandatlas,aK-Ldistancerecognitionalgorithmisemployedtoobtainthespatialpositioninformationofthemid-sagittalplaneinclinicalbraindata,thecoordinatespacesofthemarethenstandardized;secondly,sincetheanatomyofclinicaldataisdifferentwiththatoftheatlasinoverallsize,aFCMclusteringalgorithmisservedastoextractthemarkersofcerebralcortex,thenapiecewiseaffineregistrationmethodisexecutedtoachievethecoarseregistrationoftheclinicalimagesandtheIII Abstractatlasbymappingthebrainblocksandcalculatingthecorrespondingscalingratios;finally,anon-rigidregistrationalgorithmbasedonB-splinefunction,alongwithamulti-scaleB-splinefreedeformationoptimizedingridparameters,isemployedtoachievethefineregistrationresultsbetweentheclinicaldataandtheatlas.TheresultsofthestudyabovehavesuccessfullyworkedontheregistrationoftheICBM-152atlasandclinicalMRimages.However,duetotheimagingconditionsintheDBSsurgery,thepreoperativeT1WIimagescontainwholebraindatawithoutnucleusinformation;thepreoperativeT2WIimagesandpostoperativeimagescontainnucleusandelectrodeinformationbutarelackofcorticalmarkers,theinformationrequiredfortheregistrationalgorithmisdistributedamongdataindifferentmodalitiesandimagequality.Tomakeuseoftheinformationeffectivelyandefficiently,acompletedataprocessingpipelineofsurgicalplanningandperformanceassessmentisdesigned,whichisappliedtoachievenucleuslocationforthepreoperativeplanofDBSsurgery,aswellasthenucleusandelectrodepositionforthepostoperativeeffectassessment.Thetestofthepipelineusingtheactualdatashowsgoodresults.Keywords:Parkinson'sdisease,deepbrainstimulationsurgery,subthalamicnucleus,brainatlas,registrationIV 目录目录摘要..................................................................................................................................................................................IAbstract........................................................................................................................................................................III目录................................................................................................................................................................................V第一章绪论..........................................................................................................................................................11.1研究目的及意义............................................................................................................................................11.1.1帕金森病..................................................................................................................................................11.1.2帕金森病脑深部电刺激手术.............................................................................................................11.1.3DBS手术流程..........................................................................................................................................31.1.4帕金森病DBS手术的难点..................................................................................................................51.2DBS术前靶点定位及术后电极靶点定位研究状况.............................................................................71.2.1DBS术前靶点识别技术研究现状.....................................................................................................71.2.2DBS术后电极靶点定位技术研究现状............................................................................................81.3本文研究内容及难点...................................................................................................................................91.4本文组织结构...............................................................................................................................................12第二章基于ICBM-152图谱的临床数据空间标准化........................................................................142.1数字化人脑图谱研究现状........................................................................................................................142.1.1Brainweb图谱.......................................................................................................................................142.1.2WholeBrain图谱.................................................................................................................................152.1.3Talairach图谱........................................................................................................................................152.1.4MNI脑模板............................................................................................................................................162.2ICBM-152图谱空间....................................................................................................................................172.3基于K-L距离的脑中膈面识别算法.....................................................................................................172.3.1K-L距离的定义....................................................................................................................................182.3.2问题描述................................................................................................................................................182.3.3算法设计................................................................................................................................................202.4临床数据校准至图谱空间........................................................................................................................212.5实验结果与分析..........................................................................................................................................232.5.1模拟数据结果.......................................................................................................................................242.5.2实际数据结果.......................................................................................................................................26V 东南大学硕士学位论文2.6本章小结........................................................................................................................................................29第三章基于脑皮质标识点的分段仿射配准............................................................................................303.1大脑分区以及脑皮层标识点...................................................................................................................303.1.1大脑分区................................................................................................................................................303.1.2脑皮质层标识点..................................................................................................................................313.2脑皮层标识点的提取.................................................................................................................................323.2.1FCM聚类算法.......................................................................................................................................333.2.2FCM聚类以及图像二值化................................................................................................................363.2.3形态学图像处理..................................................................................................................................373.2.4实验结果与分析..................................................................................................................................403.3基于标识点的大脑分段仿射变换匹配.................................................................................................433.3.1仿射变换的原理与性质....................................................................................................................433.3.2计算分块比例.......................................................................................................................................443.3.3实验结果与分析..................................................................................................................................453.5本章小结........................................................................................................................................................47第四章基于B样条函数的非刚性配准.....................................................................................................484.1基于B样条函数的非刚性配准算法.....................................................................................................484.1.1互信息测度...........................................................................................................................................484.1.2B样条函数自由形变模型..................................................................................................................504.1.3L-BFGS优化算法.................................................................................................................................514.1.4实验结果与分析..................................................................................................................................524.2多层B样条自由形变配准........................................................................................................................594.2.1不同网格尺度的形变特点................................................................................................................594.2.2多层B样条自由形变优化................................................................................................................604.2.3实验结果与分析..................................................................................................................................604.3本章小结........................................................................................................................................................65第五章面向帕金森深部脑刺激手术的系统应用...................................................................................675.1临床数据及扫描参数.................................................................................................................................675.2数据处理流程设计......................................................................................................................................685.3数据处理流程的应用..................................................................................................................................695.4本章小结........................................................................................................................................................74第六章总结与展望..........................................................................................................................................75VI 目录6.1总结.................................................................................................................................................................756.2展望.................................................................................................................................................................76致谢..............................................................................................................................................................................77参考文献.....................................................................................................................................................................78作者简介.....................................................................................................................................................................81VII 第一章绪论第一章绪论1.1研究目的及意义1.1.1帕金森病帕金森病(Parkinson’sDisease,PD),又称震颤麻痹,是一种基底神经节功能障碍引起的神经病症,因黑质多巴胺能神经元退行性变,纹状体多巴胺能递质减少而发病[1]。在1817年,由英国医生JamesParkinson首次完整描述了震颤麻痹的临床特征而得名。帕金森病患者的临床症状通常表现为动作迟缓,静止性震颤,姿势平衡障碍和肌张力增高等,同时伴有精神方面的抑郁、焦虑、幻觉等问题[2-4]。帕金森病发病较为隐蔽,进展缓慢,在漫长的发展过程中会表现出上述不同的症状。目前,中国有近230万帕金森病患者,占全球患者约50%,且患者的数量正在以每年10万的速度增加[5],随着我国人口老龄化程度的加剧,该数字还在急剧上升。到2030年,我国帕金森病人预计将会达到五百万人,占全球患者数的57%[6]。因此,我国帕金森病的防治将面临严峻考验。1.1.2帕金森病脑深部电刺激手术左旋多巴(L-Dopa)替代治疗一直是帕金森病的最有效的药物治疗方案[7],核团捣毁手术也是治疗帕金森病的方式之一。左旋多巴可以增加纹状体中的多巴胺递质浓度,从而明显改善患者的症状。大部分患者都会进行左旋多巴替代治疗,并将其作为控制帕金森病情的主要手段,是帕金森病治疗的金标准[8]。但是随着症状的加重,用药量也需要逐渐增加,而疗效会逐渐降低,并伴有药物过量带来的“开-关现象”、“剂末现象”、“晨僵少动”、“运动徐缓”等并发症[9],部分帕金森后期病人会实施立体定向苍白球丘脑损毁手术进行治疗。然而损毁手术会有30%~40%的几率导致偏瘫、言语障碍等并发症[10],此外,手术造成的创伤及副作用较大,这些因素导致该手术的应用较为受限。脑深部电刺激手术(DeepBrainStimulation,DBS)是帕金森病外科手术上的巨大突破,对于长期服用帕金森病治疗药物(如左旋多巴等)出现疗效减退的患者是有效的治疗手段,该手术被用于改善包括震颤、僵硬、运动迟缓、姿势不稳、特发性震颤(ET)和肌张力障碍等帕金森病症状,是目前帕金森病后期首选的治疗方法。DBS系统通过将颅内电极植入大脑深部的目标区域,并使用皮下导线连接到植入脉冲发生器(IPG),该1 东南大学硕士学位论文脉冲发生器通过电极提供电刺激,从而达到治疗效果。颅内电极为四级电极,电极分为间距1.5毫米的3387型,间距为0.5毫米的3389型等。植入脉冲发生器可以被外部编程器控制,从而释放不同电压、脉宽以及频率的脉冲信号。图1-1Medtronic3387以及3389型电极示意图[11]目前DBS系统可以选择不同的手术植入靶点解决相应神经障碍疾病。其中,以内侧苍白球(GPi)作为植入靶点可以治疗肌张力障碍;以丘脑腹侧中间核(Vim)作为植入靶点可以治疗特发性震颤;本文所探讨帕金森病则通常选择丘脑底核(STN)作为植入靶点,其解剖示意图见图1-2。图1-2脑深部核团解剖示意图[12]DBS手术的疗效通常会与核团损毁术进行对比,核团损毁术的手术靶点主要为Vim核以及Gpi核。在手术效果上,Vim核损毁术可以显著改善患者的震颤以及僵直症状,但是对于运动减少症状的治疗效果较差;Gpi核损毁术对于僵直、运动减少症状的治疗效果较好,而对于震颤的治疗效果较差;DBS对于震颤、僵直以及运动减少这三种症状均具有较好的效果[13]。从中远期效果来看,双侧DBS治疗可以维持5年的治疗效果[14],其疗效更加持久稳定,相比损毁术有明显优势。在手术安全性上,DBS损伤的范围较小,2 第一章绪论且造成的损伤是可逆的,安全性较高。对于手术便捷性,DBS电极可以同时实施双侧植入,而核团损毁术通常只能进行单侧手术。此外,DBS治疗的可控性较强,所植入的电极的多个放电点可以任意组合,且可以调节放电电压、频率、脉宽等参数,达到控制症状且不产生副作用的目的[15,16]。因此,DBS手术是目前帕金森后期患者控制帕金森症状最有效的治疗手段。虽然DBS对于帕金森病以及特发性震颤的治疗有良好的效果,但是临床上对脑深部电刺激手术(DBS)作用机制的认识还相对匮乏[17]。1.1.3DBS手术流程DBS手术流程图如图1-3所示,主要分为术前评估、MR全脑扫描、安装立体定向头架及其后续CT成像、DBS植入手术和术后程控等5个过程。术前评估主要用于评估患者的症状严重程度以及其身体状况是否满足实施手术的条件;对患者进行术前MR全脑扫描以及安装立体定向头架后的CT扫描,用以准备手术计划以及确定电极植入路径所需的数据;DBS植入手术后,待患者经过一段时间的静养,对整个系统进行开机以及刺激参数的调节,并实施术后程控,达到治疗的目的。各个过程的实施细节如下:图1-3DBS手术流程图(1)术前评估脑深部电刺激系统在植入前,需要对患者进行术前评估,全面地评价患者症状的严重程度,判断其是否满足进行DBS手术治疗的条件。评估内容对DBS术后的疗效预测,3 东南大学硕士学位论文以及术中和术后程控过程中刺激参数调节等也有参考作用。一般采用帕金森综合评分量表[18](UnifiedParkinson’sDiseaseRatingScale,UPDRS)评分及多巴胺冲击实验。(2)获取患者头部MR信息:在术前对患者进行MR全脑扫描,获得患者T1以及T2加权的全脑矢状位以及横断位MR图像。这些图像主要用于获取患者头部的标识线,前连合与后连合(Anteriorcommissure-posteriorcommissureline,AC-PC)之间的连线,以及术前植入部位的初步评估,确定头架的基准线以及切口线。图1-4ACPC解剖示意图[19](左)以及临床数据中的AC-PC连线(右)(3)立体定向头架的安装立体定向头架是排除术中抖动等干扰因素,确保电极稳定并精确植入的金属框架。在其安装前,首先需要使用乙醇对患者头部消毒,再对患者进行局部麻醉,随后戴上立体定向头架。在安装时,头架要与患者的头部标识以及基准线重合,前方螺钉位于AC-PC线上方,后方螺钉位于AC-PC线下方。随后对带立体定向头架[20,21]的患者进行CT扫描,将MR与CT图像进行融合,制定精确的手术计划,确认植入靶点在立体定向头架下的坐标。(4)DBS植入手术对患者实施局部麻醉并开颅,使用单极电极套针管刺破硬脑膜,以明胶海绵及棉片保护骨孔。随后植入微电极,固定并连接导线,从植入靶点的上方10mm至下方5mm推进微电极,并记录不同位置神经元细胞的放电情况。根据记录的结果确定电极最终的植入位置,将微电极退回至植入点位置。随后将刺激电极植入最终的记录位置并固定,选择适当的刺激参数观察患者上肢以及下肢的震颤、运动灵活性以及肌张力的改善情况,同时嘱咐患者进行一些规定动作,询问患者是否出现身体麻痹等不良反应,如患者需要进行双侧植入,则在另一侧进行相同的操作。完成植入后关闭双侧切口,对患者进行全4 第一章绪论身麻醉并进行脉冲发射器植入手术。图1-5DBS手术植入效果示意[22](5)系统开机以及术后程控术后对患者进行CT或者低剂量MR检查,评估电极位置,结合术中电生理记录以及术中病人反馈判断电极的最佳触点。在DBS系统开机前一天,患者需要停用抗帕金森病药物,防止药物对程控过程的干扰。在程控过程中,通常对症状较重的身体一侧优先开机,并逐步提升电压,观察患者症状改善情况、是否出现副反应等,最后选择对症状控制好且副反应的电压值作为最高电压阈值点。1.1.4帕金森病DBS手术的难点DBS治疗成功的关键是术中电极植入的精确度,以及术后根据电极的植入情况,调节刺激参数的准确度。为了保证以上两个过程的精度,需要克服一系列的难题。1.1.4.1DBS手术靶点定位DBS手术植入的靶点为STN核,成年人的丘脑底核大小约为80𝑚𝑚𝑚𝑚3,直径约5~8mm,体积较小。在MR图像中,只能通过T2WI加权或SWI加权成像观察,且在图像中其灰度与周围组织相近,成像结果通常没有明确的边界,因此很难直接使用算法将其从周围组织中分割出来。而帕金森患者对身体的控制力较差,在MR成像过程中容易晃动,更降低了图像的质量。5 东南大学硕士学位论文图1-6术前T2WI图像中的丘脑底核在手术进行时,由于患者颅腔被打开,部分空气进入导致颅内压力变化,引起脑部形变;在双侧电极植入过程中,第一根电极在植入过程中同样会使脑组织发生形变,影像第二根电极的植入精度。因此,这两个因素容易使得术前所规划的电极植入路径失效,从而影响手术的效果。1.1.4.2DBS术后靶点与电极位置关系的确定脑深部电刺激手术后,为了防止患者在高场强磁共振下被颅内的金属电极灼伤,通常会对患者头部使用CT或者1.5T场强MR进行成像,获得金属电极在患者脑部的位置。术后MR为了减少成像时间,通常不会对患者进行全脑扫描,只对丘脑底核部分进行薄层扫描。这导致术后数据在分辨率的下降及大量的信息缺失,只能在术后图像中获得电极在颅内的位置,无法直接观察核团靶点与电极之间的相互关系。图1-7术后T1WI图像中的电极轨迹手术所植入的电极通常是Medtronic所生产的3387型或3389型电极(图1-1)。在术后程控过程中,医生依据电极与核团的相互关系,调节电极尖端四个金属触点的刺激6 第一章绪论参数,来缓解PD患者的症状。如果术后的MR无法提供二者相互位置关系,就会给术后的程控带来较大的困难。一旦选择错误的刺激参数,往往会给病人带来痛苦,严重时还会导致失明、失语、行动障碍等并发症。1.2DBS术前靶点定位及术后电极靶点定位研究状况综上可知,利用术前图像对DBS的植入点进行定位是整个DBS手术的基础,也是决定手术效果的根本环节。因此,如何在图像中获得患者适合的植入靶点信息成为了被广泛关注的问题。而术后程控的效果直接决定了DBS手术的成败,因此对于术后的靶点与电极相互位置关系的精确计算也是一个重要的研究课题。国内外研究者已经基于DBS手术前后获得的影像数据,对上述两个问题进开展了大量的研究工作。1.2.1DBS术前靶点识别技术研究现状在国际上,Villeger等人[21]通过基于模糊逻辑的图谱数据融合方法可以识别MRI中的大脑结构,其中每个信息源(即患者数据和专家知识)在相同的可能性框架中被模型化。模糊集合和可能性理论框架适合解释在活组织和解剖结构中观察到的个体间变异性,特别是对丘脑底核的分割获得了初步成果。YuanLiu等人[23]提出了一种基于多模式学习的方法,使用回归森林来自动定位术前MR脑部扫描中的目标。通过利用患者的大数据集,该方法通过将原始估计值调整到患者的个体解剖结构来改善间接定位结果,并且利用多对比度信息来结合不同直接定位方法,该方法与现有的基于统计图谱方法的精度相似。Castro等人[24]比较了专家手动定位与基于图谱的算法以及基于非刚性配准的算法的精度,并进行了交叉验证,认为两种方法的精确度相近,但是由于实验采用的数据都是清晰可见的STN进行T2加权的图像,因此样本空间较小;T2图像到T1空间的刚性变换也会引入误差;基于AC-PC图谱的算法需要手动进行坐标点识别,同样会引入误差。Walter等人[25]使用经颅超声检查(TCS),实现了术中监测DBS电极,以及高分辨率的实时位置成像,高端TCS在立体定向手术中,结合X射线图像,可以精确记录电极的正确位置。Pouratian等人[26]使用了基于概率纤维束追踪成像定位DBS手术的靶点,这种新颖的丘脑靶向方法是基于识别与前运动和辅助运动皮质连接的可能性最高的丘脑区域。该方法比传统的间接靶向术前方法具有优势,提供了可提高深部脑刺激手术的精确度、功效和效率的患者特异性靶点。7 东南大学硕士学位论文在国内,南京医科大学第一附属医院使用神经导航技术以及图像融合技术实现对术前的手术计划[27],有效降低了术中微电极出血,减少了脑脊液的丢失,增加了手术的可靠性。耿馨佚等人[28]总结了临床上成熟的DBS电极靶点定位方法,基于影像解剖学定位以及微电极记录功能定位,影像学定位方法存在个体性差异较大的问题,且受选用的图谱的影响较大。功能定位法直接反映脑功能状态,以及靶点区功能信息,两者结合可以提高定位的准确性。1.2.2DBS术后电极靶点定位技术研究现状在脑深部电刺激手术过程中,由于开颅手术造成的脑转移等各方面因素的影响,电极往往无法精确的按照手术计划阶段设定的路径进行植入,从而影响治疗效果,因此对于DBS术后的植入效果评估的研究也具有深刻的意义。Horn等人[22]为了确定深部脑刺激(DBS)手术后电极的放置位置,开发了一种便于重建铅电极轨迹和接触位置的新型工具箱。该工具箱通过将术前术后影像数据向脑图谱配准,以及一些半自动算法实现DBS电极导线以及触点的精确快速重建。Videen等人[29]开发了一个识别靠近丘脑底核结构基准点的算法,基于这些基准点,将术前MR图像以及术后CT图像精确配准至Mai脑图谱,同时根据术后CT图像中的电极伪影定位电极中心。Silva等人[30]提出了一个系统,该系统使用非刚性配准算法将患者数据配准至标准图谱,基于深部大脑结构的分割和多模式成像方法,准确地获得DBS电极相对于解剖结构的三维位置,减少了临床医生的工作量。在国内,第二军医大学郝斌使用手动标记方法将丘脑底核标记,实现术后丘脑底核区域的可视化,根据术后MRI数据作为判断电极位置的依据,并根据两者的位置关系探究了丘脑底核靶点位置关系与DBS治疗帕金森病的疗效之间的关系[31]。苏州生物医学工程技术研究所张芷齐等人[32]在标准脑图谱空间下将双侧丘脑底核刺激治疗患者的高分辨率术前图像与低分辨率的术后图像进行配准并融合,同时使用阈值法确定横断面图像上电极位置,重建脑深部电刺激电极的路径,并使用图谱信息分割出丘脑底核等核团,实现核团与电极的三维空间可视化,为医生的术后程控提供辅助指导。综上所述,在脑深部电刺激靶点定位过程中,主要是使用基于医学影像学的CT/MRI配准融合算法实现,该方法是临床公认的判断DBS手术靶点位置的有效途径。部分情况下会结合超声等其他模态的临床图像,或基于神经生理学的功能定位法,比如微电极8 第一章绪论记录、局部场电位等手段进行联合定位,以增加精确度。但是这些研究中,往往使用图像质量较好的非DBS手术图像,基于临床DBS手术MR图像的研究较少。对于实现DBS术后靶点与电极的相互位置关系的确定,不但需要依靠术后图像,也需要术前图像的参与。由于目前脑深部核团成像技术的限制,在DBS术后效果评估算法中,使用脑图谱辅助进行丘脑底核以及其附近核团的识别与分割是较为常用的方法。但是脑图谱与实际人脑的差异较大,对于配准算法提出了较高的要求。上述几种基于图谱的算法都采用将患者术前术后的图像向标准的脑图谱进行配准,虽然在一定程度上解决了丘脑底核靶点与电极的定位问题,但是在配准过程中容易导致电极产生扭曲的现象,从而在结果中引入误差,并不能准确体现二者的位置关系。1.3本文研究内容及难点对于DBS手术靶点定位和术后程控治疗的良好实现,需要精确定位核团位置,然而核团在术前术后的图像中都较难直接进行分割得到。结合研究可以发现,脑图谱对于脑部组织结构具有较好的指示性作用,因此对脑图谱和术前术后图像进行配准,并利用脑图谱已标注的核团信息进行位置的精确定位,是比较可行的技术路线。在实施细节上,需要解决图谱与临床数据空间位置、实际大小以及内部纹理等的不一致性的问题,因此脑图谱与临床数据的精确配准是其中最为关键的一环,它是术前靶点识别以及靶点与电极的相互位置关系确定的基础。脑图谱只是一个标准的人脑模板,虽然在大致结构上与临床数据做到了“形似”,但是二者在很多方面还存在着差异。如图1-8所示,从全局来看,人脑图谱的大脑尺寸略大于临床数据中的大脑,图谱数据中的大脑形态相较于临床数据更为狭长,且图谱数据所摆放的位置为标准的AC-PC连线为水平的位置,而临床数据位置较为随机;从局部来看,图谱数据的结构相对于脑中膈面是对称的,而临床数据的左右两侧并不对称。因此对于二者的配准过程,需要同时考虑全局以及局部的差异情况,直接使用刚性或非刚性配准算法无法实现二者的有效配准。9 东南大学硕士学位论文图1-8临床数据域图谱在三个标准方向上的比较本文通过调研确定了一种脑图谱与脑图像的配准算法,根据帕金森影像数据的特点确定了算法的执行步骤。配准流程如图1-9所示:图1-9图谱与临床数据的配准流程图首先将术前临床图像标准化至图谱的坐标空间,使临床数据与图谱的坐标空间匹配;10 第一章绪论随后使用对图谱与临床图像进行粗配准,实现二者全局大小的统一;最后对二者进行精配准,达到图谱与临床数据在细节上的统一。在实现图谱与临床图像的配准后,通过图谱所对应的子图谱可以得到脑深部核团的位置,丘脑底核在图谱数据中的标记如图1-10所示。图1-10丘脑底核在图谱中的标记图1-11DBS术前术后数据在三个标准方向上的部分断面在结合脑图谱与DBS图像数据,获得手术的电极植入效果的过程中,由于DBS手术成像条件的限制,单一模态的临床数据所包含的信息不全,信息分散在不同时间及模态的图像中。图1-11为DBS手术前后的不同模态图像数据的断面,从中可以看出,术11 东南大学硕士学位论文前T2WI以及术后图像中的信息大量缺失,而这两者中分别包含了丘脑底核靶点以及电极的位置信息,术前T1WI数据是完整的全脑数据,提供了图谱配准过程中的一些全局配准信息,但是缺少丘脑底核靶点以及电极信息,此外由于手术过程中的环境变化,术后患者脑部往往会发生形变。因此,要了解术后的丘脑底核与电极的相互位置关系,需要考虑术后患者大脑形变的因素,跨模态利用图像中各自的有效信息。然而,目前基于临床图像的丘脑底核以及DBS电极的相互位置关系的应用还不成熟。因此对于DBS手术电极与核团相互位置关系的确定,需要结合手术前后图像特点,在实现图谱与图像配准的前提下,确定相应的数据处理流程。根据对现有文献的调研,本课题的研究难点如下:1、临床数据坐标空间标准化问题。临床数据由于在成像过程中受众多随机因素影响,导致其坐标空间与图谱的坐标空间不同。一般来说,脑图谱的坐标空间是符合解剖学规范的,应作为基准的坐标空间。因此,通过一定的评判标准以及算法,实现临床数据至图谱的坐标空间标准化是本文需要解决的一个难点。2、脑图谱与临床数据配准问题。该问题的解决是实现脑组织的分割[33]以及丘脑底核靶点识别的基础。临床中的医学数据配准通常涉及到的是同一患者的单模态或多模态数据配准,如CT-MR,MR-超声等之间的配准,由于均为患者本人的影像数据,通常采用简单的全局刚性或者仿射变换就可以取得较好的效果。脑图谱是用于展示人类脑组织结构的模拟数据,与患者的脑部图像在分辨率、大小、形状以及局部细节等方面都有着较大的差异,单纯的刚性或者非刚性配准算法都无法满足脑图谱与临床数据配准的要求。因此,实脑图谱与DBS临床数据精确配准具有较大难度。3、临床数据中有效信息的使用效率问题。由于成像技术限制,术前对患者进行T2WI加权的MR图像为核团部位的薄层数据,导致图像中缺失全脑信息。在脑深部电刺激手术后,为了防止患者在高场强磁共振下被灼伤,医生会只对核团部分进行薄层扫描以减少成像时间,导致图像中缺失靶点信息。术前T2WI与术后数据在分辨率以及包含的信息存在大量缺失,因此如何结合术前T1WI数据的全脑信息,对术前术后的临床数据做到跨模态高效使用也是本课题的重点。1.4本文组织结构论文共分为六个章节,具体内容安排如下:12 第一章绪论论文第一章为绪论,主要介绍帕金森病深部脑刺激手术中配准问题研究的背景和意义,以及国内外的文献调研情况,梳理本文的难点以及章节安排。论文第二章介绍所使用的脑图谱以及对临床数据进行空间标准化,是后续配准工作开展的前提。论文第三章主要介绍脑皮质层标识点的提取以及在此基础上脑图谱与临床图像的分段仿射配准算法,实现脑图谱与脑图像的粗配准。论文第四章主要介绍基于B样条函数的非刚性配准算法,以及对算法进行的多层B样条优化。实现脑图谱与脑图像的精确配准。论文第五章主要结合以上算法及数据特点,以江苏省人民医院使用的DBS手术前后数据为例,构建了一个面向DBS手术的系统应用,实现了术前靶点识别以及术后电极与核团相互位置关系的获取。论文第六章对课题的工作进行了总结与展望。13 第二章基于ICBM-152图谱的临床数据空间标准化第二章基于ICBM-152图谱的临床数据空间标准化对于脑图谱与DBS手术图像的配准,首先需要解决在成像过程中由于患者的头部摆放位置的随机性,造成临床数据与标准图谱坐标系统不同中的问题。通常的做法是使用基于脑室的坐标系统[34],将临床数据变换至标准脑图谱的坐标空间,使二者的坐标空间统一。本章首先对数字化人脑图谱的研究进行了简单的回顾,对本文使用的ICBM-152图谱标准空间的概念进行了介绍,并使用基于Kullback-Leibler度量(K-L距离)的方法提取临床数据脑中膈面,将其校正至图谱空间矢状位垂直位置,然后将临床数据AC-PC校准至图谱空间水平位置,统一临床数据坐标空间与ICBM-152图谱的标准空间,为下一步配准工作打下基础2.1数字化人脑图谱研究现状目前国内外学者已经研究并开发出了多种数字化人脑图谱,比较具有代表性的图谱为加拿大McGill大学开发的Brainweb图谱,美国哈佛医学院开发的WholeBrain图谱Talairach-Tournoux立体序列医学脑图谱(Talairach脑图谱),蒙特利尔神经科学研究所开发的MNI脑模板等,每种图谱都各自的特点,下面将分别介绍。2.1.1Brainweb图谱Brainweb图谱[35](图2-1)是一个模糊表达解剖知识的脑模型,由10个数据集构成,该图谱对不同组织的空间灰度分布进行了定义,数据集的像素点的强度表示该像素值对应的组织在脑中的成分占比。该图谱为多模图谱,可以作为测试多模态数据配准算法以及分割标注算法的标准。图2-1Brainweb图谱14 东南大学硕士学位论文2.1.2WholeBrain图谱WholeBrain图谱[36](图2-2)由美国哈佛医学院开发,它不但包含了正常人脑的图谱,同时也包含了脑血管疾病、肿瘤、变性疾病、炎症以及传染病相关的病例。目前整个数据库中已经有超过13000个脑图像,其中包括了许多临床的病例,每个病例包含了完整的体数据,使用多模影像(MR、CT和PET)表达。同一病例的不同模态图像间已经完成了配准,可以实现不同模态、不同时间点之间的时域及空域比较。图2-2WholeBrain图谱2.1.3Talairach图谱Talairach脑图谱[37]由两位法国研究员TalairachJ和TournouxP共同开发。使用的数据来源为一个59岁法国妇女尸体的脑标本切片。图谱通过将脑标本沿矢状位方向间隔2~5mm进行切片,并对切面进行摄影,随后在照片的基础上勾勒出各个脑组织的轮廓线,并用色块进行填充得到。该图谱建立在大脑前联合(AnteriorCommissure,AC)、后联合(PosteriorCommissure,PC)和脑中矢状面为基准的坐标系中。目前有数字化的Talairach图谱提供给研究人员使用。图谱为单模图谱,仅有36个矢状层,因此使用时需要进行插值。该图谱具有丰富的解剖位置和形状信息,经常应用于立体神经外科手术以及计算机辅助神经放射学中。15 第二章基于ICBM-152图谱的临床数据空间标准化图2-3Talairach图谱2.1.4MNI脑模板MNI脑模板是由加拿大蒙特利尔神经研究所(MontrealNeurologicalInstitute,MNI)制作的一系列脑图谱。在临床实践中,功能影像的分析常叠加于MNI脑模板中。MNI305是第一个MNI模板,首先使用241组正常的MR数据,手动定义了各种标志点,每组数据都被缩放至Talairach图谱上;随后再次使用了305组正常的MR数据,将其缩放至上述241组缩放后的数据的均值上,再对所得到的305组缩放后的数据取平均值,最后得到MNI305图谱。目前最常用的标准是ICBM-152(ICBM-InternationalConsortiumforBrainMapping)图谱,它是由152名23岁左右的健康成年人的MR头部扫描所得序列配准至MNI305上所得的平均值。MNI采用了同样的方法,推出了比ICBM-152图谱分辨率更高的ICBM-452图谱,但是该图谱用途尚不广泛。图2-4ICBM-152图谱在这一系列MNI脑模板中,ICBM-152图谱当前使用最为广泛,在许多研究者的努力下,基于ICBM-152图谱空间已经延伸出了众多以该图谱为基准的脑组织图谱,如穹窿模板、皮质层及皮质下结构图谱、基于DTI的白质图谱等,其中也包括了丘脑下核团的图谱。因此本文选用ICBM-152图谱作为后续算法的标准图谱。16 东南大学硕士学位论文2.2ICBM-152图谱空间图谱空间指在脑部选取一些不变的部位作为界标,脑中其余组织可以根据这些界标来定位,从而对脑图像进行规格化的一种方法。在脑部规格化之后,可以实现不同个体的脑图像比较。由于ICBM-152图谱是在Talairach图谱的基础上构建的,因此二者图谱空间的定义十分相似。ICBM-152图谱空间以大脑前连合(AC)、后连合(PC)以及脑中膈面作为基准(如图2-5所示),将前连合与后连合之间的连线(AC-PC线)作为水平轴,将脑中膈面作为垂直面,该规格化方法在神经外科立体定向中被广泛使用。图2-5前连合、后连合以及ICBM坐标系该图谱空间建立的重要标识是通过大脑中央回的脑中膈面。该平面具有最为显著的解剖学特征,且构成图谱空间的AC、PC两个基准点都处于这一平面上,因此对于临床数据空间标准化的第一步是在体数据中识别其脑中膈面。2.3基于K-L距离的脑中膈面识别算法脑中膈面是分割左右大脑的平面,确定脑中隔面是将临床数据坐标空间统一至图谱坐标空间的重要步骤。本课题使用的是Volkau等人所提出的基于Kullback-Leibler度量(K-L距离)的脑中膈面提取法,该方法是一个基于统计学的方法,使用K-L距离表征两个分布之间的差异,不需要对数据做任何如重新格式化、颅骨剥离等预处理,在不同模态的数据(T1WI,T2WI,FLAIR,CT)中都有较好的鲁棒性[38]。17 第二章基于ICBM-152图谱的临床数据空间标准化2.3.1K-L距离的定义假设一个离散随机变量X的概率分布为p={𝑝𝑝𝑖𝑖},𝑝𝑝𝑖𝑖表示系统第i个状态发生的可能性,则表达式log(1/𝑝𝑝𝑖𝑖)的值被称为意外性(unexpectednessorsurprise)。如果𝑝𝑝𝑖𝑖=1,表示事件i必定会发生,那么事件的意外性为0。如果事件几乎不可能发生,即𝑝𝑝𝑖𝑖≈0时,则意味着事件i发生的意外性无穷大。假设两个离散分布p={𝑝𝑝𝑖𝑖}和q={𝑞𝑞𝑖𝑖},𝑝𝑝𝑖𝑖和𝑞𝑞𝑖𝑖表示他们在系统第i个状态时发生的可能性。则变量p关于变量q的意外性可以被定义为log(1𝑞𝑞𝑖𝑖)−log(1𝑝𝑝𝑖𝑖)。通过对𝑝𝑝𝑖𝑖的意外性做整体平均可以得到:I(𝑝𝑝𝑞𝑞)=∑𝑖𝑖𝑝𝑝𝑖𝑖log𝑝𝑝𝑖𝑖−∑𝑖𝑖𝑝𝑝𝑖𝑖log𝑞𝑞𝑖𝑖=∑𝑖𝑖𝑝𝑝𝑖𝑖log(𝑝𝑝𝑖𝑖/𝑞𝑞𝑖𝑖)(2.1)表达式I(𝑝𝑝𝑞𝑞)被称为两个概率分布p与q的Kullback-Leibler度量(K-L距离),表示了两种分布之间的差异程度。我们可以通过K-L距离定义两个矢状面的强度概率变化。2.3.2问题描述假设p={𝑝𝑝𝑖𝑖}和q={𝑞𝑞𝑖𝑖}代表两个矢状面的灰度分布状态。𝑝𝑝𝑖𝑖和𝑞𝑞𝑖𝑖的值代表了在图像中灰度值为i的像素出现的可能性,该值可以通过灰度直方图来计算。对公式(2.1)做如下约束:(1)0log(0𝑝𝑝𝑖𝑖)=0;(2)分母为0的项不参与统计。脑中膈面也被称为中矢状面(midsagittalplane,MSP),被定义为穿过半球间裂隙(interhemisphericfissure,IF)的平面,其中包含了前连合、后连合、脑室、丘脑和胼胝体等脑组织,其表现出近似的双侧对称性。在正常数据集中,脑中膈面可以被认为是平面,其在所有矢状切片中具有最大量的脑脊液(cerebrospinalfluid,CSF)。在理想的情况下,如果半球间裂隙不是弯曲的,且体素的矢状方向宽度小于半球间裂隙的宽度,则脑中膈面通常包含以下大脑结构:前连合和后连合第三以及第四脑室,脑沟以及半脑间裂隙的脑脊液胼胝体、脑干、松果体、小脑和中间层结构的横截面非皮层区域的横断面18 东南大学硕士学位论文图2-6脑中膈面空间及解剖示意图2-7ICBM-152图谱脑中膈面其两侧距其2cm、4cm处矢状面图像以ICBM-152图谱为例,其脑中膈面及其两侧距其2cm、4cm处的矢状面图像如图2-7所示,从中可以看到脑中膈面的平均灰度相比其他矢状面更小。解剖学上的脑中膈面与其他矢状面的区分在矢状面的MR图像中很容易使用灰度直方图来表达(如图2-8所示)。典型的脑中膈面图像的灰度直方图的峰值出现在脑脊液所对应的灰度区。与之相反,非脑中膈面图像的灰度直方图的峰值则出现在白质所对应的灰度区。图2-8ICBM-152图谱脑中膈面及距其2cm处灰度直方图当我们将脑中膈面的像素强度作为随机事件考虑时,脑脊液所对应的像素出现概率较高,白质所对应的像素出现的概率较低,而非脑中隔面的图像则反之。两张图像之间19 第二章基于ICBM-152图谱的临床数据空间标准化的K-L距离表示了对应的像素灰度出现概率的变化情况。若选择非脑中膈面作为基准图像与目标图像计算K-L距离时,若该图像为非脑中隔面图像,则他们像素出现的概率相似,二者K-L距离较小;若该图像为脑中膈图像时,他们像素出现的概率相反,二者K-L距离较大。因此K-L距离可以作为判断脑中膈面的度量。2.3.3算法设计首先对算法使用的图像三维坐标空间进行定义,设(𝑥𝑥𝑚𝑚,𝑦𝑦𝑚𝑚,𝑧𝑧𝑚𝑚)代表了整数体素坐标。体数据三维空间坐标示意如图2-9所示,其中x代表了体素坐标的左右方向,y代表了体素坐标的前后方向,z代表了体素坐标的上下方向,坐标原点为体数据的几何中心点。图2-9体数据三维坐标空间示意该算法分为两个阶段。在第一个阶段,通过分析给出的大脑体数据的矢状位图像,获得一个粗略的脑中膈面近似值(cMSP)。在第二个阶段,算法会将cMSP值进行优化,得到脑中膈面的最佳位置。在第一阶段,算法主要进行以下工作:1)读取大脑的数据,以及体素的大小等信息;2)若原始的数据方向是横断位或冠状位方向,需要将其插值生成矢状位方向数据;3)在矢状位方向中心切片周围±2cm范围内定义一个感兴趣体积(VOI);4)选择VOI中最边缘的图像作为初始的基准图像;5)计算VOI中所有其他矢状位图像与基准图像的K-L距离;20 东南大学硕士学位论文6)将其中K-L距离最大的图像作为cMSP。假设当K-L距离达到最大时所对应的图像的x方向坐标为s,则初始的cMSP的四个顶点坐标分别为(s,𝑦𝑦1,𝑧𝑧1),(s,𝑦𝑦2,𝑧𝑧2),(s,𝑦𝑦3,𝑧𝑧3),(s,𝑦𝑦4,𝑧𝑧4)。在第二阶段中,使用第一阶段获得的cMSP坐标,经过精确查找过程,通过改变cMSP的顶点,找到最大K-L距离值,所得的结果即为最优的脑中膈面。算法的细节如下:1)选择cMSP中任意三个顶点构成一个平面,使用(s,𝑦𝑦𝑖𝑖,𝑧𝑧𝑖𝑖)表示,其中i=1,2,3;2)设置矢状位方向上距离cMSP平面2cm处的平面作为新的基准图像,后续过程中基准图像不再变化;3)设L为当前VOL中矢状位面图像的总数,以𝛥𝛥0=min(⌈𝐿𝐿8⌉,8)为初始步长沿矢状位方向搜索,搜索只在VOL内进行;4)设每一步平面移动的顶点为(a,𝑦𝑦1,𝑧𝑧1),(b,𝑦𝑦2,𝑧𝑧2),(c,𝑦𝑦3,𝑧𝑧3),其中a,b,c∈[s±𝛥𝛥0](选择三个点在VOL边界上构成了一个新的平面,其中y,z为常数,只有x坐标被改变,这样就定义了沿矢状位方向的VOL中待测平面的新位置),根据上述点的定义计算其与基准图像间的K-L距离;5)当前K-L距离最大值的一组顶点所对应的平面为新的感兴趣区域中的初始搜索平面,新的步长𝛥𝛥=⌈𝛥𝛥2⌉,新的感兴趣区域的的大小在X方向上压缩到4𝛥𝛥,计算新的VOL中的每个平面与基准面的K-L距离;6)当步长大于等于1时重复上一步;7)算法停止时与基准面的K-L距离最大的一组顶点构成的面就是所求的脑中膈面。2.4临床数据校准至图谱空间根据以上算法就可以获得最终构成脑中膈面的三个顶点坐标a=(𝑥𝑥1,𝑦𝑦1,𝑧𝑧1),b=(𝑥𝑥2,𝑦𝑦2,𝑧𝑧2),c=(𝑥𝑥3,𝑦𝑦3,𝑧𝑧3),计算出该平面所对应的法向量N=(𝑛𝑛𝑥𝑥,𝑛𝑛𝑦𝑦,𝑛𝑛𝑧𝑧)以及平面中心o=(𝑜𝑜𝑥𝑥,𝑜𝑜𝑦𝑦,𝑜𝑜𝑧𝑧),将整个体数据绕着平面中心旋转,旋转度数为平面法向量与单位向量I=(0,1,0)的夹角θ,并使用旋转公式将数据旋转至脑中膈面垂直于y轴的方向。其示意图如图2-10所示,其中红色平面为初始MSP平面,蓝色平面为目标平面。21 第二章基于ICBM-152图谱的临床数据空间标准化图2-10图像初次旋转示意随后在旋转后的脑中膈面上标出AC、PC点,计算出AC-PC连线与x轴的夹角φ,再次对数据进行旋转,最终所得的体数据就是校准至ICBM-152图谱标准空间的数据了。图2-11图像AC-PC旋转示意对于体数据的旋转,涉及到空间向量旋转的概念,其方法如下:𝑇𝑇考虑对一个向量做旋转,其中𝑣𝑣⃗是原向量,三维单位向量𝑘𝑘�⃗=�𝑘𝑘𝑥𝑥𝑘𝑘𝑦𝑦𝑘𝑘𝑧𝑧�是旋转轴,θ是旋转角,𝑣𝑣���𝑟𝑟𝑟𝑟𝑟𝑟����⃗是旋转后的向量。首先通过点积得到𝑣𝑣⃗在𝑘𝑘�⃗的方向上的平行分量𝑣𝑣���‖⃗。�𝑣𝑣��‖⃗=�𝑣𝑣⃗·𝑘𝑘�⃗�𝑘𝑘�⃗(2.2)再通过叉乘得到与𝑘𝑘�⃗正交的两个向量𝑣𝑣����⊥⃗和𝑤𝑤��⃗。𝑣𝑣����⊥⃗=𝑣𝑣⃗−𝑣𝑣���‖⃗=𝑣𝑣⃗−�𝑣𝑣⃗·𝑘𝑘�⃗�𝑘𝑘�⃗=−𝑘𝑘�⃗×(𝑘𝑘�⃗×𝑣𝑣⃗)(2.3)22 东南大学硕士学位论文𝑤𝑤��⃗=𝑘𝑘�⃗×𝑣𝑣⃗(2.4)这三个向量互相正交,不难得出:𝑣𝑣���𝑟𝑟𝑟𝑟𝑟𝑟����⃗=𝑣𝑣���‖⃗+cos𝜃𝜃𝑣𝑣����⊥⃗+𝑠𝑠𝑠𝑠𝑛𝑛𝜃𝜃𝑤𝑤��⃗(2.5)𝑇𝑇设K为𝑘𝑘�⃗=�𝑘𝑘𝑥𝑥𝑘𝑘𝑦𝑦𝑘𝑘𝑧𝑧�的叉积矩阵,显然K是一个反对称矩阵。0−𝑘𝑘𝑧𝑧𝑘𝑘𝑦𝑦𝐾𝐾=�𝑘𝑘𝑧𝑧0−𝑘𝑘𝑥𝑥�(2.6)−𝑘𝑘𝑦𝑦𝑘𝑘𝑥𝑥0它满足如下性质:𝑘𝑘�⃗×𝑣𝑣⃗=𝐾𝐾𝑣𝑣⃗(2.7)利用该性质,需要将𝑣𝑣���𝑟𝑟𝑟𝑟𝑟𝑟����⃗代换为𝑣𝑣⃗与𝑘𝑘�⃗的叉积可得:𝑣𝑣���‖⃗=𝑣𝑣⃗+𝑘𝑘�⃗×(𝑘𝑘�⃗×𝑣𝑣⃗)(2.8)随后可以推知:𝑣𝑣���𝑟𝑟𝑟𝑟𝑟𝑟����⃗=𝑣𝑣⃗+𝑘𝑘�⃗×�𝑘𝑘�⃗×𝑣𝑣⃗�−cos𝜃𝜃𝑘𝑘�⃗×�𝑘𝑘�⃗×𝑣𝑣⃗�+𝑠𝑠𝑠𝑠𝑛𝑛𝜃𝜃𝑘𝑘�⃗×𝑣𝑣⃗(2.9)根据叉乘矩阵的性质:𝑣𝑣�������⃗=𝑣𝑣⃗+(1−cos𝜃𝜃)𝐾𝐾2+𝑠𝑠𝑠𝑠𝑛𝑛𝜃𝜃K𝑣𝑣⃗(2.10)𝑟𝑟𝑟𝑟𝑟𝑟𝑣𝑣�������⃗=(𝐼𝐼⃗+(1−cos𝜃𝜃)𝐾𝐾2+𝑠𝑠𝑠𝑠𝑛𝑛𝜃𝜃K)𝑣𝑣⃗(2.11)𝑟𝑟𝑟𝑟𝑟𝑟最后将𝑣𝑣⃗、𝑣𝑣�������⃗替换为B、C,就得到了罗德里格斯公式[39]的标准形式。𝑟𝑟𝑟𝑟𝑟𝑟B=�𝐼𝐼⃗+(1−cos𝜃𝜃)𝐾𝐾2+𝑠𝑠𝑠𝑠𝑛𝑛𝜃𝜃K�C⇔R=𝐼𝐼⃗+(1−cos𝜃𝜃)𝐾𝐾2+𝑠𝑠𝑠𝑠𝑛𝑛𝜃𝜃K(2.12)其中R为最终的旋转矩阵。假设需要围绕的单位向量为�𝑛𝑛𝑥𝑥,𝑛𝑛𝑦𝑦,𝑛𝑛𝑧𝑧�,需要旋转的角度为θ,则根据公式,所求得的旋转矩阵为:𝑛𝑛𝑥𝑥2(1−cosθ)+cosθ𝑛𝑛𝑥𝑥𝑛𝑛𝑦𝑦(1−cosθ)+𝑛𝑛𝑧𝑧sinθ𝑛𝑛𝑥𝑥𝑛𝑛𝑧𝑧(1−cosθ)−𝑛𝑛𝑦𝑦sinθ�𝑛𝑛𝑥𝑥𝑛𝑛𝑦𝑦(1−cosθ)−𝑛𝑛𝑧𝑧sinθ𝑛𝑛𝑦𝑦2(1−cosθ)+cosθ𝑛𝑛𝑦𝑦𝑛𝑛𝑧𝑧(1−cosθ)+𝑛𝑛𝑥𝑥sinθ�𝑛𝑛𝑥𝑥𝑛𝑛𝑧𝑧(1−cosθ)+𝑛𝑛𝑦𝑦sinθ𝑛𝑛𝑦𝑦𝑛𝑛𝑧𝑧(1−cosθ)−𝑛𝑛𝑥𝑥sinθ𝑛𝑛𝑧𝑧2(1−cosθ)+cosθ(2.13)2.5实验结果与分析本文在VC++编程平台上对上述算法进行了实现,分别使用ICBM-152图谱作为模拟数据以及临床图像作为实际数据对算法进行了验证。23 第二章基于ICBM-152图谱的临床数据空间标准化2.5.1模拟数据结果首先,本文通过对ICBM-152图谱进行绕Z轴旋转,模拟患者在接受MR扫描时的不同角度的头部摆放效果,其初始图像如图2-12所示:图2-12绕z轴旋转各个角度(5°、10°、15°和20°)的模拟数据随后使用本章所介绍的脑中膈面矫正算法对各组模拟数据进行矫正,结果如图2-13:图2-13矫正算法对各个角度(5°、10°、15°和20°)模拟数据的恢复效果24 东南大学硕士学位论文由矫正后的图像我们可以看出,算法将模拟数据矫正至了中膈面垂直的状态(由于在生成模拟数据的过程中15°以及20°的数据发生数据越界,因此这两幅数据在矫正后产生了右下角的部分图像缺失)。将每幅图像旋转度数的计算结果与原始的图像旋转度数进行对比,其结果如表2-1所示:表2-1原始旋转度数与算法计算结果数据编号原始旋转度数(°)算法结果(°)误差(°)154.970.032109.880.1231514.630.3742019.810.19由上述结果可知,算法的平均误差为0.18°,因此该算法具有较好的鲁棒性。随后本文对算法的搜索规则进行了验证,通过记录每次搜索时的当前最大K-L距离进行实现,对上述四组模拟数据的搜索次数与最大K-L距离绘制了曲线,其结果见图2-14:图2-14算法搜索次数与最大K-L距离的关系曲线由上述曲线我们可以发现,随着搜索次数的增大,当前对应的基准面与搜索平面的最大K-L距离也在增大,当算法停止搜索时,最大K-L距离所对应的平面即为中膈面,这与我们的预期相符,因此本章所使用的算法是有效的。25 第二章基于ICBM-152图谱的临床数据空间标准化2.5.2实际数据结果在对本章所使用的算法使用模拟数据进行分析验证之后,本文使用了该算法对实际临床数据进行了处理,其中部分处理结果如下:图2-15算法对测试数据1的矫正效果原始数据(左)算法结果(右)图2-16算法对测试数据2的矫正效果原始数据(左)算法结果(右)观察图2-16可以看到,矫正前图像临床数据中的脑中膈面存在着一定角度的偏转,在经过算法处理后脑中膈面基本处于垂直状态,由于临床数据为三维数据,因此我们可以之间观察过体数据坐标原点的矢状面来评估算法的处理结果。由于患者临床数据存在的偏转角度通常较小,因此通常不会出现数据越界的情况,若数据偏转角度较大,可以通过对数据边缘进行补零,再使用算法对数据进行处理来避免该问题。26 东南大学硕士学位论文图2-17测试数据1过体数据坐标原点矢状面(左)算法处理结果(右)图2-18测试数据2过体数据坐标原点矢状面(左)算法处理结果(右)图2-19测试数据1(a)算法所得中矢状面前一帧图像(b)算法所得中矢状面图像(c)算法所得中矢状面后一帧图像27 第二章基于ICBM-152图谱的临床数据空间标准化图2-20测试数据2(a)算法所得中矢状面前一帧图像(b)算法所得中矢状面图像(c)算法所得中矢状面后一帧图像观察以上两组数据(图2-17及图2-18)中的红圈标识部分,我们仍然可以在初始的过体数据坐标原点的矢状面中清楚的看到脑沟等结构,说明该平面并不是解剖学上的脑中膈面,而在算法处理后,图像中的多数脑沟边缘变得模糊甚至不可见,说明该平面中对应脑脊液的部分较多,这与解剖学上的脑中膈面定义相符,从图2-19以及图2-20观察,算法所得中矢状面与其前后一帧的矢状面对比同样符合上述性质。分别比较两组数据的灰度直方图(图2-21)可以看到,算法执行后的中膈面中灰度值为200左右的像素数相比于执行前增加了,而灰度值为300左右的像素数将较于算法执行后减少了,这与模型建立时的预测相同。图2-21算法执行前后数据中膈面灰度直方图比较由上述结果可知,算法对实际数据同样有较好的处理效果,在中膈面被校准后,在算法所得的中膈面上手动选取AC-PC点,计算AC-PC连线y轴夹角度数,并根据所得度数旋转体数据,使图像的数据空间与ICBM-152图谱的数据空间达到最终的匹配,其结果如图2-22所示:28 东南大学硕士学位论文图2-22AC-PC旋转前(左)AC-PC旋转后(中)脑图谱中膈面(右)由上图,我们可以观察到,旋转后的图像与脑图谱在视觉上已经有了一定的相似之处。算法处理后,临床数据中脑中膈面上的AC及PC点清晰可见,且二者连线已处于水平位置。因此本章所介绍的图像预处理算法有效地解决了临床图像与脑图谱的坐标空间标准化问题。2.6本章小结本章对标准图谱空间进行了介绍。选择ICBM-152图谱,通过使用基于K-L距离作为测度来实现对原始数据的脑中膈面提取,在脑中膈面上手动定位AC-PC点,并对数据进行旋转,使原始数据坐标空间统一到标准脑图谱的坐标空间,实现了对原始临床数据的坐标空间标准化。实验结果表明,本文采用的方法对于解决临床数据与图谱数据坐标空间不一致的问题是有效的,而且效果显著。29 第三章基于脑皮质标识点的分段仿射配准第三章基于脑皮质标识点的分段仿射配准上一章统一了临床数据与标准脑图谱坐标空间,但临床脑图像大小形状千差万别,存在显著的个体差异,而ICBM-152图谱是用于展示人类脑组织结构的模拟数据,与患者的脑部图像在分辨率、大小、形状以及内部纹理等方面都有着较大的差异,因此一步完成二者精确配准较为困难。本章将Talairach脑图谱所包含的脑皮层标识点概念套用至ICBM-152图谱中,并使用分段仿射变换完成了图谱与人脑的粗配准问题[40],缩小二者之间的大小差异,为下一步的精确配准过程提供一个良好的初始状态。3.1大脑分区以及脑皮层标识点3.1.1大脑分区Talairach图谱采用了三维栅格系统,根据该系统依据AP、MSP、VAC以及VPC这四个平面将脑部划分成12个区间,其示意图如图3-1所示。其中AP平面表示过AC、PC的横断面,VAC平面为过AC的冠状面,VPC平面为过PC点的冠状面。图3-1四个基准面在脑部体数据中的示意图MR全脑体数据的12个区间如图3-2所示:30 东南大学硕士学位论文图3-2体数据12分区示意图3.1.2脑皮质层标识点脑皮层标识点均位于大脑皮质层的外边缘,一共包括6个,分别记为A、P、R、L、S以及I,他们的原始定义如下:A:额前皮质的最前点;P:枕侧皮质的最后点;R:右脑半球顶颞皮质的最边缘点;L:左脑半球顶颞皮质的最边缘点;S:颅顶皮质的最高点;I:颞颥皮质的最低点。文献[41]研究了Talairach图谱中的标识点,并指出了原始的脑皮质标记点在Talairach图谱中的缺陷:首先这些标识点存在冗余,在建立Talairach图谱空间时,部分标识点在脑图谱中的位置和其定义产生了矛盾,是不必要的;其次,该定义下皮层标识点的定位方法十分困难,需要十分强的医学背景知识。因此Nowinski对原始的Talairach皮层标识点的定义进行了改进,重新定义了一组更符合实际算法运用的标识点——Talairach-Nowinski标识点,其定义如下:A:AP平面、MSP平面以及通过额前皮质最前端点的冠状面的这三个平面的交界点;P:AP平面、MSP平面以及通过枕侧皮质最后端点的冠状面的这三个平面的交界点;R:AP平面、过PC的冠状面以及通过右脑半球顶颞皮质层最边缘点的矢状面的这三个平面的交界点;L:AP平面、过PC的冠状面以及通过左脑半球顶颞皮质层最边缘点的矢状面的这三个平面的交界点;31 第三章基于脑皮质标识点的分段仿射配准S:通过PC的冠状面、MSP平面以及通过颅顶皮质最高点的横断面的这三个平面的交界点;I:通过AC的冠状面、MSP平面以及通过颞颥皮质最低点的横断面的这三个平面的交界点。改进后的标识点定义简洁明了,应用意义更大,在三个平面上的位置如图3-3所示:图3-36个脑皮质标识点在三个基准面上的示意图ICBM-152图谱与Talairach图谱在空间上有着一定的相似性,因此本文将Talairach图谱中的大脑分区以及Talairach-Nowinski标识点的概念引用到ICBM-152图谱中,用于解决临床图像与脑图谱在全局上的粗配准问题。3.2脑皮层标识点的提取本文使用的粗配准方式为分段线性配准,因此需要计算图谱中各个脑分区与临床数据中对应脑分区的缩放比例,这些比例系数是通过脑皮层标识点所提供的,因此对于患者脑皮质层标识点的提取就显得尤为重要。AP平面上的标识点数量最多,因此对于AP平面脑组织的准确分割,是获取脑皮质标识点的关键步骤。由于原始脑数据中存在较多的边缘,且图像的质量层次不齐,因此本文采用聚类算法处理图像[42],获取图像中代表背景、脑脊液、白质、灰质以及颅骨的像素强度分布范围;随后将像素强度较小的脑脊液以及背景的像素强度范围投射至低值区间,将像素强度较大的灰质、白质以及颅骨投射至高值区间,从而起到数据二值化的效果;随后使用取得图像中的最大连通域,对其进行形态学的操作,得到脑组织的分割结果,并根据定义提取脑皮层标识点。32 东南大学硕士学位论文3.2.1FCM聚类算法FCM算法是模糊C均值(FuzzyC-means)算法的简称,是一种基于目标函数的模糊聚类算法。算法通过对目标函数进行迭代优化,获取数据集的模糊分类,在图像的特征分析、模式识别以及智能信息处理等方面有较多应用。FCM算法的主要思想是使被划分至同一类的对象之间的相似度达到最大,而使不同类的对象之间的相似度达到最小。FCM算法建立在模糊集合的概念之上。算法中的隶属度是模糊集合的重要概念,隶属度函数用于描述一个对象x隶属于某一集合M的程度,通常记作𝜇𝜇𝑀𝑀(𝑥𝑥),函数的自变量范围为集合M所在空间的所有点,值域为[0,1]。𝜇𝜇𝑀𝑀(𝑥𝑥)=1表示x完全隶属于集合M,相当于传统集合中的x∈M。一个模糊集合M可以用一个位于空间X={x}中的隶属度函数来定义。对于有限个对象𝑥𝑥1,𝑥𝑥2,…,𝑥𝑥𝑛𝑛,模糊集合M表示为:M={(𝜇𝜇𝑀𝑀(𝑥𝑥𝑖𝑖),𝑥𝑥𝑖𝑖)|𝑥𝑥𝑖𝑖∈M}(3.1)在模糊集合的概念之上某个元素隶属于模糊集合不再是硬性的。在聚类问题中,若将聚类所生成的簇看作模糊集合,那么每个样本点对于簇的隶属度便是[0,1]区间中的值,该值表示一个样本相似于不同结果的程度。假设数据集为X,并把这些数据划分成c类,那么对应的有c个类的中心集合C,每一个样本j属于某一类i的隶属度为𝑢𝑢𝑖𝑖𝑖𝑖,那么定义一个FCM目标函数及其约束条件如下:𝐶𝐶𝑛𝑛𝑚𝑚2J=∑𝑖𝑖=1∑𝑖𝑖=1𝑢𝑢𝑖𝑖𝑖𝑖�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�(3.2)∑𝑐𝑐𝑢𝑢=1,𝑗𝑗=1,2,…,𝑛𝑛(3.3)𝑖𝑖=1𝑖𝑖𝑖𝑖由目标函数(3.2)可知,该函数是由相应样本的隶属度与该样本到各个类中心的距离相乘组成的,指数m是隶属度因子,属于样本的轻缓程度。式(3.3)为约束条件,表示一个样本与其属于所有类的样本的隶属度之和必须为1。根据约束条件,通过迭代运算求式(3.2)的极值就是聚类结果。使用拉格朗日乘数法将约束条件代入目标函数中,在前方加上系数,并把约束条件中的j项展开,式(3.2)就变换成如下形式:𝐶𝐶𝑛𝑛𝑚𝑚2ccJ=∑𝑖𝑖=1∑𝑖𝑖=1𝑢𝑢𝑖𝑖𝑖𝑖�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�+𝜆𝜆1(∑i=1𝑢𝑢𝑖𝑖1−1)+⋯+𝜆𝜆𝑛𝑛(∑i=1𝑢𝑢𝑖𝑖𝑛𝑛−1)(3.4)要求式(3.4)的极值,需要对其中的变量𝑢𝑢𝑖𝑖𝑖𝑖和𝑐𝑐𝑖𝑖进行求导,首先对𝑢𝑢𝑖𝑖𝑖𝑖进行求导,分析式(3.4)可以将前部分求和公式展开成如下矩阵:33 第三章基于脑皮质标识点的分段仿射配准𝑢𝑢11𝑚𝑚‖𝑥𝑥1−𝑐𝑐1‖2⋯𝑢𝑢1𝑚𝑚𝑛𝑛‖𝑥𝑥𝑛𝑛−𝑐𝑐1‖2�⋮⋱⋮�(3.5)𝑢𝑢𝑐𝑐𝑚𝑚1‖𝑥𝑥1−𝑐𝑐𝑐𝑐‖2⋯𝑢𝑢𝑐𝑐𝑚𝑚𝑛𝑛‖𝑥𝑥𝑛𝑛−𝑐𝑐𝑐𝑐‖2𝑚𝑚2该矩阵需要对𝑢𝑢𝑖𝑖𝑖𝑖进行求导,即只需要对有𝑢𝑢𝑖𝑖𝑖𝑖对应的𝑢𝑢𝑖𝑖𝑖𝑖�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�项进行保留,其他𝑚𝑚2所有项中因为不包含𝑢𝑢𝑖𝑖𝑖𝑖因此求导结果都为0。𝑢𝑢𝑖𝑖𝑖𝑖�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�对𝑢𝑢𝑖𝑖𝑖𝑖的求导结果为m�𝑥𝑥𝑖𝑖−2𝑚𝑚−1𝑐𝑐𝑖𝑖�𝑢𝑢𝑖𝑖𝑖𝑖。对于后半部分多项式的求导,同样将求和展开,再除去与𝑢𝑢𝑖𝑖𝑖𝑖不相关的部分,只剩下了𝜆𝜆𝑖𝑖(𝑢𝑢𝑖𝑖𝑖𝑖−1),对于𝑢𝑢𝑖𝑖𝑖𝑖求导的结果为𝜆𝜆𝑖𝑖。最终J对𝑢𝑢𝑖𝑖𝑖𝑖的求导结果等于0时,就是J的极值,即𝜕𝜕𝜕𝜕2𝑚𝑚−1=m�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�𝑢𝑢𝑖𝑖𝑖𝑖+𝜆𝜆𝑖𝑖=0(3.6)𝜕𝜕𝑢𝑢𝑖𝑖𝑖𝑖将该式化简,并将𝑢𝑢𝑖𝑖𝑖𝑖解出可得:1−𝜆𝜆𝑖𝑖𝑚𝑚−11𝑢𝑢𝑖𝑖𝑖𝑖=��(2)(3.7)𝑚𝑚()�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�𝑚𝑚−1解出𝑢𝑢𝑖𝑖𝑖𝑖需要将𝜆𝜆𝑖𝑖进行替换,将式(3.3)的约束条件再次使用,将计算出的𝑢𝑢𝑖𝑖𝑖𝑖代入式中可得:11𝑐𝑐𝑐𝑐−𝜆𝜆𝑖𝑖𝑚𝑚−11−𝜆𝜆𝑖𝑖𝑚𝑚−1𝑐𝑐11=∑𝑖𝑖=1𝑢𝑢𝑖𝑖𝑖𝑖=∑𝑖𝑖=1�𝑚𝑚�(2)=��∑𝑖𝑖=1(2)(3.8)()𝑚𝑚()�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�𝑚𝑚−1�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�𝑚𝑚−1这样就有:1−𝜆𝜆𝑖𝑖𝑚𝑚−111�𝑚𝑚�=𝑐𝑐1=𝑐𝑐1(3.9)∑𝑖𝑖=1(2)∑𝑘𝑘=1(2)�����𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�𝑚𝑚−1�𝑥𝑥𝑖𝑖−𝑐𝑐𝑘𝑘�𝑚𝑚−1将上式重新代入式(3.7)中可得:⎛⎞111𝑢𝑢𝑖𝑖𝑖𝑖=⎜⎟�2�=2(3.10)⎜⎟����1�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�𝑚𝑚−1�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�𝑚𝑚−1∑𝑐𝑐��∑𝑐𝑐��𝑘𝑘=12𝑘𝑘=1�𝑚𝑚−1��𝑥𝑥𝑖𝑖−𝑐𝑐𝑘𝑘�⎝�𝑥𝑥𝑖𝑖−𝑐𝑐𝑘𝑘�⎠式(3.10)就是最终的𝑢𝑢𝑖𝑖𝑖𝑖的迭代公式。随后需要求J对于𝑐𝑐𝑖𝑖的导数,根据式(3.3)可知只有𝐶𝐶𝑛𝑛𝑚𝑚2∑𝑖𝑖=1∑𝑖𝑖=1𝑢𝑢𝑖𝑖𝑖𝑖�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖�部分含有𝑐𝑐𝑖𝑖,对其二级求和展开,则其对𝑐𝑐𝑖𝑖的导数为:𝜕𝜕𝜕𝜕𝑛𝑛𝑚𝑚𝜕𝜕𝑐𝑐=∑𝑖𝑖=1�−𝑢𝑢𝑖𝑖𝑖𝑖×2×�𝑥𝑥𝑖𝑖−𝑐𝑐𝑖𝑖��=0(3.11)𝑖𝑖将其解出可得:∑𝑛𝑛(𝑥𝑥𝑖𝑖𝑢𝑢𝑖𝑖𝑖𝑖𝑚𝑚)𝑖𝑖=1𝑐𝑐𝑖𝑖=𝑛𝑛𝑚𝑚(3.12)∑𝑢𝑢𝑖𝑖=1𝑖𝑖𝑖𝑖34 东南大学硕士学位论文式(3.12)为类中心的迭代公式。由式(3.10)及式(3.12)可知,𝑢𝑢𝑖𝑖𝑖𝑖与𝑐𝑐𝑖𝑖是相互关联的,彼此包含对方。在算法初始时,只需将初值任意赋予其中一个,使其数值满足条件即可。随后便可以开始迭代,循环计算𝑢𝑢𝑖𝑖𝑖𝑖与𝑐𝑐𝑖𝑖的值,使得目标函数J的值逐渐趋向稳定。当J不再变化时算法便收敛到了一个较为稳定的解。式(3.10)及(3.12)是算法的关键,现假设样本集中有三个类(即c=3),同时假设指数m=3,根据式(3.10)可求得样本j到三个类中心的隶属度为:(𝑥𝑥𝑖𝑖−𝑐𝑐2)(𝑥𝑥𝑖𝑖−𝑐𝑐3)𝑢𝑢1𝑖𝑖=(3.13)�𝑥𝑥𝑖𝑖−𝑐𝑐1��𝑥𝑥𝑖𝑖−𝑐𝑐2�+�𝑥𝑥𝑖𝑖−𝑐𝑐1��𝑥𝑥𝑖𝑖−𝑐𝑐3�+(𝑥𝑥𝑖𝑖−𝑐𝑐2)(𝑥𝑥𝑖𝑖−𝑐𝑐3)(𝑥𝑥𝑖𝑖−𝑐𝑐1)(𝑥𝑥𝑖𝑖−𝑐𝑐3)𝑢𝑢2𝑖𝑖=(3.14)�𝑥𝑥𝑖𝑖−𝑐𝑐1��𝑥𝑥𝑖𝑖−𝑐𝑐2�+�𝑥𝑥𝑖𝑖−𝑐𝑐1��𝑥𝑥𝑖𝑖−𝑐𝑐3�+(𝑥𝑥𝑖𝑖−𝑐𝑐2)(𝑥𝑥𝑖𝑖−𝑐𝑐3)(𝑥𝑥𝑖𝑖−𝑐𝑐1)(𝑥𝑥𝑖𝑖−𝑐𝑐2)𝑢𝑢3𝑖𝑖=(3.15)�𝑥𝑥𝑖𝑖−𝑐𝑐1��𝑥𝑥𝑖𝑖−𝑐𝑐2�+�𝑥𝑥𝑖𝑖−𝑐𝑐1��𝑥𝑥𝑖𝑖−𝑐𝑐3�+(𝑥𝑥𝑖𝑖−𝑐𝑐2)(𝑥𝑥𝑖𝑖−𝑐𝑐3)可以看出𝑢𝑢𝑖𝑖𝑖𝑖的分母相同,当𝑢𝑢𝑖𝑖𝑖𝑖越大时,表示到i之外的类的距离的乘积越大,即越靠近i类,且越远离其他类。对于式(3.12),当样本所属的类别i被确定后,该式中的分母求和变为一个常数,则该式可以化简为:∑𝑛𝑛(𝑥𝑥𝑖𝑖𝑢𝑢𝑖𝑖𝑖𝑖𝑚𝑚)𝑢𝑢𝑚𝑚𝑖𝑖=1𝑛𝑛𝑖𝑖𝑖𝑖𝑐𝑐𝑖𝑖=∑𝑛𝑛𝑢𝑢𝑚𝑚=∑𝑖𝑖=1∑𝑛𝑛𝑢𝑢𝑚𝑚𝑥𝑥𝑖𝑖(3.16)𝑖𝑖=1𝑖𝑖𝑖𝑖𝑘𝑘=1𝑖𝑖𝑘𝑘式(3.16)展现了类中心的更新法则,其本质是对样本点的加权平均,在所属的类i确定后,首先需要将所有样本点对该类的隶属度求和,对每个点的隶属度除以该和便是该点所占的比重,将该比重乘以𝑥𝑥𝑖𝑖就是该样本点对于这个类的贡献值。以上分析为式(3.10)及(3.12)对于在迭代时的作用,在实际算法实现时,FCM聚类算法使用如下步骤确定各聚类中心𝑐𝑐𝑖𝑖以及隶属度,隶属度一般使用矩阵存储,记为U:1)确定对样本分类的数目,指数m的值,以及算法的迭代次数(作为终止条件,可以设置迭代次数以外的多种终止条件);2)使用在[0,1]之间的随机数对隶属度矩阵U进行初始化,使其满足式(3.3)中的约束条件;3)根据隶属度矩阵U计算各个类的聚类中心;4)根据当前的聚类中心,更新隶属度矩阵U;5)反复执行步骤(3)(4)直至迭代次数到达上限或者满足其他的终止条件。当算法结束时,每个样本对于每一个类都有一个对应的隶属度u,根据模糊集合中35 第三章基于脑皮质标识点的分段仿射配准最大隶属原则可以确定每个样本所属的类别。而算法所得聚类中心,则表示每个类的平均特征,也可以认为是每个类的代表值。3.2.2FCM聚类以及图像二值化对于实际数据进行处理,根据图像中的脑组织以及其灰度特征,将确定5个FCM聚类中心,分别表示图像中的背景、脑脊液、白质、灰质以及颅骨,随后对图像进行模糊聚类分析,当FCM聚类算法结束时,图像中的每一个像素都会被分类至预先分类的五个类别中。如果分别以1~5大小的值来区分每个类,就可以直观的观察聚类的结果。如图3-4所示:图3-4(a)原始图像(b)聚类后投射的结果随后,去除对于图像中对皮层标识点的定位有干扰的分类(脑脊液、背景以及头骨),保留灰质和白质,即对图像中的像素点根据其类型进行二值化。图3-5图像二值化结果36 东南大学硕士学位论文由于在颅骨边缘软组织的灰度与脑组织相同,同时图像中也存在着噪声,因此二值化的结果并不完全与我们的要求相符,其中包含了许多噪声以及颅骨中的组织,因此需要对二值图像进行形态学图像处理,获得脑组织的边缘信息。3.2.3形态学图像处理数学形态学[43]在图像处理中引用十分广泛,对于二值化后的临床数据,通过使用一定的形态学运算,可以得到较好的脑组织分割结果。3.2.3.1开闭运算与最大连通域二值化临床图像中,脑组织与颅骨之间较为容易发生粘连,若在此基础上提取图像中的最大连通域,往往会导致脑组织连同部分颅骨一起被提取,如图3-6所示,开运算较好的去除了脑皮质与颅骨的粘连情况,由此提取的脑部二值图像更加准确,去除了其他组织的干扰。图3-6(a)原始二值化效果(b)对a提取最大连通区域(c)对a进行开运算(d)对c提取最大连通区域37 第三章基于脑皮质标识点的分段仿射配准3.2.3.2孔洞填充一个孔洞被定义为由前景像素相连接所包围的一个背景区域[44],孔洞填充是一个基于集合膨胀,求补以及交集的算法。该算法可表示为:𝑋𝑋=(𝑋𝑋⊕B)∩𝐴𝐴𝑐𝑐𝑘𝑘=1,2,3…(3.17)𝑘𝑘𝑘𝑘−1其中B为结构元。若𝑋𝑋𝑘𝑘=𝑋𝑋𝑘𝑘−1,则算法在迭代第k步时结束,随后,集合𝑋𝑋𝑘𝑘包含所有的被填充的孔洞,𝑋𝑋𝑘𝑘和A的并集包含所有填充的孔洞及这些孔洞的边缘。使用孔洞填充去除闭运算无法去除的孔洞,并对图像进行边缘提取。结果见图3-7。图3-7(a)对脑组织进行闭运算以及孔洞填充后的结果(b)对a的结果进行边缘提取获取脑皮质边缘的分割结果后,便可以结合AP平面上AC、PC点的坐标进行对A、P、L、R根据其定义进行提取。图3-8(a)边缘信息提取的4个皮层标记点(b)与原始图像进行对应对于VAC及VPC平面,同样可以通过上述算法流程进行处理,算法的部分中间过程及最终结果如图3-9及图3-10:38 东南大学硕士学位论文图3-9VAC平面经过算法处理的中间过程以及分割结果图3-10VPC平面经过算法处理的中间过程以及分割结果随后基于算法的处理结果进行对剩余的S、I这两个标记点的提取。对于Z方向分39 第三章基于脑皮质标识点的分段仿射配准辨率较低的临床图像,在分割结果不理想的情况下,也可以使用手动标记的方法获得这两个标记点。3.2.4实验结果与分析本文在Matlab平台实现了基于FCM聚类的脑皮层标识点提取算法,并使用医院提供的DBS术前MR图像进行了测试,下面给出部分数据的运算结果。对于AP平面中的脑组织的分割问题,所使用的四组测试数据的分割结果及标识点如下:图3-11测试数据1(左)与测试数据2(右)分割结果及标识点图3-12测试数据3(左)与测试数据4(右)分割结果及标识点从分割结果可以看出,本章中所提出的方法较好的解决了AP平面中的脑组织分割问题,分割所得到的脑皮层边缘与实际边缘较为接近。40 东南大学硕士学位论文对于S、I这两个脑皮层标识点则需要对VAC以及VPC平面进行分割,本文同样列出了4组数据的分割结果以及标识点:图3-13测试数据1VAC平面(左)与VPC平面(右)分割结果及标识点图3-14测试数据2VAC平面(左)与VPC平面(右)分割结果及标识点图3-15测试数据3VAC平面(左)与VPC平面(右)分割结果及标识点41 第三章基于脑皮质标识点的分段仿射配准图3-16测试数据4VAC平面(左)与VPC平面(右)分割结果及标识点随后根据分割结果得到皮层标识点的坐标,与手动标记的皮层坐标进行对比,其结果如下表所示。表3-1测试数据1各个标志点结果及误差标识点手动标注结果(mm)算法结果(mm)误差(mm)A206.05203.122.93P38.0939.06-0.97L157.23155.671.56R13.1814.160.98S204.59203.121.47I102.54101.560.98表3-2测试数据2各个标志点结果及误差标识点手动标注结果(mm)算法结果(mm)误差(mm)A189.94186.043.90P16.1517.58-1.43L153.32152.340.98R11.7212.21-0.49S205.08203.121.96I120.12119.630.49表3-3测试数据3各个标志点结果及误差标识点手动标注结果(mm)算法结果(mm)误差(mm)A196.48193.363.12P43.9545.41-1.46L152.34150.411.93R19.5317.581.95S214.84211.912.93I126.46125.980.4842 东南大学硕士学位论文表3-4测试数据4各个标志点结果及误差标识点手动标注结果(mm)算法结果(mm)误差(mm)A192.37188.963.41P35.6636.13-0.47L149.41148.930.48R17.0918.55-1.46S201.17199.221.95I117.68118.16-0.484.543.53数据12.5数据22数据3误差(mm)1.5数据410.50A误差P误差L误差R误差S误差I误差图3-17各组数据误差统计图从各组数据的误差可以得出A点的误差均值为3.27mm,P点的误差均值为1.08mm,L点的误差均值为1.23mm,R点的误差均值为1.22mm,S点的误差均值为2.08mm,I点的误差均值为0.61mm。相比其余点,A点存在的误差较大,分析其原因为大脑前端的脑组织附近存在的脑脊液较多,且脑组织的灰度值与其他区域的脑组织相比较小,导致FCM聚类算法将一部分的脑组织划分为脑脊液类中。3.3基于标识点的大脑分段仿射变换匹配在实现脑皮质层标记点的分割提取后,便可以计算出各个大脑分块与图谱中对应分块的大小差异,并使用分段仿射变换实现脑图谱与图像的粗配准。3.3.1仿射变换的原理与性质一个在两个仿射空间之间的仿射变换,是在向量上呈线性的坐标点的变换。一个对向量𝑥𝑥⃗平移𝑏𝑏�⃗,并旋转、放大、缩小A的仿射变换可以表达为:43 第三章基于脑皮质标识点的分段仿射配准𝑦𝑦⃗=𝐴𝐴𝑥𝑥⃗+𝑏𝑏�⃗(3.18)在齐次坐标系中,该式等价于:𝑦𝑦⃗𝐴𝐴𝑏𝑏�⃗𝑥𝑥⃗��=����(3.19)10,…,011仿射变换保留空间点之间的共线性,即通过同一直线的点在变换后仍然在同一直线上,而且沿着一线的比例也保持不变。根据仿射不变性可知,对于分块后的每一块图谱,在经过仿射变换后其脑组织的相对位置不会发生变化。3.3.2计算分块比例对于3.2节中提取所得的A、P、L、R、S、I六个脑皮层标识点以及AC、PC两点,可以结合3.1节的大脑分区,得到用于计算各块用于仿射变换缩放比例的七个标识点间的距离:A点至AC点的距离𝐿𝐿𝐴𝐴−𝐴𝐴𝐶𝐶;P点至PC点的距离𝐿𝐿𝑃𝑃−𝑃𝑃𝐶𝐶;L点至PC点的距离𝐿𝐿𝐿𝐿−𝑃𝑃𝐶𝐶;R点至PC点的距离𝐿𝐿𝑅𝑅−𝑃𝑃𝐶𝐶;S点至PC点的距离𝐿𝐿𝑆𝑆−𝑃𝑃𝐶𝐶;I点至AC点的距离𝐿𝐿𝐼𝐼−𝐴𝐴𝐶𝐶;AC点至PC点的距离𝐿𝐿𝐴𝐴𝐶𝐶−𝑃𝑃𝐶𝐶。得到MR图像中的上述参数后,我们可以计算12个分块中图像与图谱的比例系数。表3-5图谱与图像的缩放比例系数参数名图谱中的距离(mm)临床图像中的距离(mm)图谱与图像距离的比值𝐿𝐿𝐴𝐴−𝐴𝐴𝐶𝐶69.0060.551.14𝐿𝐿𝑃𝑃−𝑃𝑃𝐶𝐶74.0067.381.10𝐿𝐿𝐴𝐴𝐶𝐶−𝑃𝑃𝐶𝐶28.0025.921.08𝐿𝐿𝐿𝐿−𝑃𝑃𝐶𝐶72.0068.371.05𝐿𝐿𝑅𝑅−𝑃𝑃𝐶𝐶72.0069.821.03𝐿𝐿𝐼𝐼−𝐴𝐴𝐶𝐶45.0027.831.62𝐿𝐿𝑆𝑆−𝑃𝑃𝐶𝐶64.0057.131.12在x,y,z方向上一共有7个缩放参数,而在12个分块中,对于每两个相邻的分块,他们互相接触的两个维度上的缩放比例相同,因此,对于这12个分块进行缩放后的边界长度是相等的,整个图谱数据在缩放后是连续的。对图像进行三线性插值[45],最44 东南大学硕士学位论文后计算出图谱与图像中AC、PC点的平移量,将图谱与图像中的AC、PC点平移至重合,达到初步的配准效果。3.3.3实验结果与分析在计算出各个脑皮层标志点后,利用本章介绍的分段仿射变换,计算出图谱与临床数据的粗略配准结果。本文在visualstudio2010平台上对该算法进行了实现,以下为部分实验结果。图3-18测试数据1:原始图谱(左)临床图像(中)配准结果(右)图3-19测试数据2:原始图谱(左)临床图像(中)配准结果(右)45 第三章基于脑皮质标识点的分段仿射配准图3-20测试数据3:原始图谱(左)临床图像(中)配准结果(右)图3-21测试数据4:原始图谱(左)临床图像(中)配准结果(右)由各组图像的配准结果可以看出,ICBM-152图谱的脑部较为狭长,经过配准点的图谱与临床数据的脑组织在大小上更加一致,为了便于观察,将配准前后的图谱与图像进行融合,其结果如图3-22所示:图3-22配准前(左)及配准后(右)图谱与图像融合效果上图对于分段仿射配准的效果进行了直观的展示,为了对配准的效果进行量化,本46 东南大学硕士学位论文文计算了配准前后图谱与图像之间的归一化互信息值(其概念将在第四章进行详细介绍),归一化互信息值的区间为[1,2],互信息值越大代表两幅图像间的配准效果越好,其结果如表3-6所示:表3-6配准前后图谱与图像的归一化互信息值数据编号配准前归一化互信息值配准后归一化互信息值差值11.06111.11170.050621.05471.10370.049031.06131.12300.061741.06551.11840.0529根据上表可知,对于本实验所使用的四组数据,配准前的图谱与图像的归一化互信息的均值为1.0607,而在配准后二者的归一化互信息均值为1.1142,通过本章所提出的基于脑皮层标识点的分段仿射配准算法,实现了图谱与图像之间的初步配准,使二者之间的互信息增加了5.04%,取得了较好的结果。3.4本章小结本章将Talairach图谱的大脑分区以及脑皮层坐标点的概念使用在ICBM-152图谱空间中,将空间划分为12个子空间。对临床图像的AP、VAC、VPC平面,使用基于FCM的聚类算法,对图像中的组织根据其灰度不同进行聚类,将聚类结果进行二值化投射,并对二值图像使用形态学图像处理提取脑组织边界,实现了对A、P、L、R、S、I这六个脑皮层标识点的提取。最后根据上述皮层标识点,计算出图谱与临床图像对于12个子空间的缩放比例系数,应用分段仿射配准算法,将图谱数据配准至临床数据中。该过程实现了脑图谱与临床图像的粗配准,为后续的精确配准算法提供一个良好的初值。实验结果显示,本章所提出的算法可以较好的将脑图谱与图像进行配准,且方法简单有效。根据互信息指标显示配准后临床图像与图谱归一化互信息相比配准前上升了5.04%。47 第四章基于B样条函数的非刚性配准第四章基于B样条函数的非刚性配准经过第三章的处理过程,图谱在整体上已经与临床图像达到了一致,但是由于人脑存在个体差异,仅对图谱与空间矫正后的临床图像进行分段仿射配准,二者之间仍然存在着较多局部失配的现象。因此需要对它们进行更加精确的非刚性配准,使图谱与图像之间的差异进一步缩小。本章介绍了基于B样条函数的非刚性配准算法,通过该方法实现了临床数据与粗配准图谱的精确配准,并针对单层配准结果中的失配之处,对算法进行了多层B样条自由形变优化。4.1基于B样条函数的非刚性配准算法非刚性图像配准过程可以定义为一个寻找两个图像在空间以及像素强度方面映射的任务[46],该任务可以分解成以下四个部分:特征空间、相似性度量、搜索空间、搜索策略。由于本文所研究的图谱与临床图像的配准涉及多模态的影像数据,因此对于相似性度量本文选择了广泛用于多模态数据配准的互信息测度。B样条自由形变模型可以利用图像的整体信息,精度较高,较多应用于图谱于图像的非刚性配准中,因此本文使用该模型作为搜索空间。针对B样条自由形变模型的计算占用较多资源的问题,L-BFGS优化算法通常被用作搜索策略。下面对各部分展开介绍。4.1.1互信息测度熵与互信息是信息论中的重要概念。在信息系统中,熵代表了事物的不确定性,是系统包含的信息总量的反应,而互信息可以表示两个信息之间的关系,是用于统计两个随机变量相关性的测度。由于图谱与图像之间的配准属于多模影像配准,而且二者之间的特征点较少,因此本文使用互信息作为相似性测度,该测度是由信息熵[47],联合熵以及互信息这三个概念所引出的。4.1.1.1互信息互信息[48]是反映两个随机变量之间的相关性的变量,互信息能够运用到图像配准中主要是因为,物体虽然在两幅图像上的所呈现的灰度可能不同,但是物体的空间结构是一致的,因此在两幅图像的灰度上有彼此之间的对应关系[49]。48 东南大学硕士学位论文设P(𝑋𝑋)为每个事件𝑋𝑋发生的概率,则随机变量X的信息熵的公式[50]如下:𝑖𝑖𝑖𝑖H(X)=−∑𝑖𝑖𝑃𝑃(𝑋𝑋𝑖𝑖)log𝑝𝑝(𝑋𝑋𝑖𝑖)(4.1)设𝑃𝑃𝑋𝑋𝑋𝑋(x,y)两个随机变量X及Y的联合概率分布函数,则他们的联合熵的表达式为:H(X,Y)=−∑𝑥𝑥,𝑦𝑦𝑃𝑃𝑋𝑋𝑋𝑋(x,y)log𝑃𝑃𝑋𝑋𝑋𝑋(x,y)(4.2)假设存在两幅图像,令参考图像为R,浮动图像为F,则二者的互信息计算公式为:MI(R,F)=H(R)+H(F)−H(R,F)(4.3)式中,MI(R,F)为图像R与F的互信息,其中H(R)与H(F)为两幅图像的信息熵,H(R,F)表示两幅图像的联合熵。当两幅图像实现配准时,两幅图像的相关性达到最大,即H(R,F)的值达到最小,两幅图像的互信息值最大,式(4.3)表明最小联合熵与最大互信息在一定程度上是一致的,通常,配准结果越精确,它们之间的联合熵越小,互信息值越大。但是当配准图像中的物体较小,大部分为背景时,如果出现了图像失配,联合熵也会是一个较小的值,而这种情况下,互信息可以有效避免问题的发生,如果重叠区域中的背景区域占比较大,两个边缘熵也会比较小,在此基础上减去较小的联合熵所得的结果不一定为最大值。互信息的也可以使用如式(4.4)进行表达:MI(R,F)=H(R)−H(R|F)=𝐻𝐻(𝐹𝐹)−𝐻𝐻(𝐹𝐹|𝑅𝑅)(4.4)其中,H(R|F)与𝐻𝐻(𝐹𝐹|𝑅𝑅)表示条件熵,当图像F已知时,互信息就是图像F中包含图像R的信息量。4.1.1.2归一化互信息互信息的计算依赖图像重叠区域的信息,当重叠区域较小时,参与统计的像素对也相应较少,因此容易出现局部极值,导致配准失败。对于此问题,Maes等人提出了熵相关系数的概念[51],其计算公式为:2𝑀𝑀𝐼𝐼(𝑅𝑅,𝐹𝐹)ECC(R,F)=(4.5)𝐻𝐻(𝑅𝑅)+𝐻𝐻(𝐹𝐹)Studholme等人在1999年提出了归一化互信息的概念[52],其计算公式为:𝐻𝐻(𝑅𝑅)+𝐻𝐻(𝐹𝐹)NMI(R,F)=(4.6)𝐻𝐻(𝑅𝑅,𝐹𝐹)由式(4.6)可知,两幅图像的归一化互信息是他们的边缘熵与联合熵的比值。根据二者的计算公式,边缘熵可以通过联合熵计算得出,因此该值收到联合熵的约束。当重叠区域较小,且所需要配准的图像没有很好的匹配时,联合直方图的分布将会十分分散,49 第四章基于B样条函数的非刚性配准因此联合熵将会很大,而边缘熵同样很大,因此保证归一化互信息不会扩大,而当图像很好地配准后,联合熵将会减小,从而导致归一化互信息变大。因此归一化互信息与传统互信息相比可以更加准确地反应两幅图像配准结果是否达到最优,其鲁棒性也比传统互信息更好。4.1.2B样条函数自由形变模型样条配准技术假设在浮动图像以及固定图像上存在一组称为控制点的对应点,利用样条函数对从浮动图像到参考图像对应点的位移进行描述,在控制点之间,样条函数提供一个平滑的形变位移空间,通过样条插值,可以获得所有坐标点的位移。目前基于B样条函数的自由形变模型是进行非刚性配准的有效方法,自由形变的形式可以对物体的形变进行有效的拟合。基于B样条函数的配准方法可以利用图像的整体信息,因此其精度较高。B样条曲线的最初是基于差商进行定义的,但是这种方法较为复杂而且结果不稳定。Carl等人[53]在1970提出了使用通过B样条的递推关系计算B样条函数的方法,该方法十分高效,且在节点发生重合时不需要进行调整,本文采用的是三次均匀样条函数。B样条自由形变的目的是构造一个基本函数f,根据其控制点的位置,可以通过函数求得任意位置的值,通过控制点的位移,经过B样条函数可以拟合出图像中所有点的位移。若图像为二维图像,需要求得图像中所有点的位移在X,Y上的运动分量。而在本课题中,由于需要配准的数据为三维数据,因此需要求所有点在X,Y,Z三个方向上的运动分量。假设(x,y,z)是目标图像上的点,(X,Y,Z)是体数据在三个方向上的尺寸,则图像域可以表示为Ω={(x,y,z)|0≤x≤X,0≤y≤Y,0≤z≤Z},令Φ表示一个由𝑛𝑛𝑥𝑥×𝑛𝑛𝑦𝑦×𝑛𝑛𝑧𝑧个控制点𝜙𝜙𝑖𝑖,𝑖𝑖,𝑘𝑘所组成的网格,若使用某种方法将控制点进行移动,从而使图像发生形变,那么图像上的任意一点可以使用式(4.7)进行定义:f(x,y,z)=∑3∑3∑3𝜃𝜃(𝑢𝑢)𝜃𝜃(𝑣𝑣)𝜃𝜃(𝑤𝑤)𝜙𝜙(4.7)𝑙𝑙=0𝑚𝑚=0𝑛𝑛=0𝑙𝑙𝑚𝑚𝑛𝑛𝑖𝑖+𝑙𝑙,𝑖𝑖+𝑚𝑚,𝑘𝑘+𝑛𝑛其中,i=⌊𝑥𝑥/𝛿𝛿𝑥𝑥⌋−1,j=�𝑦𝑦/𝛿𝛿𝑦𝑦�−1,k=⌊𝑧𝑧/𝛿𝛿𝑧𝑧⌋−1,u=𝑥𝑥𝛿𝛿𝑥𝑥−⌊𝑥𝑥/𝛿𝛿𝑥𝑥⌋,v=𝑦𝑦𝛿𝛿𝑦𝑦−�𝑦𝑦/𝛿𝛿𝑦𝑦�,w=𝑧𝑧𝛿𝛿𝑧𝑧−⌊𝑧𝑧/𝛿𝛿𝑧𝑧⌋,⌊⌋表示向下取整,𝜃𝜃𝑖𝑖表示i阶B样条函数的基函数,其定义为:50 东南大学硕士学位论文𝜃𝜃0=(1−𝑝𝑝3)/6𝜃𝜃3−6𝑝𝑝2+4)/61=(3𝑝𝑝�320≤𝑝𝑝<1(4.8)𝜃𝜃2=(−3𝑝𝑝+3𝑝𝑝+3𝑝𝑝+1)/6𝜃𝜃3/63=𝑝𝑝𝛿𝛿𝑥𝑥,𝛿𝛿𝑦𝑦,𝛿𝛿𝑧𝑧为三个维度上网格的间距,则三个维度上控制点的尺寸𝑛𝑛𝑥𝑥=𝑋𝑋𝛿𝛿𝑥𝑥+3,𝑛𝑛𝑦𝑦=(0)𝑌𝑌𝛿𝛿𝑦𝑦+3,𝑛𝑛𝑧𝑧=𝑍𝑍𝛿𝛿𝑧𝑧+3,设初始的控制点为𝜙𝜙𝑖𝑖,𝑖𝑖,𝑘𝑘,则其初始值为式(4.9):(0)𝜙𝜙𝑖𝑖,𝑖𝑖,𝑘𝑘=(i𝛿𝛿𝑥𝑥,j𝛿𝛿𝑦𝑦,k𝛿𝛿𝑧𝑧)(4.9)其中,−1≤i<𝑛𝑛𝑥𝑥−1,−1≤j<𝑛𝑛𝑦𝑦−1,−1≤k<𝑛𝑛𝑧𝑧−1。控制点𝜙𝜙𝑖𝑖,𝑖𝑖,𝑘𝑘作为B样条自由形变模型的参数,控制着非刚性形变的程度。当网格之间的距离确定时,f的值仅取决于网格的控制点,由于每个控制点只能影响其相邻的区域,因此移动其中的一个控制点只能改变该点所控制区域。改变的范围只与控制点间的距离有关,当间距较大时,该形变模型所产生的效果类似更接近于一个全局的非刚性变换,当间距较小时,形变的效果则用于描述局部的形变状态。根据这个特性,我们可以通过控制网格的大小来改变形变模型的刻画精度。4.1.3L-BFGS优化算法在基于B样条函数的自由形变模型中,每个控制点的位移都是通过一定的优化算法迭代使得互信息度量达到最大值所得到的。令我们在每个维度上的控制点数量分别为𝑛𝑛𝑥𝑥、𝑛𝑛𝑦𝑦、𝑛𝑛𝑧𝑧,所求的控制点在X方向,Y方向以及Z方向上都有分量,因此,优化参数的总数为3×𝑛𝑛𝑥𝑥×𝑛𝑛𝑦𝑦×𝑛𝑛𝑧𝑧。假设𝑛𝑛𝑥𝑥=𝑛𝑛𝑦𝑦=𝑛𝑛𝑧𝑧=10,那么所有优化参数的总数将会达到3000个,面对如此多的优化参数,选择适当的优化算法十分重要。迭代是解决优化问题的最优解的通常做法,牛顿法通过在目标函数f(x)现有极小值的附近对函数进行二阶泰勒展开,从而得到下一个极小值点的估计,其迭代公式如下:𝑥𝑥=𝑥𝑥−(𝐻𝐻)−1𝑔𝑔𝑘𝑘=0,1,…,𝑔𝑔=𝑓𝑓′(𝑥𝑥)(4.10)𝑘𝑘+1𝑘𝑘𝑘𝑘𝑘𝑘𝑘𝑘𝑘𝑘牛顿法的核心计算出Hessian矩阵所提供的曲率信息,该算法收敛速度快,但在运算过程中需要计算Hessian矩阵以及它的逆矩阵,需要较大的运算量。L-BFGS优化算法[54]是在BFGS算法基础上衍生出的算法,对于大量参数的优化问题有着较强的解决能力。该算法不需要对Hessian矩阵进行存储,而是通过计算搜索方向,并保留近m次的迭代信息来降低计算对内存的要求,从而提高算法的速度。L-BFGS算法是一种基于拟牛顿法的有限内存优化算法,它基于目标函数f(x)在点𝑥𝑥𝑘𝑘处的二次模型:51 第四章基于B样条函数的非刚性配准f(𝑥𝑥𝑘𝑘+1)=f(𝑥𝑥𝑘𝑘+𝑎𝑎𝑘𝑘𝑑𝑑𝑘𝑘)=f(𝑥𝑥𝑘𝑘)+𝑎𝑎𝑘𝑘𝑔𝑔𝑘𝑘𝑇𝑇𝑑𝑑𝑘𝑘+𝑜𝑜(𝑎𝑎𝑘𝑘‖𝑑𝑑𝑘𝑘‖)(4.11)其中,𝑥𝑥𝑘𝑘为第k次迭代所得的结果,𝑔𝑔𝑘𝑘为函数f(𝑥𝑥)在𝑥𝑥𝑘𝑘处的一阶导数向量,即𝑔𝑔𝑘𝑘=∇f(𝑥𝑥𝑘𝑘),𝑎𝑎为搜索步长,𝑑𝑑为搜索方向。L-BFGS算法的执行步骤如下[55]:𝑘𝑘𝑘𝑘1)选择初始点𝑥𝑥∈𝑅𝑅𝑛𝑛以及非负整数m,初始正定阵𝐻𝐻∈𝑅𝑅𝑛𝑛×𝑛𝑛,误差ε>0,并令k=0002)计算𝑔𝑔𝑘𝑘=∇f(𝑥𝑥𝑘𝑘),若‖𝑔𝑔𝑘𝑘‖≤𝜀𝜀,则停止迭代3)计算𝑑𝑑𝑘𝑘=−𝐻𝐻𝑘𝑘𝑔𝑔𝑘𝑘4)计算𝑎𝑎𝑘𝑘,以𝑥𝑥𝑘𝑘为起点,沿𝑑𝑑𝑘𝑘进行搜索,若满足f(𝑥𝑥𝑘𝑘+𝑎𝑎𝑘𝑘𝑑𝑑𝑘𝑘)=min𝑓𝑓(𝑥𝑥𝑘𝑘+𝑎𝑎𝑘𝑘𝑑𝑑𝑘𝑘),𝑎𝑎≥0则使𝑥𝑥𝑘𝑘+1=𝑥𝑥𝑘𝑘+𝑎𝑎𝑘𝑘𝑑𝑑𝑘𝑘5)如果k>m,那么只保留前m次的向量对(𝑠𝑠𝑘𝑘−𝑖𝑖+1,𝑦𝑦𝑘𝑘−𝑖𝑖+1)i=1,2,…,m。计算并保存𝑠𝑠𝑘𝑘=𝑥𝑥𝑘𝑘+1−𝑥𝑥𝑘𝑘,𝑦𝑦𝑘𝑘=𝑔𝑔𝑘𝑘+1−𝑔𝑔𝑘𝑘,随后更新𝐻𝐻𝑘𝑘+1:𝐻𝐻𝑘𝑘+1=(𝐼𝐼−𝜌𝜌𝑘𝑘𝑠𝑠𝑘𝑘𝑦𝑦𝑘𝑘𝑇𝑇)𝐻𝐻𝑘𝑘(𝐼𝐼−𝜌𝜌𝑘𝑘𝑦𝑦𝑘𝑘𝑠𝑠𝑘𝑘𝑇𝑇)+𝜌𝜌𝑘𝑘𝑠𝑠𝑘𝑘𝑠𝑠𝑘𝑘𝑇𝑇=𝑉𝑉𝑘𝑘𝑇𝑇𝐻𝐻𝑘𝑘𝑉𝑉𝑘𝑘+𝜌𝜌𝑘𝑘𝑠𝑠𝑘𝑘𝑠𝑠𝑘𝑘𝑇𝑇其中,𝜌𝜌𝑘𝑘=1𝑦𝑦𝑇𝑇𝑠𝑠𝑘𝑘,𝑉𝑉𝑘𝑘为n×n的矩阵。令k=k+1,返回2)。𝑘𝑘通常𝐻𝐻0被初始化为单位阵,由于𝐻𝐻𝑘𝑘具有正定性,因此算法是收敛的。4.1.4实验结果与分析本文将上述基于B样条函数的非刚性配准算法在VC++平台上进行了编程实现,并在上一章对于脑图谱进行分段仿射配准的基础上对其进行了非刚性配准。为方便起见,本章中使用“配准前图谱”指代使用第三章的算法所得到的粗配准图谱,使用“配准后图谱”指代粗配准图谱经过本章算法所得的结果。实验中,所采用的三维网格为5×5×5网格,算法迭代次数为500次,代价函数收敛因子为1×107。52 东南大学硕士学位论文4.1.4.1非刚性配准效果本文使用四组临床数据对分别与对应的粗配准脑图谱进行非刚性配准,结果如下:图4-1测试数据1:临床图像(左)配准前图谱(中)配准后图谱(右)图4-2测试数据2:临床图像(左)配准前图谱(中)配准后图谱(右)图4-3测试数据3:临床图像(左)配准前图谱(中)配准后图谱(右)53 第四章基于B样条函数的非刚性配准图4-4测试数据4:临床图像(左)配准前图谱(中)配准后图谱(右)根据通过观察上述图像可以看到,通常患者大脑的左右部分以及颅骨都不是对于脑中膈面完全对称的,但是患者的𝐿𝐿𝐿𝐿−𝑃𝑃𝐶𝐶以及𝐿𝐿𝑅𝑅−𝑃𝑃𝐶𝐶两个距离通常不会有很大的差异,导致通过分段仿射变换所得到的粗配准脑图谱通常对于中膈面是近似对称的。在对粗配准图像进行非刚性配准后,从图谱的颅骨外轮廓以及左右半脑的对称情况可以观察到配准后的图像发生了形变,且其趋势与患者图像相同。为了进一步观察配准算法的效果,本文分别将临床数据与其对应的配准前后图谱进行了分块融合,同时对算法执行前后的脑图谱作差,其结果如下图所示:图4-5测试数据1:图像与配准前(左)及配准后(中)图谱融合结果算法执行前后脑图谱之差(右)54 东南大学硕士学位论文图4-6测试数据2:图像与配准前(左)及配准后(中)图谱融合结果算法执行前后脑图谱之差(右)图4-7测试数据3:图像与配准前(左)及配准后(中)图谱融合结果算法执行前后脑图谱之差(右)图4-8测试数据4:图像与配准前(左)及配准后(中)图谱融合结果算法执行前后脑图谱之差(右)根据以上融合结果可以看出,在粗配准后的图谱数据与临床数据在细节上仍有局部55 第四章基于B样条函数的非刚性配准失配的现象存在,这种现象在颅骨部分以及脑皮质边缘可以明显的观察到,而在经过本章所使用的非刚性配准算法后,图谱与临床图像在颅骨的走势以及内部的纹理上更加趋于一致。本文同时对配准前后的图谱作差,所得图像中灰度值越高的部分代表二者之间的差异越大(图4-5~图4-8红色箭头所指处),此方法可以直观的看出非刚性算法的形变效果,其结果显示,灰度值较高的部分重点集中于颅骨、脑皮质边缘以及脑室部分。因此,基于B样条函数的非刚性配准算法对于这部分具有明显边界的区域具有较大的影响。随后,本文计算了粗配准后及精配准后图谱与图像之间的归一化互信息测度,对配准效果进行量化分析,其结果如表4-1所示:表4-1粗配准后及精配准后图谱与图像的归一化互信息值数据编号粗配准后归一化互信息精配准后归一化互信息差值11.11201.12700.015021.10591.11660.010731.11711.13840.021341.12111.13490.0138根据此表可知,对于配准前,图谱与图像的归一化互信息均值为1.1140,而配准后二者的归一化互信息值为1.1292,配准前后二者的归一化互信息值增长了约1.36%。4.1.4.2不同加权条件下MR图像的配准效果DBS手术的MR图像通常会涉及到多种加权方式。本文所使用的图像数据中包括了三种加权,分别为:T1WI、T2WI以及SWI。为了验证本章算法对于不同加权条件下的数据的配准效果,本文将粗配准后的图谱数据与同一患者的三种加权方式下的MR图像数据进行了配准。结果如下:56 东南大学硕士学位论文图4-9T1WI下原始图像(左)图谱配准结果(右)图4-10T2WI下原始图像(左)图谱配准结果(右)57 第四章基于B样条函数的非刚性配准图4-11SWI下原始图像(左)图谱配准结果(右)图4-12不同加权方式下图像与其配准结果分块融合效果分别对T1WI及SWI加权方式的配准结果的脑部区域进行反白,并对三者进行调窗后与原图像进行分块融合,其结果见图4-12,从中可以看出,本章所使用的算法对于不同加权条件下的MR图像均有较好的效果。随后计算三种加权条件下的配准前后图谱与图像的归一化互信息值,其结果如表4-2所示:表4-2不同加权方式下配准前后图谱与图像的归一化互信息值加权方式配准前归一化互信息配准后归一化互信息差值T1WI1.1121.1270.015T2WI1.12911.14950.0204SWI1.12011.12640.0063根据上表可知,在本章算法对于不同加权方式图像的配准均有一定效果。由于上述58 东南大学硕士学位论文测试使用的图谱为T2加权的ICBM-152图谱,因此该算法对于T2加权的图像配准结果较好。然而,在实际的图谱与临床数据的对比中仍然可以发现处于脑皮质的边缘部分有一定程度的失配现象(图4-12中标识部分),因此需要对算法进行一定程度的改进。4.2多层B样条自由形变配准由于三维图像的三阶B样条非刚性配准算法配准参数较多,运算所需的内存较大,为了平衡运算效率与最终结果,本文使用了5×5×5的网格对图谱与临床数据进行了配准。由上一节的实验结果可以看到,三维非刚性配准对于图谱与临床图像具有一定的效果,但是仍有部分失配现象,其主要原因是B样条配准中所使用的空间网格数较少导致的。本文为解决该问题使用多层B样条自由形变优化。4.2.1不同网格尺度的形变特点以二维数据为例,B样条非刚性配准算法主要使用控制网格来改变图像的形变状态。当控制网格的尺度较大时,其效果近似于一个全局配准;当控制网格的尺度较小时,接近于在控制网格附近的局部配准。如图4-13所示,在二维配准过程也可以近似看作一个拟合过程,由其拟合结果可以看出,网格较大时,所得结果误差较大,引起形变的范围较大(图4-13(b)),曲面较为光滑;随着网格的缩小,拟合结果的误差越小,引起形变的范围也被限制在了较小的区域内(图4-13(d))。图4-13不同尺度下B样条函数拟合结果[56](a)原始数据(b)8×8网格(c)16×16网格(d)32×32网格59 第四章基于B样条函数的非刚性配准4.2.2多层B样条自由形变优化假设Φ0,Φ1,Φ2,…,Φ𝑛𝑛是区域Ω中的控制网格,且每一层控制网格的控制点间距𝐿𝐿Φ𝑛𝑛1为上一层的一半,即𝐿𝐿Φ𝑘𝑘+1=2𝐿𝐿Φ𝑘𝑘。若使用第一层控制网格对图像进行非刚性配准,可以得到一个较为平滑的配准结果。若令对应的B样条函数为𝑓𝑓0,假设理想状态下图像的形变结果为𝑅𝑅𝐼𝐼,其误差可表示为:∆1=𝑅𝑅𝐼𝐼−𝑓𝑓0(4.12)随后使用更小尺度的控制网格Φ1对第一层所得的配准结果进行进一步的优化,令所得的B样条函数为𝑓𝑓1,其误差可表示为:∆2=𝑅𝑅𝐼𝐼−𝑓𝑓0−𝑓𝑓1(4.13)同理,使用尺度为𝐿𝐿Φ𝑘𝑘的控制网络Φ𝑘𝑘,对第k层所得的结果进行非刚性配准,设所得的B样条函数为𝑓𝑓𝑘𝑘,其误差可表示为:∆=𝑅𝑅−∑𝑘𝑘𝑓𝑓=∆−𝑓𝑓(4.14)𝑘𝑘+1𝐼𝐼𝑖𝑖=0𝑖𝑖𝑘𝑘𝑘𝑘由上式可以看出,随着非刚性配准所使用的网格层数的增加,配准所带来的误差在每一层中都被缩小,最终非刚性配准的结果可以表示为所有层级的样条插值函数之和,即:f=∑𝑛𝑛𝑓𝑓(4.15)𝑘𝑘=0𝑘𝑘因此多层B样条自由形变整个流程为,首先使用大尺度的控制网格得到一个全局近似解,随后对于整幅图像使用更小尺度的网格对其进行优化,消除残差,其优点是可以通过这个过程达到较高的精度,同时,防止直接使用小尺度网格导致算法陷入局部极值内。4.2.3实验结果与分析本文对于上一节三维B样条自由形变结果进行了多层B样条的优化配准。在实际操作中,由于三维非刚性配准算法需要占用较大的计算内存,初始的三维B样条空间控制网格数为5×5×5,其内存占用约为16GB,耗费了较大的系统资源,因此本文在初始的三维B样条自由形变配准的基础上,对三维数据逐层使用二维小尺度控制网格的B样条自由形变模型进行优化。对三维非刚性配准后的临床数据作为初始数据,分别对数据进行对其继续进行8×8,10×10,12×12,14×14,16×16及32×32的网格优化,其结果如下图所示:60 东南大学硕士学位论文图4-14测试数据1不同网格数下的配准结果(a)临床数据(b)三维配准图谱(c)8×8(d)10×10(e)12×12(f)14×14(g)16×16(f)32×32图4-15测试数据2不同网格数下的配准结果(a)临床数据(b)三维配准图谱(c)8×8(d)10×10(e)12×12(f)14×14(g)16×16(f)32×3261 第四章基于B样条函数的非刚性配准图4-16测试数据3不同网格数下的配准结果(a)临床数据(b)三维配准图谱(c)8×8(d)10×10(e)12×12(f)14×14(g)16×16(f)32×32图4-17测试数据4不同网格数下的配准结果(a)临床数据(b)三维配准图谱(c)8×8(d)10×10(e)12×12(f)14×14(g)16×16(f)32×3262 东南大学硕士学位论文观察配准后的图像可以发现,当网格数较高时,配准的结果反而使图谱脑白质部分呈现丝状的扭曲,且在32×32网格的配准结果中最为明显(各组图像中红色箭头标识)。这与实际临床数据的特征不符,但二者的归一化互信息值却在一个较高的水平。这是因为使用较多的网格数进行配准时,每个网格所控制的区域较少,配准结果陷入了局部极值。对于各组数据的配准结果,分别计算其与临床数据的归一化互信息值,其结果如表4-3所示:表4-3不同网格数下图谱与图像配准归一化互信息值数据粗配三维8×810×1012×1214×1416×1632×32序号准配准网格网格网格网格网格网格11.12901.14981.15731.16061.16971.17271.17321.186421.11061.12541.13151.14161.14461.15771.16111.173331.15451.16721.16981.17521.17831.18291.19001.202041.12041.13881.14211.14851.15431.15971.16471.1726根据不同网格数下的归一化互信息值,可以得到如图4-18的变化曲线。由该曲线可以看出,当网格数较少时,略微增加网格数量就可以使得图像的归一化互信息值得到较大提升,而当网格数提高到一定数值之后(≥14~16),即使对空间网格的数量进行翻倍,归一化互信息的值的增加数量也较为有限。1.211.21.191.181.171.161.15归一化互信息值1.141.131.12581012141632B样条网格数数据1数据2数据3数据4图4-18不同网格数下归一化互信息变化曲线同时,在实验过程中,记录程序所占用的内存情况,如表4-4所示。可以看到二维网格相比于三维网格占用内存小,但是当网格数量达到32×32时,其内存占用会达到7554MB,会占用较多的系统资源。63 第四章基于B样条函数的非刚性配准表4-4不同网格数下程序所占内存网格数81012141632所占内存(MB)400583860138115417554结合以上两个因素,以8×8的空间网格数作为多层B样条自由形变配准的第二层网格参数较为合适。为了验证算法的效果将分别将配准前、第一层配准后以及第二层配准后的图谱与临床数据进行融合,并图谱进行直方图均衡,结果如下。图4-19测试数据1与临床数据融合结果(a)配准前(b)第一层配准后(b)第二层配准后图4-20测试数据2与临床数据融合结果(a)配准前(b)第一层配准后(b)第二层配准后64 东南大学硕士学位论文图4-21测试数据3与临床数据融合结果(a)配准前(b)第一层配准后(b)第二层配准后图4-22测试数据4与临床数据融合结果(a)配准前(b)第一层配准后(b)第二层配准后通过观察上述图像可以发现,配准前临床数据与原始图谱有较多失配之处(红色箭头标识),而在进行第一层配准后,图像中仍有部分不连续之处(黄色箭头标识),其中较为明显的是颅骨以及脑皮质边缘部分,而经过第二层配准之后二者的融合图像在边缘上更为连续(白色箭头标识)。因此本文所使用的基于三维5×5×5网格+二维8×8网格双层结构的多层B样条自由形变配准相比单层的三维B样条非刚性配准算法可以更好的解决图谱与临床图像的精确配准问题。4.3本章小结本章主要对图谱与临床图像的非刚性配准进行了实现。首先介绍了非刚性图像配准的基本框架,基于此框架,使用图像灰度域作为特征空间、互信息作为相似性测度、B样条形变模型作为搜索空间、L-BFGS优化算法作为搜索策略来进行精细配准工作。通65 第四章基于B样条函数的非刚性配准过多组实验数据测试,基于B样条自由形变模型的配准算法可以较好解决图谱与临床图像的精确配准问题,其对于多种加权方式的MR数据均有较好的鲁棒性。针对算法结果中仍然存在的失配问题,基于三维非刚性形变算法的结果,对配准模型进行了多层B样条自由形变配准优化。通过实验验证了不同二维网格参数对输出的影响,确定了三维5×5×5网格+二维8×8网格双层结构的多层B样条自由形变配准算法具有较好的效果该结构是较适合解决本文问题层次结构。66 第五章面向帕金森深部脑刺激手术的系统应用第五章面向帕金森深部脑刺激手术的系统应用前面的章节介绍了DBS术前图像与ICBM-152图谱的配准方法,较好的解决了二者的配准问题。本章根据深部脑刺激手术不同阶段数据的特点,将上述方法应用于帕金森病DBS术前靶点识别以及术后效果评估的环节中,并设计了一套标准的配准流程,将配准获得DBS术前靶点位置以及DBS术后电极与丘脑底核相互位置关系进行了可视化。5.1临床数据及扫描参数目前临床上一般使用高分辨率核磁共振在T1WI加权扫描方式下获得患者术前的全脑扫描数据,并使用T2加权扫描方式得到患者脑深部的核团位置信息,用于患者大脑的总体情况评估,判断电极的植入路径等。在电极完成植入后,需要对植入的情况进行评估时,考虑到高场强、长时间的MR扫描会导致金属材质的刺激电极发热,造成患者植入路径附近脑组织发生水肿等不良反应,因此在术后成像过程中通常使用CT或低分辨MRI对患者头部进行扫描。本文临床数据源于江苏省人民医院帕金森病患者的双侧丘脑底核脑深部电刺激手术,植入电极为美国美敦力公司3389电极。手术前,使用德国西门子公司3.0TMRI系统,在快速自旋回波序列T1相重复时间(timeofrepetition,TR)为1600ms,回波时间(timeofecho,TE)2.48ms,可视范围(fieldofview,FOV)250mm×250mm,像素512×512,层厚1.5mm条件下以矢状位无间隔扫描获取全脑数据;随后在快速自旋回波序列T2相TR5000ms,TE84.00ms,FOV229mm×229mm,像素384×384,层厚2mm条件下以横断位扫描获取核团部位数据。手术后,使用美国通用电气公司1.5TMRI系统,在快速自旋回波T1相TR4000ms,TE105.62ms,FOV250mm×20mm,像素512×512,层厚2mm条件下以横断位扫描获得带有电极的MRI图像。67 东南大学硕士学位论文图5-1(a)术前T1矢状位影像(b)术前T2横断位影像(c)术后T2横断位影像由于实际成像条件限制,术前T2加权MR序列的扫描时间较长,通常只对患者的植入核团区域进行薄层扫描代替全脑扫描。术后扫描时为了减少扫描时间,同样只对患者核团附近的脑部区域进行薄层扫描。因此确定核团位置的T2加权MR数据以及用于确定电极位置的术后MR数据都是部分脑部区域图像,仅术前T1加权序列是全脑数据。因此利用前文配准方法对脑部术前术后的核团和电极位置定位较为困难,需要基于上述DBS手术影像数据特点,设计一套有针对性的数据处理流程。5.2数据处理流程设计由第四章结果可知,在不同的加权条件下,MR对大脑的不同组织的增强效果不同,配准的效果也不同,因此在理想状况下,使脑图谱与包含丘脑底核信息的T2WI加权图像进行配准,会取得更好的靶点识别效果。然而该图像并不是全脑图像,该图像中缺少完整的AC、PC以及脑皮质层标识点信息,因此无法使用第二章的算法对图像进行空间标准化,也无法使用第三章中的算法实现图谱与脑图像的粗配准。而术前的T1WI全脑图像提供了完整的AC、PC以及脑皮质层标识点信息,弥补了这个不足。因此,我们可以使用第二章以及第三章的算法处理T1WI数据,得到与临床图像粗配准的图谱数据,随后利用该数据与T2WI图像进行第四章所介绍的非刚性配准过程,得到精确配准的图谱数据。然而,T1WI与T2WI数据是在不同时域以及空域得到的,缺少相应的空间对应关系,因此二者在粗配准到精确配准之间还存在着一个数据空间统一问题。刚性配准算法相对于非刚性配准算法起步较早,发展已经相对成熟[57,58],对于患者不同模态间的脑图像配准问题,可以使用刚性配准算法解决[59]。对于术后数据,同样存在着与术前数据的数据空间转换问题,同样可以使用刚性配准算法与术前T1进行配准。随后,考虑到术后大脑发生形变的因素,因此需要对术后数据与术前数据进行非刚性配准,计算形68 第五章面向帕金森深部脑刺激手术的系统应用变的变化量,并将其应用到术前精确配准的图谱上,从而得到术后图像中的核团估计。根据上述分析,本文设计了如图5-2的数据处理流程:1)将术前T2数据以及术后数据与术前T1数据进行刚性配准,将使三者的坐标空间达到一致;2)利用术前T1数据为全脑数据的特性,通过使用第二章所介绍的预处理过程,计算空间转换矩阵,将其同时应用至刚性配准后的术前T2数据以及术后数据,使三种数据的坐标空间与标准的ICBM-152坐标空间相同;3)找出术前T1数据的脑皮层标识点将ICBM-152图谱与其进行分段仿射配准;4)截取术前T2中的有效部分,将分段仿射配准后的脑图谱与T2图像进行非刚性配准,得到患者的个性化图谱;5)根据ICBM-152图谱的核团标记,获得术前图像中的核团预测;6)手术的开颅以及电极植入过程会导致患者大脑发生形变,因此,术前术后的数据需要进行进一步非刚性配准。将术前的核团位置使用相同的B样条自由形变参数变换至术后数据中;7)提取术后数据中的电极路径,就可以获得核团与电极的相互位置关系;图5-2数据处理流程总体框图5.3数据处理流程的应用本文在VisualStudio2010开发平台上,结合ITK的基本配准框架[60]对上述的流程进行了实现,计算机的硬件配准为:IntelXeonE3-1231v3CPU,16GB内存,GTX970M69 东南大学硕士学位论文显卡,操作系统为Windows1064位。以下为影像处理管线的各阶段结果。DBS术前以及术后数据的断层如图5-3所示。从图中我们可以看到,(a)为T1全脑图像,初始位置为矢状位;(b)为T2核团区域图像,其初始位置为横断位,且半球间裂隙出现了略微的倾斜;(c)为术后核团区域数据,其初始位置为横断位,半球间裂隙较呈标准的垂直状态。根据这三幅断层数据可知三份原始体数据的初始空间位置不同,其主要原因是患者在不同时空条件下很难保持其头部在图像中的相对位置不变。图5-3DBS术前术后数据的所在空间比较(a)术前T1加权图像(b)术前T2加权图像(c)术后图像图5-4为将术前T2加权数据以及术后T1加权数据向术前T1全脑数据进行刚性配准后的效果,其中(a)为术前数据的配准效果,(b)为术后数据的配准效果,从中间的融合效果可以看出,二者在空间位置上基本达到了一致。图5-4(a)术前T2加权MR与术前T1加权MR刚性配准结果(b)术后MR图像与术前T1加权MR刚性配准结果70 第五章面向帕金森深部脑刺激手术的系统应用通过计算术前T1-术前T2,术前T1-术后图像以及术前T2与术后图像的各组数据在配准前后的归一化互信息值(表5-1),从结果可以看出,在配准后,各组图像的归一化互信息值均有所上升。其中,术前T1-术前T2以及术前T1-术后的配准后归一化互信息值增加幅度较小,其原因为术前T2数据以及术后数据的数据量较少,背景占了图像中的很大一部分,从配准后的情况可以看到,二者的重合部分比例较低,这导致了图像的联合灰度直方图分布十分分散,因此归一化互信息值普遍较低。而对于二者之间的相互配准,由于二者的重合部分比例较高,联合直方图分布较为集中,因此,归一化互信息的增长值较大。表5-1配准前后各组数据的归一化互信息值组别配准前归一化互信息值配准后归一化互信息值差值术前T1-术前T21.02671.03280.0061术前T1-术后1.01681.04250.0257术前T2-术后1.0311.29660.2656对术前图像及脑图谱使用上述数据流程进行处理,所得到的手术靶点识别结果如图5-5所示。图5-5术前核团识别结果(a)专家手动标记的结果(b)配准流程所得结果表5-2为二者核团重心坐标计算结果。对于左侧核团,算法所得左侧核团与金标准的误差距离为0.22mm,右侧核团与金标准的误差距离为1.29mm,算法的平均误差为0.76mm。71 东南大学硕士学位论文表5-2专家预测结果与配准管线所得术前核团重心坐标左侧核团重心坐标(mm)右侧核团重心坐标(mm)专家手动标记(-9.04,126.30,9.66)(14.47,126.20,9.70)配准管线结果(-8.72,126.38,9.24)(13.42,123.79,10.10)在电极植入后,大脑会发生形变。术前大脑额叶与颅骨之间的间隙较小,使用术后刚性配准后的同一断面数据进行对比(图5-6),可以很容易的观察到,额叶与颅骨之间的间隙明显的增大了。因此,对于将术前图像中的核团预测直接用于术后数据中会引入误差。图5-6术前(a)术后(b)大脑形变对比因此本文对于手术前后的图像进行了非刚性配准,并对二者作差,查看形变情况。其结果如图5-7所示,观察图像的差值可以发现其脑组织前部的差值较大,而脑组织内部也有一定程度的形变,因此,配准结果与我们的预期相符。图5-7术前图像非刚性配准结果分析(a)原始图像(b)配准结果(c)二者之差随后将配准所得形变模型使用至图谱配准结果,利用其子图谱得到术后核团位置的预测结果,并于专家手动标记的结果进行对比,见图5-8:72 第五章面向帕金森深部脑刺激手术的系统应用图5-8术后核团位置估计结果(a)专家手动标记结果(b)配准流程所得的结果计算二者核团的重心坐标,其结果见表5-3。算法所得左侧核团与专家手动标记结果的误差距离为0.97mm,右侧核团与金标准的误差距离为2.57mm,算法的平均误差为1.77mm。表5-3专家预测结果与配准管线所得术后核团重心坐标左侧核团重心坐标(mm)右侧核团重心坐标(mm)专家手动标记(-8.97,126.26,9.69)(14.45,126.16,9.18)配准管线结果(-8.62,125.81,8.91)(13.16,123.98,9.65)最后对电极路径进行手动提取,并进行三维可视化,其结果如图5-9:图5-9核团与电极的三维可视化效果73 东南大学硕士学位论文5.4本章小结本章利用前面各个章节所提出的脑图谱与DBS手术图像的配准算法,结合江苏省人民医院提供的术前全脑T1WI、薄层T2WI图像以及术后薄层MR图像,利用术前T1WI数据的全脑信息,将图谱与其进行粗配准。将含核团信息的T2WI数据与粗配准后所得图谱进行非刚性配准,得到精配准图谱以及术前核团分割结果。对于术后数据,考虑大脑发生形变因素,计算了术前术后数据的形变模型参数,将其应用至精配准图谱中,得到术后核团预测,并根据术后图谱中的电极信息,将结果进行了可视化。通过将数据处理流程所得结果与专家手动勾画的金标准进行对比,算法术前靶点平均误差为0.76mm,术后靶点平均误差为1.77mm,对于临床判断手术前后靶点位置具有一定的指导意义。74 第六章总结与展望第六章总结与展望6.1总结随着我国人口老龄化问题的日益严重,我国帕金森患者的数量也在不断增加,深部脑刺激手术对于帕金森后期患者来说是最有效的治疗手段。其中,术前对电极植入靶点的精确导航,以及术后对于电极与靶点之间的相互位置关系的准确判断都是影响DBS手术效果的关键因素。由于目前成像方面的技术条件限制,DBS手术图像中的靶点难以识别,图像模态较多,且信息不全,通常的做法是引入脑图谱来辅助解决术前导航以及术后判断电极与核团位置关系的问题。图谱与临床图像的配准问题是这一过程中最为关键的一环,本文以ICBM-152图谱作为标准脑图谱,就其与DBS手术中的临床数据的配准问题展开了相关研究,主要完成工作如下:(1)在成像过程中,患者的头部摆放存在着随机性,临床数据与标准图谱通常不在同一坐标系中,因此需要将临床数据的坐标空间进行标准化。在此过程中,首先需要对坐标空间中最重要的基准面——脑中膈面进行识别。本文实现了基于K-L距离的脑中膈面识别算法,在此基础上将中膈面旋转至标准的矢状位,并使AC-PC连线呈水平状态的位置。实验结果表明,脑中膈面算法对于不同倾斜角度的模拟数据有较好的校准效果,应用于临床数据的效果也比较理想,满足对临床数据坐标空间标准化的需求。(2)由于脑图谱与临床数据在大小以及形状上存在着一定的差异,很难简单使用现有算法将二者配准。本文借助脑图谱的大脑分块概念,将临床数据通过脑皮层标识点分为12个区间。随后采用基于FCM的聚类算法,根据图像中不同组织的灰度进行聚类,将聚类结果进行二值化投射,并使用形态学图像处理提取脑组织边界,得到标识点坐标。在此基础上,使用分段仿射变换作为配准算法,计算脑图谱对应区间三个维度上与临床数据的缩放比例,将图谱粗配准至与临床数据外形大小相同的状态。实验结果表明,FCM聚类算法对于多组临床数据的脑皮质层标识点提取都具有良好的效果,分段仿射配准算法最后的配准效果较为理想,配准前后图谱与临床图像的互信息值得到了提高。(3)在进行图谱与数据的粗配准后,二者之间仍然存在着局部失配的现象,本文使用基于图像的灰度域为特征空间;互信息为相似性测度;B样条形变模型为搜索空间;L-BFGS优化算法为搜索策略的非刚性配准算法,进一步提高二者的匹配程度。实验结果表明,基于B样条自由形变模型的配准算法可以解决图谱与临床图像的精确配准问75 东南大学硕士学位论文题,且对于多种加权方式的MR数据均具有较好的鲁棒性。针对算法结果中仍然存在的失配问题,对配准模型进行了多层B样条自由形变配准优化,验证了不同二维网格参数对输出的影响。针对DBS手术图像,确定了三维5×5×5网格+二维8×8网格双层的最优B样条自由形变配准结构。(4)结合DBS手术不同阶段的临床图像设计了针对DBS手术的数据处理流程。首先,利用术前T1WI数据的全脑信息,将图谱于其进行粗配准;随后,将含核团信息的T2WI数据与粗配准后所得图谱进行非刚性配准,得到精配准图谱以及术前核团分割结果;最后,对于术后数据,考虑大脑发生形变因素,计算了术前术后数据的形变模型参数,将其应用至精配准图谱中,得到术后靶点预测,并根据术后图谱中的电极信息,将结果进行了可视化。实验结果表明,该流程所得的术前术后核团位置估计结果与专家提供的估计结果相近,可以作为临床诊断的辅助参考。6.2展望本文针对ICBM-152图谱与DBS手术图像配准的问题进行了讨论研究,设计实现了相应的配准方案,并利用实验验证了算法的有效性。由于时间以及实验条件的限制,与本研究相关的工作没有完全地展开并实施,主要包括以下几个方面:(1)在临床数据进行空间标准化的过程中,由于整个图谱空间是建立在AC、PC两点的基础上的,这两点相当于整个坐标空间的基准点,因此AC以及PC点的确定是一个十分重要的过程。而本文在具体的实施过程中采用的是在MSP平面上通过手动标记来确定AC以及PC点,在这个过程中可能因为操作者的主观因素引入误差,因此可以考虑使用算法识别数据中的AC、PC点,实现整个处理管线的完全自动化。(2)本论文采用的三维非刚性配准算法,存在运算速度过慢以及对于内存需求较大的弊端,导致所生成的个性化图谱局限在目标核团周围的区域。因此,对于算法的运算时间以及空间复杂度的优化也是后续工作开展的重点之一,可以考虑使用CUDA等并行化技术对算法进行优化。(3)本文所得到的电极与核团的相互位置关系仅仅是一个初步的应用,术后程控过程中,不同的电极型号以及刺激参数的选择对于植入靶点周围的刺激程度以及范围也不同。后续的研究工作将基于手术图像,开展对电极刺激范围等方面的研究工作。76 致谢致谢首先感谢我的导师罗守华老师,本论文是在罗老师的悉心指导下完成的,他深厚的学术功底、丰富的工程经验令人钦佩。他的治学方法、工作态度以及对学生的感化教育让我受益匪浅,促进了我的成长,为我今后的学习和工作带来积极的影响。罗老师不仅在学习上给予我诸多指导,在生活上对我也非常关心,对此我十分感谢。感谢江苏省脑科医院的郑慧芬主任和江苏省人民医院的曹胜武主任,他们在我完成课题的过程中提供了很多帮助,并用丰富的临床经验为本课题的研究指明了方向。同时,对于我所提出的问题,他们总是耐心的给予回答,在此对他们表示深深地感谢。感谢师兄董歌、李光、张奎以及马孚骁,平时在一间办公室工作,他们为我的毕业设计工作提供了很好的建议,耐心解答我不懂的疑问,并在我刚进实验室时给了我耐心的指导。感谢我的女朋友朱珠,熬夜协助我完成了毕业设计中展示效果图的制作。感谢赵海桐师妹,在我完成毕业设计的过程中给了我很多有用的建议,感谢实验室的李中源师兄、吴沛泽师兄、孙翌师兄、沈涛师兄、高晟凯师兄、吴华珍师姐和姚佳丽师姐,他们给我传授了很多经验,并且给予我很多帮助。感谢骆爽、贺凯、许海涛、顾晓卉、蒋建慧等实验室同学,大家在学习生活中相互交流,分享经验,促进共同进步,感谢大家的支持与配合。感谢李昊、贾浩然、陈良等同学,感谢他们在学习和生活上对我的关心,让我在繁忙的学习工作之余能够感受到集体生活的温暖。感谢父母和家人对我无私的关爱、理解和支持,他们是我坚强的后盾。感谢所有关心、帮助过我的人。77 参考文献参考文献[1]ZhangZX,RomanGC,HongZ,etal.Parkinson'sdiseaseinChina:prevalenceinBeijing,Xian,andShanghai[J].Lancet,2005,365(9459):595-597.[2]EisensehrI,LinkeR,NoachtarS,etal.Reducedstriataldopaminetransportersinidiopathicrapideyemovementsleepbehaviourdisorder.ComparisonwithParkinson'sdiseaseandcontrols[J].Brain,2000,123(6):1155-1160.[3]KayeJ,GageH,KimberA,etal.ExcessburdenofconstipationinParkinson'sdisease:apilotstudy[J].MovementDisorders,2006,21(8):1270–1273.[4]SiddiquiMF,RastS,LynnMJ,etal.AutonomicdysfunctioninParkinson'sdisease:acomprehensivesymptomsurvey[J].Parkinsonism&RelatedDisorders,2002,8(4):277-284.[5]KlugerBM,KlepitskayaO,OkunMS.Surgicaltreatmentofmovementdisorders[J].NeurolClinics,2009,27(3):633-677.[6]高翔.老龄化时代:来势汹汹的帕金森病[J].健康管理,2013,7:116-119.[7]钱浩.丘脑底核脑深部电刺激(DBS)治疗帕金森病的临床初探和基础研究[D]:[硕士学位论文].广州:中山大学,2010.[8]王琼,韩丁,陈彤,等.帕金森病运动并发症的调查分析[J].中华老年心脑血管病杂志,2013,15(4):390-392.[9]KopellBH,RezaiAR,ChangJW,etal.Anatomyandphysiologyofthebasalganglia:ImplicationsfordeepbrainstimulationforParkinson'sdisease[J].MovementDisorders,2006,21(14):238–246.[10]GabrielEM,NasholdBS.Evolutionofneuroablativesurgeryforinvoluntarymovementdisorders:anhistoricalreview[J].Neurosurgery,1998,42(3):590-591.[11]RezaiAR,KopellBH,GrossRE,etal.DeepbrainstimulationforParkinson'sdisease:Surgicalissues[J].MovementDisorders,2006,21(14):197–218.[12]底丘脑.百度百科.[EB/OL].[2016-10-27]https://baike.baidu.com/item/底丘脑/1966279.[13]程振国,孙来广,高俊红,等.脑深部核团电刺激与毁损术治疗帕金森病疗效比较[J].中国实用医刊,2007,34(18):34-35.[14]GuridiJ,RodriguezorozMC,ManriqueM.SurgicaltreatmentforParkinson'sdisease[J].Neurocirugia,2004,15(1):5-16.[15]AshbyP,RothwellJC.Neurophysiologicaspectsofdeepbrainstimulation[J].Neurology,2000,55(6):17-20.[16]BenabidAL,KrackPP,BenazzozA,etal.DeepbrainstimulationofthesubthalamicnucleusforParkinson'sdisease:methodologicaspectsandclinicalcriteria[J].Neurology,2000,55(6):40-44.[17]GrillWM,MclntyreCC.Extracellularexcitationofcentralneurons:implicationsforthemechanismsofdeepbrainstimulation[J].Thalamus&RelatedSystems,2001,1(3):269-277.[18]统一帕金森氏病评分量表.临床医学[EB/OL].http://www.1mpi.com/doc/d54ae1ac36a4b8-7f612e947d.[19]CuiX.BrainAC-PCline[EB/OL].[2008-11-21]http://www.alivelearn.net/?p=492.[20]DormontD,CornuP,PidouxB,etal.Chronicthalamicstimulationwiththree-dimensionalMR78 东南大学硕士学位论文stereotacticguidance[J].AmericanJournalofNeuroradiology,1997,18(6):1093-1107.[21]VillegerA,OuchchaneL,LemaireJJ,etal.Assistancetoplanningindeepbrainstimulation:datafusionmethodforlocatinganatomicaltargetsinMRI[C].ProceedingsoftheInternationalConferenceoftheIEEEEngineeringinMedicine&BiologySociety.2006:144.[22]HornA,KühnAA.Lead-DBS:atoolboxfordeepbrainstimulationelectrodelocalizationsandvisualizations[J].Neuroimage,2014(107):127-135.[23]LiuY,DawantBM.Multi-modallearning-basedpre-operativetargetingindeepbrainstimulationprocedures[C].ProceedingsoftheIEEE-EMBSInternationalConferenceonBiomedical&HealthInformatics,2016:17.[24]CastroFJ,PolloC,MeuliR,etal.Acrossvalidationstudyofdeepbrainstimulationtargeting:fromexpertstoatlas-based,segmentation-basedandautomaticregistrationalgorithms[J].IEEETransactionsonMedicalImaging,2006,25(11):1440-1450.[25]WalterU,WoltersA,WittstockM,etal.Deepbrainstimulationindystonia:sonographicmonitoringofelectrodeplacementintotheglobuspallidusinternus[J].MovementDisorders,2009,24(10):1538–1541.[26]PouratianN,ZhengZ,BariAA,etal.Multi-institutionalevaluationofdeepbrainstimulationtargetingusingprobabilisticconnectivity-basedthalamicsegmentation[J].JournalofNeurosurgery,2011,115(5):995-1004.[27]叶宇阳,曹纹平,曹胜武,等.神经导航技术和图像融合技术在帕金森病脑深部电刺激术中的应用[J].江苏医药,2017,43(1):59-61.[28]耿馨佚,何江弘,王守岩.脑深部电刺激术靶点定位方法的现状及进展[J].中华神经外科杂志,2015,31(9):964-967.[29]VideenTO,CampbellMC,TabbalSD,etal.Validationofafiducial-basedatlaslocalizationmethodfordeepbrainstimulationcontactsintheareaofthesubthalamicnucleus[J].JournalofNeuroscienceMethods,2008,168(2):275-281.[30]SilvaNMD,RozanskiVE,CunhaJPS.A3DmultimodalapproachtopreciselylocateDBSelectrodesinthebasalgangliabrainregion[C].ProceedingsoftheInternationalIEEE-EMBSConferenceonNeuralEngineering,2015:292-295.[31]郝斌.丘脑底核核磁共振三维重建与脑深部电刺激治疗帕金森病的相关研究[D]:[博士学位论文].上海:第二军医大学,2008.[32]张芷齐,耿馨佚,徐欣,etal.丘脑底核深部脑刺激三维空间结构可视化[J].生物医学工程学杂志,2016(3):405-412.[33]潘红.脑部MRI海马体三维分割算法研究[D]:[硕士学位论文].成都:电子科技大学,2015.[34]李屹.人脑计算机图谱区配与三维重建问题研究[D]:[硕士学位论文].南京:东南大学,1999.[35]CocoscoC,KollokianV,KwanR,etal.BrainWeb:OnlineInterfacetoa3DMRISimulatedBrainDatabase[J].NeuroImage,1997,5(4):425.[36]KeithA.TheWholeBrainAtlas[EB/OL].http://www.med.harvard.edu/aanlib/.[37]TalairachJT,TournouxL.Co-PlanarStereotaxicAtlasOfTheHumanBrain:AnApproachToMedicalCerebralImaging[J].Thieme,1988(122).[38]VolkauI,AzizA,NowinskiW,etal.ExtractionofthemidsagittalplanefrommorphologicalneuroimagesusingtheKullback-Leibler'smeasure[J].MedicalImageAnalysis,2006,10(6):863-874.79 参考文献[39]GongBJ.罗德里格斯公式理解、推导[EB/OL].[2017-12-29].https://blog.csdn.net/q58395-6932/article/details/78933245.[40]安淑卿.基于B样条的Talairach图谱与脑图像的配准与融合方法研究[D]:[硕士学位论文].南京:东南大学,2008.[41]NowinskiWL.ModifiedTalairachlandmarks[J].Neurochir,2001,143(10):1045-1057.[42]常秋月.基于T2加权MR像的脑皮质标志点的定位方法研究[D]:[硕士学位论文].哈尔滨:哈尔滨工业大学,2007.[43]SerraJ.Imageanalysisandmathematicalmorphology[M].AcademicPress,1982,536-538.[44]GonzalezRC,WintzP.Digitalimageprocessing[M].Addison-Wesley,1977,484-486.[45]诸葛斌,冯焕清,周荷琴.医学图像体绘制中的快速三线性插值算法[J].航天医学与医学工程,2003,16(3):206-209.[46]BrownLG.Asurveyofimageregistrationtechniques[J].AcmComputingSurveys,1992,24(4):325-376.[47]ShannonCE.Amathematicaltheoryofcommunication[J].BellSystemTechnicalJournal,1948,27(3),379-423.[48]朱雪龙.应用信息论基础[M].北京:清华大学出版社,2001,36-52.[49]ViolaP,WellsWM.Alignmentbymaximizationofmutualinformation[J].InternationalJournalofComputerVision,1997,24(2):137-154.[50]李勇.基于互信息的图像拼接算法研究[D]:[硕士学位论文].哈尔滨:哈尔滨工业大学,2016.[51]MaesF,CollignonA,VandermeulenD,etal.Multimodalityimageregistrationbymaximizationofmutualinformation[J].IEEETransMedImaging,1997,16(2):187-198.[52]StudholmeC,HillDLG,HawkesDJ.Anoverlapinvariantentropymeasureof3Dmedicalimagealignment[J].PatternRecognition,1999,32(1):71-86.[53]BoorCD.OncalculatingwithB-splines[J].JournalofApproximationTheory,1972,6(1):50-62.[54]NocedalJ.Updatingquasi-Newtonmatriceswithlimitedstorage[J].MathematicsofComputation,1980,35(151):773-782.[55]李琦.多模态医学影像鲁棒配准方法研究[D]:[博士学位论文].西安:西安电子科技大学,2013.[56]李宇森.基于多层B样条和L2正则化的图像配准方法研究[D]:[硕士学位论文].济南:山东大学,2017.[57]VanDE,PolEJD,ViergeverMA.Medicalimagematching-areviewwithclassification[J].IEEEEngineeringinMedicine&Biology,1993,12(1):26-39.[58]陈珊.医学图像刚性配准与非刚性配准技术研究[D]:[硕士学位论文].杭州:浙江工业大学,2015.[59]李谭.多模态医学图像配准算法研究[D]:[硕士学位论文].天津:天津大学,2007.[60]IbanezL,SchroederW,NgL,etal.TheITKsoftwareguide[J].ComputationalStatistics&DataAnalysis,2005(21):231–256.80 作者简介作者简介倪杨阳(1993—),男,汉族,上海嘉定人,现为东南大学信号与图像处理实验室硕士研究生,研究方向为医学图像处理。攻读硕士学位期间发表的论文[1]倪杨阳,郑慧芬,曹胜武,罗守华,一种DBS术前图像与ICBM-152图谱的配准算法[J].中国医疗设备.(已录用)81

当前文档最多预览五页,下载文档查看全文

此文档下载收益归作者所有

当前文档最多预览五页,下载文档查看全文
温馨提示:
1. 部分包含数学公式或PPT动画的文件,查看预览时可能会显示错乱或异常,文件下载后无此问题,请放心下载。
2. 本文档由用户上传,版权归属用户,天天文库负责整理代发布。如果您对本文档版权有争议请及时联系客服。
3. 下载前请仔细阅读文档内容,确认文档内容符合您的需求后进行下载,若出现内容与标题不符可向本站投诉处理。
4. 下载文档时可能由于网络波动等原因无法下载或下载错误,付费完成后未能成功下载的用户请联系客服处理。
大家都在看
近期热门
关闭