CN102282587B - 用于处理图像的系统和方法 - Google Patents
用于处理图像的系统和方法 Download PDFInfo
- Publication number
- CN102282587B CN102282587B CN200980147491.XA CN200980147491A CN102282587B CN 102282587 B CN102282587 B CN 102282587B CN 200980147491 A CN200980147491 A CN 200980147491A CN 102282587 B CN102282587 B CN 102282587B
- Authority
- CN
- China
- Prior art keywords
- centerdot
- image
- generating unit
- time
- time domain
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
- G06T7/0012—Biomedical image inspection
- G06T7/0014—Biomedical image inspection using an image reference approach
- G06T7/0016—Biomedical image inspection using an image reference approach involving temporal comparison
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/0033—Features or image-related aspects of imaging apparatus classified in A61B5/00, e.g. for MRI, optical tomography or impedance tomography apparatus; arrangements of imaging apparatus in a room
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/05—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves
- A61B5/055—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves involving electronic [EMR] or nuclear [NMR] magnetic resonance, e.g. magnetic resonance imaging
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
- A61B6/02—Devices for diagnosis sequentially in different planes; Stereoscopic radiation diagnosis
- A61B6/03—Computerised tomographs
- A61B6/032—Transmission computed tomography [CT]
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
- A61B6/02—Devices for diagnosis sequentially in different planes; Stereoscopic radiation diagnosis
- A61B6/03—Computerised tomographs
- A61B6/037—Emission tomography
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
- A61B6/48—Diagnostic techniques
- A61B6/481—Diagnostic techniques involving the use of contrast agents
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
- A61B6/48—Diagnostic techniques
- A61B6/486—Diagnostic techniques involving generating temporal series of image data
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
- A61B6/50—Clinical applications
- A61B6/504—Clinical applications involving diagnosis of blood vessels, e.g. by angiography
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
- A61B6/50—Clinical applications
- A61B6/507—Clinical applications involving determination of haemodynamic parameters, e.g. perfusion CT
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B6/00—Apparatus for radiation diagnosis, e.g. combined with radiation therapy equipment
- A61B6/52—Devices using data or image processing specially adapted for radiation diagnosis
- A61B6/5258—Devices using data or image processing specially adapted for radiation diagnosis involving detection or reduction of artifacts or noise
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/22—Matching criteria, e.g. proximity measures
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2200/00—Indexing scheme for image data processing or generation, in general
- G06T2200/24—Indexing scheme for image data processing or generation, in general involving graphical user interfaces [GUIs]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10081—Computed x-ray tomography [CT]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10088—Magnetic resonance imaging [MRI]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10104—Positron emission tomography [PET]
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Medical Informatics (AREA)
- Physics & Mathematics (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- General Health & Medical Sciences (AREA)
- Radiology & Medical Imaging (AREA)
- Public Health (AREA)
- Pathology (AREA)
- Veterinary Medicine (AREA)
- Animal Behavior & Ethology (AREA)
- Surgery (AREA)
- Molecular Biology (AREA)
- Biophysics (AREA)
- Heart & Thoracic Surgery (AREA)
- Biomedical Technology (AREA)
- High Energy & Nuclear Physics (AREA)
- Optics & Photonics (AREA)
- Theoretical Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Dentistry (AREA)
- Oral & Maxillofacial Surgery (AREA)
- General Physics & Mathematics (AREA)
- Quality & Reliability (AREA)
- Vascular Medicine (AREA)
- Data Mining & Analysis (AREA)
- Pulmonology (AREA)
- Evolutionary Computation (AREA)
- General Engineering & Computer Science (AREA)
- Artificial Intelligence (AREA)
- Evolutionary Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Apparatus For Radiation Diagnosis (AREA)
- Image Processing (AREA)
Abstract
一种处理多个时间分离的图像的方法,包括:在各个图像中选择多个成像单元;测量在各个成像单元中的时域差;以及选择在阈限以上的时域差。
Description
相关申请的交叉引用
本申请要求于2008年10月2日提交的、名为“用于处理图像的系统和方法(systemandmethodforprocessingimages)”的、李(Lee)的美国临时申请No.61/102,206的权益,该美国临时申请的内容通过引用而并入本申请中。
技术领域
本发明总体上涉及图像处理,并且更具体地,涉及处理通过使用动态成像技术而获得的图像数据的系统和方法。
背景技术
医学成像涵盖用于为了临床目的创建人体的图像的技术和过程,包括用于诊断或监视疾病的医疗程序。医学成像技术已经发展为涵盖很多图像记录技术,所述图像记录技术包括电子显微镜法、荧光透视法、磁共振成像(MRI)、核医学、光声成像、正电子发射断层摄影术(PET)、投影射线摄影术、热成像法、计算机断层摄影术(CT)和超声。
医学成像可以包含,使用被称为造影剂或造影材料的复合物,以提高图像中体内结构的可见度。
医学成像原本局限于采集单个静态图像以捕获身体的器官/区域的解剖结构。目前,更复杂精密的成像技术的使用允许进行动态研究,所述动态研究提供可以表征生理或病理生理信息的图像的时域序列。
动态医学成像涉及获取过程,即,随时间的推移而拍摄受关注的器官/区域/身体的许多“快照(snapshot)”,以便捕捉随时间变化的行为,例如,造影剂的分布,且因而捕捉具体的生物状态(疾病、不适、生理现象等)。随着医学成像的速度和数字性的发展,此获取数据会具有极大的瞬时清晰度且会导致大量数据。
医学成像技术已广泛用于改善对诸如癌症、心脏疾病、大脑失调和心血管疾病的病情的诊断和护理。大部分评估表明,数以百万的生命由于这些医学成像技术而得以挽救或者显著得到改善。然而,必须对这类医学成像技术给患者带来的辐射暴露风险加以考虑。
在这点上,由于CT中所递送的辐射剂量比更为常见的常规x射线成像过程大,医学诊断中CT越来越多的使用已经对被暴露的患者增加的罹患癌症风险进行了突出的关注。对于两相CT灌注成像方案而言,则尤其如此。为了说明目的,将剂量-长度乘积用于由可市购的CT扫描仪(GEHealthcare)提供的方案,CT中风系列的有效剂量为10mSv,其中4.9mSv是由两相CT灌注成像方案贡献的,其中所述CT中风系列由以下组成:1)非增强CT扫描,以排除出血;2)CT血管造影术,定位导致中风的梗塞;和3)两相CT灌注成像方案,以用血流量和血容量来定义缺血区,并用血脑屏障渗透表面积乘积(BBB-PS)来预测出血性转变(HT)。预测CT中风系列对于每暴露的100,000个男性和女性急性缺血性中风(AIS)患者分别会引发80和131例额外癌症,并且仅两相CT灌注成像方案就被预测为导致了一半的额外癌症(HealthRisksfromExposuretoLowLevelsofIonizingRadiation:BEIRVII.TheNationalAcademiesPress,WashingtonDC,2006(暴露于低水平电离辐射的健康风险BEIRVII.国家学术出版社华盛顿特区2006))。
除非减少两相CT灌注成像方案的有效剂量,否则医学成像的益处将会由于对癌症引发的担忧而被破坏。此外,对辐射暴露风险的担忧延伸到其他医学成像技术,尤其是对于儿科患者或正经受重复暴露的那些患者。
例如,在2007年3月6日授权的第7,187,794号美国专利中,描述了使用统计滤波技术以增大低剂量辐射成像中的信噪比。然而,目前,统计滤波的应用限于图像构建之前的投影数据。若干不利之处与投影数据的统计滤波相关。例如,因为来自动态序列的各个图像一般是由~1,000个投影重建而成,所以使用投影数据工作时计算负担高。另外,滤波过程必须被“硬布线”至扫描仪的图像重建管线中。如此一来,投影数据的统计滤波器由于其通常专用于扫描仪平台而缺乏灵活性。
不论是否已被觉察或证实,对辐射暴露风险的担忧将必须得到解决,以使患者对医学成像技术的长期安全放心。据此,需要允许辐射剂量的减少的医学成像技术。
因此,本发明的目的是要提供一种新颖的用于处理图像的系统和方法。
发明内容
在一方面,提供有一种处理多个时间分离的图像的方法,包括:在各个图像中选择多个成像单元;为各个成像单元确定时域差;以及选择在阈限以上的时域差。
在另一方面,提供有一种处理多个时间分离的图像的系统,包括:接口,用于接收多个时间分离的图像;和处理器,用于在各个图像中选择多个成像单元,为各个成像单元确定时域差,以及选择在阈限以上的时域差。
在又一方面,提供有一种包含处理多个时间分离的图像的计算机程序的计算机可读介质,所述计算机程序包括:用于在各个图像中选择多个成像单元的计算机程序代码;用于为各个成像单元确定时域差的计算机程序代码;和用于选择在阈限以上的时域差的计算机程序代码。
附图说明
现在,将仅通过示例的方式、参照附图来描述实施例,在附图中:
图1是图像处理方法的流程图;
图2示意性图示了多个时间分离的图像中的成像单元的定位;
图3A至3E示出了用统计滤波技术处理过的、来自CT图像的脑血流量图;
图4图示性示出了使用图3E中勾画出的区域为图3A至3D中所示的图像所计算的血流量品质因数(FigureofMerit);以及
图5是用于实施图1中所示图像处理方法的系统的图示。
具体实施方式
动态医学成像涉及采集过程,即,随时间的推移而拍摄受关注的器官/区域/身体(即靶区)的许多“快照”,以便捕捉随时间变化的行为,例如造影剂的吸取和/或消除。此采集过程导致多个时间分离的图像的产生。本申请中所描述的方法和系统涉及这样的时间分离的图像的处理。
现在转到图1,图示出了根据该方法在时间分离的图像的处理期间所执行的步骤。多个时间分离的图像可以在使用或者不使用注射的造影剂的情况下从患者的动态医学成像获得。造影剂可以例如根据靶区的结构或生理状态,来选择性地增加图像中靶区的对比度。例如,可以将对特定器官、疾病、状态或生理过程具有生物物理、分子、基因或细胞亲和性的复合物注射到患者体内。选择这样的造影剂,使之具有下述属性:通过改变成像条件,例如,通过改变图像对比度,对给定成像技术提供增强的信息,以反映体内复合物的行为。这可以通过局域性部位处增大的X射线衰减(例如对于CT/X射线)、被改变的顺磁特性(例如对于MRI)或者对于核医疗/PET使用放射性同位素来实现。用于许多成像技术的造影剂是公知的。在一些情况下,对比增强并不依赖于所注入的造影剂,所述一些情况例如:在磁共振成像(MRI)中使用“黑血”或“白血”,其中,使用特定脉冲序列以改变血的磁饱和且因而改变其在图像中的呈现;或者使用带标记的MRI序列,所述带标记的MRI序列改变具体组织或流体的磁性能。所注入的造影剂或者具有改变了的性能的组织或流体可以全部被认为是“成像剂”。
多个时间分离的图像不限于任何具体的成像技术,并且包括例如使用磁共振成像(MRI)、计算机断层摄影术(CT)、核医学(NM)或正电子发射断层摄影术(PET)的动态医学成像。多个时间分离的图像也不限于具体图像类型。例如,可以处理灰度或彩色图像。
在该方法期间,在各个时间分离的图像内,选择多个成像单元(步骤110)。“成像单元”可以为用于将图像分成等同的部分的任何所期望的单元,包括例如像素、多个像素、部分像素、体素、多个体素、部分体素等。
成像单元数据由例如数字值的值来表示,并且因而可以量化和测量。一般而言,为时间分离的图像中的对应的成像单元测量图像强度。对于使用造影剂的动态成像而言,可以关于造影剂浓度确定时域差。造影剂浓度的时域差可以被表示为作为时间的函数的对比增强的图示。
图1中示出了用于确定造影剂浓度的时域差并选择时域差的两种备选方案。在图1所示的一个备选方案中,一旦在步骤110从各个时间分离的图像中选择了多个成像单元,就为各个成像单元确定时域差(步骤115)。然后根据所观察到的时域差的程度来对所确定的时域差或其表示进行排序(步骤120)。然后选择在阈限以上的时域差(步骤125)。然后基于所选的时域差来构建处理后的图像(步骤150)。
在图1所示的另一备选方案中,不对时域差排序。而是,利用统计技术,且更具体地,利用盲信号分离统计技术,来分析成像单元的时域数据,以选择在阈限以上的时域差。例如,关于使用造影剂的动态CT成像,一旦在步骤110已经从各个时间分离的图像中选出了多个成像单元,就建立与平均曲线相关的对比增强图示的所有对的协变的协方差矩阵(步骤130)并使用盲信号分离统计技术进行分析(步骤135),诸如对协方差矩阵应用主分量分析。基于统计分析来计算特征向量/特征值对(步骤140)。然后所得到的在阈限以上的特征向量/特征值对基于它们的特征值被选择出来(步骤145)。然后利用所选的特征向量来构建处理后的图像,并且使用所选特征向量的不同权重来构建各成像单元(步骤160)。
可选地,可以配准多个时间分离的图像,以使时间分离的图像中的对应位置彼此相关。仅需要多个时间分离的图像的配准来校正显著失配。否则,在时间分离的图像之间的轻微移动可以为本申请中所描述的处理方法所容许甚至可以被消除。
作为另一可选步骤,可以对多个时间分离的图像进行预处理,使得时间分离的图像由可量化的值或代码表示。例如,可以利用已知数字化技术而由数字值表示时间分离的图像。如果多个时间分离的图像的成像数据是以数字形式提供的,则无需数字化步骤。
现在将参照图2来进一步描述造影剂浓度的时域差的确定。图2示出了一系列6个按照时间顺序排列的时间分离的图像t1至t6。在各个时间分离的图像上标记有成像单元的4×6阵列。在此示例中,选择了五(5)个成像单元IU1至IU5。如从比较整个系列的时间分离的图像t1至t6清楚可知,对于各个时间分离的图像在阵列的相同位置(即第1列第2行)处选择成像单元IU1。相似地,贯穿该系列时间分离的图像在成像单元阵列的对应位置处选择成像单元IU2至IU5中的各个。据此,可以将贯穿时间分离的图像在其对应位置处的各个成像单元的值绘制为时间的函数,以获得时域曲线。绘制步骤由连接成像单元IU1、IU2和IU5的对应位置的箭头图示出。然后可以对各个成像单元的时域曲线进行分析,以确定时域差。例如,可以计算所获得的全部时域曲线的平均曲线。然后从各个时域曲线减去平均曲线,以获得各个所选成像单元的时域差曲线。将予以理解的是,所标记的阵列和成像单元IU1至IU5的选择、以及表示将成像单元IU1、IU2和IU5的值作为时间的函数进行绘图的箭头,仅用于说明目的,并且可以按照需要选择任何数目的成像单元并对其进行分析。此外,通过从成像单元的时域曲线减去平均曲线而进行的时域差的确定仅用于说明目的,并且可以想到分析时域曲线的任何其他方法,例如使用导数和/或统计技术。
现在,将使用来自CT灌注研究的时间分离的图像作为示例,提供用于多个时间分离的图像的盲信号分离统计分析的数学基础。
给定来自CT灌注研究的一系列n个动态(即时间分离的)图像,n个动态图像可以紧凑地表示为矩阵:
其中,p为CT图像中的像素数(即p=512×512、256×256等),且为n×1向量,其元素为n个图像中的像素值或来自第i个像素的时间-增强曲线。注意,的第(i,j)个元素为的第j个元素xij。各个也可以解释为单个时间-增强曲线(由n个随机变量χ1,χ2,…,χn-1,χn组成的随机过程)的重复观察。
符号解释如下:
(i)向量由具有排成一列的元素的希腊或罗马字母字符上方的颚化符号(~,tilde)来表示;
(ii)矩阵由大写的希腊或罗马字母字符上方的颚化符号来表示;
(iii)向量或矩阵的转置由符号上的上标T来表示;以及
(iv)随机变量由希腊字母字符表示,而随机变量的观察值(样本)由对应的罗马字母字符表示;希腊和罗马字母字符之间的对应关系将在下列讨论中如下所示使用:
(v)由n个随机变量组成的随机过程由n×1向量表示。
从n个动态图像中消除噪声的一个途径是,找到p×n矩阵的“展开式(expansion)”,并且然后按照某个预定标准截短展开式。常用的展开式为的奇异值分解:
其中,为p×p正交矩阵,为n×n正交矩阵,且为p×n对角矩阵并将的奇异值作为其对角元素。令l1,…,ln为的奇异值(在不失一般性的情况下假定具有满列秩n,且p>n)。则:
其中,为p×1向量且其元素除第i个元素为li以外均为零。通常,奇异值l1,…,ln按降序排列,并且小奇异值由于其是由噪声引起的因而可以被丢弃。减少中的噪声的简单方式是只保留r<n个奇异值。则:
如果截短后的奇异值序列被用以根据等式(1)来重构则:
其中,为的第i列,且为的第i列。估算展开式,
得到:
由此,在的展开式中,保留r个奇异值意味着只需要和的前r列。
然而,如这里用公式表达的的奇异值分解(SVD)(展开式)在对来自CT灌注的动态图像去噪时存在两个主要缺点,即:
1.矩阵的大小可以大至(256)2×40-60,从而使得在计算上难以找到的SVD;并且
2.虽然丢弃小奇异值的标准具有一般意义,但难以验证用以保留或丢弃奇异值的阈值。
现在将描述用以获得的SVD的替代方法。如上所述,来自各个像素(或像素块)的时间-增强曲线可以视为底层(underlying)随机过程(时间-增强曲线)的重复样本(观察值)。则随机过程的样本均值可以计算为:
按照定义为n×1向量(样本均值的时间-增强曲线),且根据等式(3),的第j个分量为:
其中为的第j个分量
来自的零均值随机过程为其中为的总体均值。对于的各个观察值(样本) 的对应的观察值(样本)由下式给出:
的样本均值(一个n×1向量)由下式给出:
其中,由于为零均值随机过程,所以如所预期,为n×1零向量。
令为p×n矩阵,其各行为 则:
其中,为ψii=1,…,n的样本方差,且为ψi和ψji,j=1,…,n的样本协方差。
或者,
为的样本协方差矩阵并且具有n×n的维数。在CT灌注中n一般为40-60。
ψi,i=1,…,n,的线性组合ζ可以定义为:
其中,为给定的n×1系数向量。ζ为随机变量(非过程),且ζ的观察值zii=1,…,p按照ζ的定义而由下式给出:
ζ的样本均值由下式给出:
由此,ζ也为零均值随机变量。ζ的样本方差由下式给出:
现在将就确定系数向量以最大化ζ的样本方差进行研究。由于可以通过按比例将放大常数倍而使与所需要的一样大,所以可以只通过施加诸如的归一化条件来进行优化。为了使满足的最大化,标准途径是使用拉格朗日乘子技术。即寻求使下式最大化:
其中,λ为拉格朗日乘子。对求微分得出:
由此,λ为的特征值,并且为对应的特征向量,且
的特征值也就是使满足的最大化,为的具有最高特征值的特征向量。
令λ1为对应特征向量为的最大特征值,则被称为的第一样本主分量。其为具有由定义的观察值的随机变量。zi1也被称为第一主分量上的第i个观察值的得分值。类似地,第二主分量在满足与ζ1不相关和归一化条件的情况下,使得ζ2的样本方差最大化。对应于ζ2的观察值为ζ1和ζ2的样本协方差为:
ζ2与ζ1不相关的要求等同于要求:
这样,下列任何条件可以被用以规定ζ2与ζ1不相关:
使满足ζ2与ζ1不相关且的条件的ζ2的样本方差最大化等同于最大化:
其中,条件被用以规定ζ2与ζ1不相关并且λ2和为拉格朗日乘子。对求微分,得到:
通过乘以得出:
由于前面两项为零且所以其简化为
因此,
所以λ2再一次是的特征值且是对应的特征向量。
再者,或者ζ2的样本方差为λ2。假定没有重复的特征值,则λ2不能等于λ1。如果等于,则会有从而违背了的约束条件。因而,λ2为第二大特征值,且为对应的特征向量。
可以看出,对第三个、第四个、……、第n个样本主分量使用相同技术,系数向量为对应于第三大、第四大、……、以及最小特征值λ3,λ4,…,λn的的特征向量。此外,主分量的样本方差也由的特征值给出。
将定义为主分量上的观察值的得分值的p×n矩阵。在此情况下,第i行由所有主分量ζ1,ζ2,…,ζn-1,ζn上的的第i个观察值(即)的得分值构成,即:
其中,为n×n矩阵,其各列为的特征向量并且正交。
最大化后的第一主分量ζ1由下式表示:
且为矩阵的第i行,而来自全部像素(像素块)的TEC中的特征向量的载荷为:
且为矩阵的第i列。因此,的第i列给出了对应于第i大特征值λi的特征向量的载荷图(loadingmap)。
为了对来自CT灌注研究的动态图像去噪,需要对特征向量的载荷图的计算求逆。需要由截短后的一系列特征向量重构原始动态图像的平滑版本
如等式(1)中所示,p×n矩阵可以通过奇异值分解而展开为:
其中,为p×p正交矩阵,为n×n正交矩阵,且为p×n对角矩阵并将的奇异值作为其对角线元素。令l1,…,ln为的奇异值(在不失一般性的情况下假定具有满列秩n,且p>n)。使用奇异值分解:
将等式两边乘以得到:
其中,为p×1向量,其元素除等于li的第i个元素以外均为零,即为且
则,
其中
令其中为n×1向量。
则,
代入等式(5)中得到:
根据对特征向量和特征值的讨论有:
或者的奇异值等于的特征值乘以(p-1)的平方根,其中
由此, 或
接下来,考虑
(注意,不是的转置矩阵)
将等式两边乘以得到:
因而,
令为如所定义的p×1向量,则
令其中为p×1向量。
则,
代入等式(6)中得到:
由此,为对应于第i个特征值的的特征向量,也碰巧为的第i个特征值,同时全部为具有零特征值的的特征向量。并且
如上所述,等式(5a)为:
将两边乘以得到:
与等式(6a)对比:
其中c为常数(7a)
类似地,将等式(6a)两边乘以得到;
与等式(5a)对比:
其中c为常数(7b)
由等式(7b):
代入等式(7a)中得到:
代入等式(6a)中得到:
因此,
在不失一般性的情况下,将c*设为:
由此,对于的奇异值分解:
其中为p×p正交矩阵:
且为p×1向量;为n×n正交矩阵:
且为n×1向量;
且为p×n对角矩阵,且将的奇异值l1,…,ln作为其对角线元素;在不失一般性的情况下还假定具有满列秩n,且p>n。
示出有:
其中为对应于特征值的的第i个特征向量;并且由等式(7c)有
如上所述,为了对除噪,可以在第r个奇异值后面截短SVD展开式:
等式(8)也提供了如下列所解释的非常简单的几何解释。
其中,为n×r矩阵,其第i列为且由等式(4a)得出:
并且的行和列具有与的行和列相同的解释。
改写等式(9):
等式(10)表明,第i个像素(像素块)的TEC为r个特征向量的加权和,且该加权的权重由第i个像素TEC在特征向量上的载荷来给定。
特征向量可以被认为在n维空间中形成正交基。在上的投影为且可以重组为
用以对去噪的算法如下:
1.计算对应于的p×n零均值矩阵,其中
2.计算n×n样本协方差矩阵
3.计算n对特征值和特征向量,
4.通过只保留最大的r个特征值来计算n×n平滑矩阵:
5.通过计算以下乘积来平滑
6.如下计算的平滑版本:
现在将就实验结果来说明本申请中所描述的方法的示例。针对实现减少CT扫描剂量的能力来对本申请中所描述的主分量分析(PCA)进行测试。
现在将基于CT灌注研究来概述PCA方法,所述CT灌注研究产生捕获通过大脑的造影剂的通路的多个时间分离的图像。对应于该图像序列中的各个像素,生成将该像素中的造影剂浓度表示为时间的函数的曲线。还生成等于该曲线与全部像素曲线的均值之差的差曲线(differencecurve)。计算所有差曲线对的协变的协方差矩阵并且其表示与平均曲线有关的全部可能偏差。主分量分析法分析协方差矩阵,以通过找出表示关于平均曲线的大部分偏差的差曲线的线性组合,而找出这些曲线之中的共有时域图案(pattern)。可以看出,这样的线性组合或者主分量,是协方差矩阵的特征向量,并且具有最大特征值的特征向量是全部差曲线之中占最大主导地位的时域特征,具有第二大特征值的特征向量为占第二大主导地位的时域特征,依此类推,则具有小特征值的特征向量来自于CT图像中的噪声。推论是,各个差曲线可以表示为全部主分量之和,并且通过丢弃具有小特征值的分量能够有效抑制噪声。为了对图像的时间序列去噪,将平均曲线加至通过使用具有大特征值的前几个主分量而重组的全部差曲线。
尽管目前可用的统计滤波器在图像重建之前消除了投影数据中的噪声,但主分量分析法消除了来自重建后的图像的噪声。下面描述的实验结果显示出执行PCA致使减少CT成像的剂量。
就CT灌注应用的辐射剂量减少,相对于传统的滤波反投影(FBP)法来测试PCA法。在健康的猪身上对PCA对比FBP在使用CT灌注来测量脑灌注上的剂量节约进行了测试,其中,对该健康的猪扫描四次,该四次扫描分别用从190mA改变到100mA、75mA和50mA的x射线管电流而保持其他扫描参数相同。
图3示出了用下述参数和方法获得的CT灌注图像且利用J-W模型而计算的猪的脑血流量图:(A)190mA且FBP重建;(B)100mA,FBP重建且PCA;(C)70mA,FBP重建且PCA;以及(D)50mA,FBP重建且PCA。图3(E)示出了用于计算图4中所示血流量图的品质因数的受关注区域(3个勾勒出的圈)。这些区域重叠在示出了通过猪的大脑的冠状断面的CT图像上。暗区域为颅骨和颌骨。
图4示出了图3(E)中所示三个区域中的血流量的标准差/均值(品质因数):70mA且PCA比190mA且FBP好,而50mAPCA比190mA且FBP差(p<0.05)。图3中血流量图的图像质量也反映了此发展(progression)。因此,PCA法能够减少用以在CT灌注研究中获得图像并且相比于仅使用FBP而言维持所得到的血流量图的质量所需要的辐射剂量。更具体来说,与仅使用FBP相比,PCA法使辐射剂量减少至将近三分之一。因为血脑屏障保持完整无缺,所以没有展示出PCA法在生成BBB-PS图时对剂量减少的效果。然而,由于血流量、血容量和BBB-PS是在J-W模型中一起计算的,所以预期利用血流量图示出的剂量减少适用于BBB-PS图。
本申请中所描述的主分量分析法已利用C++而实现于计算机可执行软件应用程序中。在诸如例如个人计算机(PC)的通用计算设备上用以处理九十五(95)个512×512的CT脑图像所需的时间少于两(2)分钟。
本申请中所描述的方法可以在诸如图5中所示的系统400的硬件上加以实现。输入510接收多个时间分离的图像,该多个时间分离的图像可以是之前存储的,从成像设备直接接收的,或者从远程位置经由有线或无线通信链接传送的。通用计算设备505接收来自输入的时间分离的图像并执行软件应用程序,从而如上所述处理时间分离的图像。通用计算设备505通过显示器520、键盘515和/或鼠标或者其他指向设备525而与用户进行接口交互。结果,例如所构建的处理后的图像,可以呈现在显示器510和/或输出至任何适当的输出设备530,例如,打印机、存储设备或者用于将结果传送给远程位置的通信设备。
本申请中所描述的软件应用程序可以作为单机应用程序运行,或者可以包含在其他可用应用程序中,以向那些应用软件提供增强的功能。所述软件应用程序可以包括包含例程、程序、对象构件、数据结构等的程序模块,并且可以实现为存储在计算机可读介质上的计算机可读程序代码。所述计算机可读介质是能够存储数据的任何数据存储设备,而所述数据之后可由计算机系统读取。计算机可读介质的示例包括,例如,只读存储器、随机存取存储器、CD-ROM、磁带和光学数据存储设备。所述计算机可读程序代码也可以分布在包含耦接的计算机系统的网络上,从而计算机可读程序代码可以分布式方式存储和执行。
上述实施例旨在作为示例,并且本领域技术人员可以在不脱离由所附的权利要求书所限定的本发明的范围的情况下对其进行改变和修改。
Claims (25)
1.一种医学图像处理方法,包括:
捕获受关注的靶区的多个时间分离的图像,从而得到所述受关注的靶区的随时间变化的行为;
将各个所述时间分离的图像分成成像单元,所述成像单元具有相等的大小;
从各个所捕获的图像选择多个成像单元,所选择的成像单元处在各个所述所捕获的图像中的相同位置;
测量各个所选择的成像单元的强度,基于所测量的强度由值来表示各个所选择的成像单元;
对于各个成像单元位置,使用表示在各个所述时间分隔的图像中所选择的成像单元的所述值,为所述选择的成像单元确定时域差;
选择在阈限以上的时域差;以及
基于所选择的时域差来构建噪声减少的图像。
2.如权利要求1所述的方法,进一步包括:在所述选择之前对各个测得的时域差排序。
3.如权利要求1所述的方法,其中,各个成像单元为一个像素和一个体素中的一种。
4.如权利要求1所述的方法,其中,各个成像单元为多个像素和多个体素中的一种。
5.如权利要求1所述的方法,其中,所述时域差基于对比增强。
6.如权利要求5所述的方法,其中,所述时域差为造影剂的差别。
7.如权利要求6所述的方法,其中,对比增强基于所述造影剂的浓度。
8.如权利要求6所述的方法,其中,对比增强基于所述造影剂的激活。
9.如权利要求1所述的方法,进一步包括:利用盲信号分离统计技术来表征所测得的时域差。
10.如权利要求1所述的方法,其中,确定所述时域差包括计算所述成像单元的时域差的协方差矩阵。
11.如权利要求10所述的方法,进一步包括:利用盲信号分离统计技术分析所述协方差矩阵,以确定所述时域差之中的共有时域图案。
12.如权利要求11所述的方法,进一步包括:分析所述协方差矩阵的特征向量。
13.如权利要求12所述的方法,其中,所述特征向量的特征值表示时域差的程度。
14.如权利要求13所述的方法,其中,所述选择时域差包括选择在所述阈限以上的特征值。
15.如权利要求1至14中任一项所述的方法,其中,所述阈限是能够调节的。
16.如权利要求11所述的方法,其中,所述盲信号分离统计技术选自由下述技术构成的组:主分量分析、独立分量分析、盲源解卷积、卡亨南-拉维变换和霍特林变换。
17.如权利要求1至14中任一项所述的方法,其中,所述多个时间分离的图像是使用动态医学成像技术而获得的。
18.如权利要求17所述的方法,其中,所述动态医学成像技术包括造影剂的使用。
19.如权利要求17所述的方法,其中,所述时间分离的图像是使用具有减少的辐射剂量的方案而获得的。
20.如权利要求17所述的方法,其中,所述时间分离的图像是使用具有减少的造影剂量的方案而获得的。
21.如权利要求17所述的方法,其中,所述动态医学成像技术为磁共振成像(MRI)、功能性MRI、正电子发射断层摄影术、核医学或计算机断层摄影术。
22.如权利要求17所述的方法,其中,所述动态医学成像技术为计算机断层摄影术(CT)。
23.如权利要求22所述的方法,其中,所述时间分离的图像是使用具有减少的辐射剂量的CT灌注方案而获得的。
24.如权利要求22所述的方法,其中,所述时间分离的图像是使用具有减少的造影剂量的CT灌注方案而获得的。
25.一种处理多个时间分离的图像的系统,包括:
接口,用于接收所捕获的受关注的靶区的多个时间分离的图像,所述多个时间分离的图像示出了所述受关注的靶区的随时间变化的行为,通过低剂量辐射成像捕获所述时间分离的图像;和
处理结构,用于与所述接口通信并且处理所述捕获的时间分隔的图像,在所述捕获的时间分隔的图像的处理期间,所述处理结构将各个所述时间分离的图像分成成像单元,所述成像单元具有相等的大小;从各个所捕获的图像选择多个成像单元,所选择的成像单元处在各个所述所捕获的图像中的相同位置;测量各个所选择的成像单元的强度,基于所测量的强度由值来表示各个所选择的成像单元;对于各个成像单元位置,使用表示在各个所述时间分隔的图像中所选择的成像单元的所述值,为所述选择的成像单元确定时域差;选择在阈限以上的时域差;以及基于所选择的时域差来构建噪声减少的图像。
Applications Claiming Priority (3)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
US10220608P | 2008-10-02 | 2008-10-02 | |
US61/102,206 | 2008-10-02 | ||
PCT/CA2009/001397 WO2010037233A1 (en) | 2008-10-02 | 2009-10-02 | System and method for processing images |
Publications (2)
Publication Number | Publication Date |
---|---|
CN102282587A CN102282587A (zh) | 2011-12-14 |
CN102282587B true CN102282587B (zh) | 2015-11-25 |
Family
ID=42073000
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN200980147491.XA Expired - Fee Related CN102282587B (zh) | 2008-10-02 | 2009-10-02 | 用于处理图像的系统和方法 |
Country Status (4)
Country | Link |
---|---|
US (4) | US20110262022A1 (zh) |
EP (1) | EP2332122A4 (zh) |
CN (1) | CN102282587B (zh) |
WO (1) | WO2010037233A1 (zh) |
Families Citing this family (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
EP2583241A1 (en) * | 2010-06-21 | 2013-04-24 | Koninklijke Philips Electronics N.V. | Method and system for noise reduction in low dose computed tomography |
EP2656241A4 (en) * | 2010-12-24 | 2017-06-07 | FEI Company | Reconstruction of dynamic multi-dimensional image data |
WO2015058320A1 (en) * | 2013-10-25 | 2015-04-30 | Acamar Corporation | Denoising raw image data using content adaptive orthonormal transformation with cycle spinning |
EP3179919B1 (en) * | 2014-08-14 | 2021-02-17 | Koninklijke Philips N.V. | Acoustic streaming for fluid pool detection and identification |
CN105809670B (zh) | 2016-02-29 | 2019-07-19 | 上海联影医疗科技有限公司 | 灌注分析方法 |
US20190150764A1 (en) * | 2016-05-02 | 2019-05-23 | The Regents Of The University Of California | System and Method for Estimating Perfusion Parameters Using Medical Imaging |
JP7086630B2 (ja) * | 2018-02-09 | 2022-06-20 | キヤノン株式会社 | 情報処理装置、情報処理方法、及びプログラム |
CN112035818B (zh) * | 2020-09-23 | 2023-08-18 | 南京航空航天大学 | 一种基于物理加密辐射成像身份认证系统 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6045775A (en) * | 1993-07-12 | 2000-04-04 | Nycomed Imaging As | Parenteral administration of positive and negative contrast agents for MRI |
WO2005120353A1 (en) * | 2004-06-14 | 2005-12-22 | Canon Kabushiki Kaisha | Image processing device and method |
Family Cites Families (17)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
AU3418700A (en) * | 1999-03-26 | 2000-10-16 | Leif Ostergaard | Method for determining haemodynamic indices by use of tomographic data |
US7187794B2 (en) | 2001-10-18 | 2007-03-06 | Research Foundation Of State University Of New York | Noise treatment of low-dose computed tomography projections and images |
US7078897B2 (en) * | 2002-01-16 | 2006-07-18 | Washington University | Magnetic resonance method and system for quantification of anisotropic diffusion |
US7206435B2 (en) * | 2002-03-26 | 2007-04-17 | Honda Giken Kogyo Kabushiki Kaisha | Real-time eye detection and tracking under various light conditions |
US20040101156A1 (en) * | 2002-11-22 | 2004-05-27 | Dhiraj Kacker | Image ranking for imaging products and services |
FR2848093B1 (fr) * | 2002-12-06 | 2005-12-30 | Ge Med Sys Global Tech Co Llc | Procede de detection du cycle cardiaque a partir d'angiogramme de vaisseaux coronaires |
GB0300921D0 (en) * | 2003-01-15 | 2003-02-12 | Mirada Solutions Ltd | Improvements in or relating to dynamic medical imaging |
US20040218794A1 (en) | 2003-05-01 | 2004-11-04 | Yi-Hsuan Kao | Method for processing perfusion images |
US20050113680A1 (en) * | 2003-10-29 | 2005-05-26 | Yoshihiro Ikeda | Cerebral ischemia diagnosis assisting apparatus, X-ray computer tomography apparatus, and apparatus for aiding diagnosis and treatment of acute cerebral infarct |
GB0326381D0 (en) * | 2003-11-12 | 2003-12-17 | Inst Of Cancer Res The | A method and means for image processing |
US7668359B2 (en) * | 2004-06-28 | 2010-02-23 | Siemens Meidcal Solutions USA, Inc. | Automatic detection of regions (such as, e.g., renal regions, including, e.g., kidney regions) in dynamic imaging studies |
EP1647937A1 (en) * | 2004-10-15 | 2006-04-19 | Sony Deutschland GmbH | Method for motion estimation |
US7826878B2 (en) * | 2004-12-07 | 2010-11-02 | Research Foundation Of City University Of New York | Optical tomography using independent component analysis for detection and localization of targets in turbid media |
JP4437071B2 (ja) | 2004-12-24 | 2010-03-24 | パナソニック株式会社 | 運転支援装置 |
WO2006114003A1 (en) * | 2005-04-27 | 2006-11-02 | The Governors Of The University Of Alberta | A method and system for automatic detection and segmentation of tumors and associated edema (swelling) in magnetic resonance (mri) images |
US7620204B2 (en) * | 2006-02-09 | 2009-11-17 | Mitsubishi Electric Research Laboratories, Inc. | Method for tracking objects in videos using covariance matrices |
US20080100473A1 (en) * | 2006-10-25 | 2008-05-01 | Siemens Corporate Research, Inc. | Spatial-temporal Image Analysis in Vehicle Detection Systems |
-
2009
- 2009-10-02 US US13/122,379 patent/US20110262022A1/en not_active Abandoned
- 2009-10-02 WO PCT/CA2009/001397 patent/WO2010037233A1/en active Application Filing
- 2009-10-02 CN CN200980147491.XA patent/CN102282587B/zh not_active Expired - Fee Related
- 2009-10-02 EP EP09817153.1A patent/EP2332122A4/en not_active Withdrawn
-
2013
- 2013-09-27 US US14/039,619 patent/US8965086B2/en not_active Expired - Fee Related
-
2015
- 2015-02-23 US US14/628,731 patent/US20150230771A1/en not_active Abandoned
-
2016
- 2016-10-28 US US15/338,000 patent/US20170109886A1/en not_active Abandoned
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6045775A (en) * | 1993-07-12 | 2000-04-04 | Nycomed Imaging As | Parenteral administration of positive and negative contrast agents for MRI |
WO2005120353A1 (en) * | 2004-06-14 | 2005-12-22 | Canon Kabushiki Kaisha | Image processing device and method |
Non-Patent Citations (2)
Title |
---|
一种视频序列的拼接算法;郭三华 等;《计算机应用》;20071130;第27卷(第11期);2786-2788 * |
应用主成分分解(PCA)法的图像融合技术;王文武;《微计算机信息》;20071231;285-287 * |
Also Published As
Publication number | Publication date |
---|---|
CN102282587A (zh) | 2011-12-14 |
US20140044334A1 (en) | 2014-02-13 |
EP2332122A1 (en) | 2011-06-15 |
EP2332122A4 (en) | 2013-11-20 |
US8965086B2 (en) | 2015-02-24 |
US20170109886A1 (en) | 2017-04-20 |
US20150230771A1 (en) | 2015-08-20 |
WO2010037233A1 (en) | 2010-04-08 |
US20110262022A1 (en) | 2011-10-27 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN102282587B (zh) | 用于处理图像的系统和方法 | |
Graff et al. | Compressive sensing in medical imaging | |
Yu et al. | Prediction of human observer performance in a 2‐alternative forced choice low‐contrast detection task using channelized Hotelling observer: Impact of radiation dose and reconstruction algorithms | |
Britten et al. | The addition of computer simulated noise to investigate radiation dose and image quality in images with spatial correlation of statistical noise: an example application to X-ray CT of the brain | |
US9924887B2 (en) | Medical image display apparatus | |
CN110462689A (zh) | 基于深度学习的断层摄影重建 | |
CN104077791B (zh) | 一种多幅动态对比度增强核磁共振图像联合重建方法 | |
CN108961237A (zh) | 一种基于卷积神经网络的低剂量ct图像分解方法 | |
Solomon et al. | Correlation between human detection accuracy and observer model-based image quality metrics in computed tomography | |
Huang et al. | Parametric image generation with the uEXPLORER total-body PET/CT system through deep learning | |
Lan et al. | SC-GAN: 3D self-attention conditional GAN with spectral normalization for multi-modal neuroimaging synthesis | |
WO2021062413A1 (en) | System and method for deep learning for inverse problems without training data | |
US11514621B2 (en) | Low-dose image reconstruction method and system based on prior anatomical structure difference | |
CN107146263B (zh) | 一种基于张量字典约束的动态pet图像重建方法 | |
WO2023134030A1 (zh) | 一种基于流模型的pet系统衰减校正方法 | |
Smith et al. | Variability in image quality and radiation dose within and across 97 medical facilities | |
Marin et al. | Numerical surrogates for human observers in myocardial motion evaluation from SPECT images | |
Tustison et al. | ANTsX: A dynamic ecosystem for quantitative biological and medical imaging | |
Clarkson et al. | A task-based approach to adaptive and multimodality imaging | |
Ba et al. | Anthropomorphic model observer performance in three-dimensional detection task for low-contrast computed tomography | |
Clark et al. | Deep learning based spectral extrapolation for dual‐source, dual‐energy x‐ray computed tomography | |
Zhang et al. | Deep learning image reconstruction in pediatric abdominal and chest computed tomography: A comparison of image quality and radiation dose | |
Zhou et al. | Residual-based convolutional-neural-network (CNN) for low-dose CT denoising: impact of multi-slice input | |
Wunderlich et al. | Confidence intervals for performance assessment of linear observers | |
Lian et al. | Spatiotemporal attention constrained deep learning framework for dual-tracer PET imaging |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20151125 Termination date: 20171002 |