CN102201038A - 脑瘤p53蛋白表达检测方法 - Google Patents

Abstract

本发明公开一种基于磁共振图像分析的脑瘤P53蛋白表达检测方法，该方法包括：（1）采集脑瘤患者的包括T1加权序列、T1增强序列、FLAIR序列的磁共振图像；（2）在前述任一种序列中截取磁共振图像中的病变区域图像，由病变区域图像构成病变区域图像集，并将病变区域图像标记为P53蛋白表达阳性或阴性；（3）对病变区域图像集进行图像分析，提取病变区域图像中的图像特征，由图像特征构成病变区域样本集；（4）从病变区域样本集中任选一部分样本作为训练样本集，其他样本作为验证样本集，使用训练样本集训练分类器；（5）使用训练后的分类器对验证样本进行分类，获得验证样本的脑瘤P53蛋白表达状况。

Description

脑瘤P53蛋白表达检测方法

技术领域

本发明属于图像处理和模式识别技术领域，具体涉及一种通过磁共振图像分析进行脑瘤P53蛋白表达状况检测的方法。

背景技术

神经胶质瘤是中枢神经系统最常见的肿瘤，具有治疗周期长、易复发、致残率和死亡率高的特点，对患者的健康和生命威胁极大。目前，神经胶质瘤可以做到早期诊断，给予及时的手术、放疗和化疗，但是恶性胶质瘤患者的预后并没有明显改善。Afshar等(GolnarAfshar，Nannette Jelluma，Yang Xiaodong et al.Radiation-Induced Caspase-8Mediates p53-IndependentApoptosis in Glioma Cells.Cancer Research.2006，66(8)：4223-4232)研究表明，成年胶质瘤患者的P53基因突变与预后的改善和放疗敏感性有关。临床通常检测phosphatase and tensinhomolog deleted on chromosome 10(PTEN)、epidermal growth factor receptor(EGFR)、theO⁶-methylguanine-DNA methyltransferase(MGMT)、tumor protein 53(P53)、monoclonal antibodyofcell proliferation associated nuclear antigen(MIB-1)等蛋白表达状况，并综合上述蛋白表达的信息来评价胶质瘤的恶性程度和患者的预后状况。因此，P53的表达状况作为一种中间结果信息，虽然不能直接评价胶质瘤的恶性程度和患者预后，但具有一定的指导意义。

当前临床使用最为广泛的P53蛋白的检测方法是免疫组织化学技术，该技术需要手术获得患者的胶质瘤病理切片后才能进行检测，因此无法指导制定术前治疗方案。同时，周晓军(周晓军.免疫组化在病理诊断中的正确应用.诊断病理学杂志.2003，10(4)：232-235)研究表明，免疫组织化学技术在标准化和结果量化方面存在不足，检测结果容易受到检测人员的主观影响。目前尚没有基于图像处理和模式识别技术进行P53蛋白检测的方法。

发明内容

本发明的目的是提供一种基于图像处理和模式识别技术的脑瘤P53蛋白表达检测方法。

为实现上述目的，本发明所采取的技术方案是：该脑瘤P53蛋白表达检测方法包括以下步骤：

(1)采集脑瘤患者的磁共振图像，所述磁共振图像包括T1加权序列、T1增强序列、FLAIR序列中的任一种或任几种；

(2)在步骤(1)的任一种所述序列中截取所述磁共振图像中的病变区域图像，由所述病变区域图像构成病变区域图像集，并将病变区域图像标记为P53蛋白表达阳性或阴性；

(3)通过对所述病变区域图像集进行图像分析，提取所述病变区域图像中的图像特征，由所述图像特征构成病变区域样本集；

(4)从所述病变区域样本集中任选n个样本作为训练样本集，利用所述训练样本集对分类器进行训练得到训练后的分类器，其中，n为正整数且

m表示所述病变区域样本集的样本个数；

(5)将所述病变区域样本集中除所述训练样本集以外的其他样本作为验证样本集，使用所述训练后的分类器对所述验证样本集进行分类，获得各验证样本的脑瘤P53蛋白的表达状况。

与现有技术相比，本发明的有益效果是：

(1)现有技术只能在通过手术获得病理组织后检测P53蛋白表达状况；而本发明所需的脑瘤患者磁共振图像可以在手术前获得，因此本发明在手术前就能够获得P53蛋白表达状况的中间结果信息，具有快速及时的特点。(2)通过采集和分析脑瘤患者的T1加权序列、T1增强序列、FLAIR序列中的任一种序列的磁共振图像来检测P53蛋白的表达状况，可以在获得肿瘤组织切片前无创地获取关于脑瘤恶性程度和患者预后的中间结果信息。(3)通过对病变区域图像集进行图像分析，提取所述病变区域图像中的图像特征，对分类器进行训练和验证，使分类器客观地获得脑瘤P53蛋白的表达状况，避免了在免疫组织化学技术中引入检测者的主观影响，也避免了检测方法标准化不足的问题。(4)通过对训练样本集进行特征优化，可以有效降低特征集的复杂程度，进一步提高检测效率。(5)由于通过图像分析获得P53蛋白表达的状况，不需要消耗化学试剂等，所以具有成本低的特点。

附图说明

图1是本发明的一个实施例中由FLAIR序列的磁共振图像中截取的病变区域图像的示意图；

图2是本发明中对训练样本集进行特征优化后，使用优化后的训练样本集训练支持向量机，然后对新的验证样本集进行分类得到的ROC曲线。

具体实施方式

本发明的基于磁共振图像分析的脑瘤P53蛋白表达检测方法包括下列步骤：

(1)采集脑瘤患者的磁共振图像，其中，磁共振图像包括T1加权序列、T1增强序列、FLAIR序列中的任一种或任几种。具体的采集方法如下：

使用磁共振扫描仪(例如GE Healthcare，1.5T)采集脑胶质瘤患者的横断位、冠状位或矢状位的磁共振图像，该磁共振图像包括T1加权序列、T1增强序列和FLAIR序列。其中，T1加权序列的成像参数优选为Repetition Time＝1966.1ms，Echo Time＝21.088ms，Inversion Time＝750ms；T1增强序列的成像参数优选为Repetition Time＝1967.25ms，Echo Time＝7.264ms，Inversion Time＝750ms；FLAIR序列图像的成像参数优选为Repetition Time＝8002ms，EchoTime＝122.904ms，Inversion Time＝2000ms。优选采集脑胶质瘤患者的横断位的磁共振图像。对患者进行扫描时，每个序列的磁共振图像的体数据分辨率均优选为512×512×16体素，即每个二维切面图像的分辨率为512×512像素，共有16个二维切面图像。磁共振图像的格式一般为DICOM。

(2)在T1加权序列、T1增强序列或FLAIR序列中截取步骤(1)所采集的磁共振图像中的病变区域图像，由截取得到的病变区域图像构成病变区域图像集。具体方法如下：

在T1加权序列、T1增强序列或FLAIR序列的横断位、冠状位或矢状位的磁共振图像的体数据中的每个二维切面图像上至多截取一个病变区域图像(即图1中的白色小方框内的区域)。病变区域图像优选为：尺寸为16×16像素、灰度级为256级、图像格式为tif。后将所截取的病变区域图像构成病变区域图像集。

在获得病变区域图像集后，使用免疫组织化学技术对P53蛋白表达进行检测，用检测结果将病变区域图像标记为P53蛋白表达阳性或阴性。例如，可临床采集22个胶质瘤患者的磁共振图像，其中包括T1加权序列、T1增强序列和FLAIR序列。使用免疫组织化学技术检测发现，22个胶质瘤患者中有10个患者的P53蛋白表达阴性，12个患者的P53蛋白表达阳性。

由于T1加权序列、T1增强序列和FLAIR序列中的部分序列上，磁共振图像的伪影比较严重，容易影响图像分析和分类的结果，因此不同序列的病变区域图像个数不尽相同。例如，在病变区域图像集中，T1加权序列上，P53蛋白表达阳性和阴性的病变区域图像个数均为66个；T1增强序列上，P53蛋白表达阳性和阴性的病变区域图像个数均为67个；FLAIR序列上，P53蛋白表达阳性和阴性的病变区域图像个数均为72个。

(3)通过对病变区域图像集进行图像分析，提取病变区域图像中的图像特征，并由图像特征构成病变区域样本集。以下优选分别对病变区域图像使用灰度共生矩阵提取13个图像特征、使用灰度-梯度共生矩阵提取15个图像特征、使用二维离散正交S变换提取34个图像特征进行进一步说明。

(a)对病变区域图像，使用灰度共生矩阵提取13个图像特征。病变区域图像通常为二维数字图像，以下以一幅二维数字图像为例进行进一步说明。此处以f(x，y)表示一幅二维数字图像，假设其大小为M×N，最高灰度级为第N_g级。在二维数字图像中，在某个方向上相隔一定距离的一对像素点的灰度出现的统计规律，从一定程度上可以反映这个二维数字图像的图像特征。这个统计规律可以用一个矩阵描述，即灰度共生矩阵。

在二维数字图像中，任意取一像素点(x，y)以及偏离它的另一像素点(x+a，y+b)形成一个点对。设该点对的灰度值为(i，j)，即像素点(x，y)的灰度值为i，像素点(x+a，y+b)的灰度值为j。固定a和b，令像素点(x，y)在整幅二维数字图像上移动，则会得到各种(i，j)值。假如二维数字图像的灰度级别为N_g，则i与j的组合共有种。在整幅二维数字图像中，统计每一种组合出现的频率为P(i，j，d，θ)，则构成大小为N_g×N_g的灰度共生矩阵，其中d是点对之间的距离

θ为点对构成的向量与坐标横轴之间的夹角，即点对的方向。灰度共生矩阵本质上就是两个像素点的联合直方图，当a和b取不同的数值组合，都可以得到二维数字图像沿一定方向θ、相隔一定距离

的灰度共生矩阵。灰度共生矩阵的数学表达式为：

P(i，j，d，θ)＝#{(x，y)，(x+a，y+b)∈M×N|f(x，y)＝i，f(x+a，y+b)＝j}

其中#{x}表示集合x中元素的个数。P为N_g×N_g的矩阵。若(x，y)与(x+a，y+b)之间距离为d，点对构成的向量与坐标横轴之间的夹角为θ，则可以得到各种间距及角度的灰度共生矩阵P(i，j，d，θ)。为表达简明起见，在下面的表述中略去d和θ，用P(i，j)表示灰度共生矩阵。

通常可以设置d＝1，分别计算0°、45°、90°、135°四个方向的图像特征，对这四个方向的特征求取均值，得到与方向无关的特征。利用灰度共生矩阵提取的13个图像特征，包括角二阶矩、对比度、相关系数、方差、逆差矩、和平均、和方差、和熵、熵、差平均、惯性、差方差和差熵。在计算上述每个特征时，设置d＝1，分别计算0°、45°、90°、135°四个方向的图像特征，对这四个方向的特征求取均值，得到与方向无关的特征，则每个特征公式如下所示：

1)角二阶矩： $f_1=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}{\left\{P(i,j)\right\}}^2$

2)对比度： $f_2=\underset{n=0}{\overset{N_g-1}{\mathrm{\Σ}}}{n}^2\left\{\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j)\right\},|i-j|=n$

3)相关系数： $f_3=\{\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j)-{\mathrm{\μ}}_x{\mathrm{\μ}}_y\}/{\mathrm{\σ}}_x{\mathrm{\σ}}_y$

式中，μ_x，σ_x分别是{P_x(i)；i＝1，2，…，N_g}的均值和标准差，μ_y，σ_y分别是{P_y(j)；j＝1，2，…，N_g}的均值和标准差。其中，

4)方差： $f_4=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}{(i-\mathrm{\μ})}^{2}P(i,j)$

式中，μ是P(i，j)的均值。

5)逆差矩： $f_5=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}\frac{1}{1+{(i-j)}^2}P(i,j)$

6)和平均： $f_6=\underset{i=2}{\overset{{2N}_g}{\mathrm{\Σ}}}{\mathrm{iP}}_{x+y}\left(i\right)$

其中， $P_{x+y}\left(k\right)=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j),k=\mathrm{2,3},\·\·\·,{2N}_g,|i+j|=k$

7)和方差： $f_7=\underset{i=2}{\overset{{2N}_g}{\mathrm{\Σ}}}{{(1-f_6)}^{2}P}_{x+y}\left(i\right)$

其中， $P_{x+y}\left(k\right)=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j),k=\mathrm{2,3},\·\·\·,{2N}_g,|i+j|=k$

8)和熵： $f_8=-\underset{i=2}{\overset{{2N}_g}{\mathrm{\Σ}}}{P}_{x+y}\left(i\right){\log }_2\left[P_{x+y}\left(i\right)\right]$

其中， $P_{x+y}\left(k\right)=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j),k=\mathrm{2,3},\·\·\·,{2N}_g,|i+j|=k$

9)熵： $f_9=-\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j){\log }_2\left[P(i,j)\right]$

10)差平均： $f_{10}=\underset{i=0}{\overset{N_g-1}{\mathrm{\Σ}}}{\mathrm{iP}}_{x-y}\left(i\right)$

其中， $P_{x-y}\left(k\right)=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j),k=0,1,\·\·\·,N_g-1,|i-j|=k$

11)惯性： $f_{11}=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}{(i-j)}^{2}P(i,j)$

12)差方差： $f_{12}=\underset{i=0}{\overset{N_g-1}{\mathrm{\Σ}}}{{(1-f_{10})}^{2}P}_{x-y}\left(i\right)$

其中， $P_{x-y}\left(k\right)=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j),k=0,1,\·\·\·,N_g-1,|i-j|=k$

13)差熵： $f_{13}=-\underset{i=0}{\overset{N_g-1}{\mathrm{\Σ}}}{P}_{x-y}\left(i\right){\log }_2\left[P_{x-y}\left(i\right)\right]$

其中， $P_{x-y}\left(k\right)=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_g}{\mathrm{\Σ}}}P(i,j),k=0,1,\·\·\·,N_g-1,|i-j|=k$

(b)对病变区域图像，使用灰度-梯度共生矩阵提取15个图像特征。病变区域图像通常为二维数字图像，以下以一幅二维数字图像为例进行进一步说明，此处用f(x，y)表示一幅二维数字图像，假设其大小为M×N，对二维数字图像的灰度矩阵进行正规化变换：

F(x，y)＝[f(x，y)×(N_g-1)/f_max]+1

式中，[x]表示x的整数部分，f_max为二维数字图像的最大灰度值，N_g是正规化变换后的灰度矩阵的最大值。

二维数字图像各像素的梯度计算使用Sobel算子，设二维数字图像的梯度矩阵为g(x，y)，其中x＝1，2，…，N，y＝1，2，…，M。为了使g(x，y)梯度值分布在更大的离散间隔Ns个等级中，我们对二维数字图像的梯度矩阵进行正规化变换：

G(x，y)＝[g(x，y)×(N_s-1)/g_max]+1

式中，[x]表示x的整数部分，g_max为二维数字图像的梯度矩阵的最大值，N_s是正规化变换后梯度矩阵的最大值。

通常，将正规化变换后的二维数字图像的灰度矩阵F(x，y)简称为正规化灰度矩阵，将正规化变换后的二维数字图像的梯度矩阵G(x，y)简称为正规化梯度矩阵，将正规化灰度矩阵和正规化梯度矩阵进行结合就可以得到灰度-梯度共生矩阵：

{H(i，j)；i＝1，2，…，N_g，j＝1，2，…，N_s}

其中，(i，j)表示灰度-梯度共生矩阵的第i行第j列元素，H(i，j)表示正规化灰度矩阵中灰度值为i，并且正规化梯度矩阵中梯度值为j的像素点的个数。N_s是正规化梯度矩阵的最大值，N_g是正规化灰度矩阵的最大值。

将灰度-梯度共生矩阵H(i，j)进行正规化变换，得到：

p(i，j)＝H(i，j)/(N_g×N_s)，i＝1，2，…，N_g，j＝1，2，…，N_s

利用灰度-梯度共生矩阵提取的图像特征主要有15个，分别是灰度平均值、梯度平均值、灰度方差、梯度方差、小梯度优势、大梯度优势、灰度分布的不均匀性、梯度分布的不均匀性、能量、相关系数、灰度熵、梯度熵、混合熵、惯性、逆差矩。上述每个特征的公式如下所示：

1)灰度平均值： $T_1={\mathrm{\μ}}_F=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}i\left[\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}p(i,j)\right]$

2)梯度平均值： $T_2={\mathrm{\μ}}_G=\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}j\left[\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}p(i,j)\right]$

3)灰度方差： $T_3={\mathrm{\σ}}_F={\left\{\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}{(i-{\mathrm{\μ}}_F)}^2\right[\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}p(i,j)\left]\right\}}^{1/2}$

4)梯度方差： $T_4={\mathrm{\σ}}_G={\left\{\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}{(i-{\mathrm{\μ}}_G)}^2\right[\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}p(i,j)\left]\right\}}^{1/2}$

5)小梯度优势： $T_5=\left[\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}\frac{H(i,j)}{j^2}\right]/H,$ 其中 $H=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}H(i,j)$

6)大梯度优势： $T_6=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}{j}^{2}H(i,j)/H$

7)灰度分布的不均匀性： $T_7=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}{\left[\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}H(i,j)\right]}^2/H$

8)梯度分布的不均匀性： $T_8=\underset{i=1}{\overset{n_g}{\mathrm{\Σ}}}\left[\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}H(i,j)\right]/H$

9)能量： $T_9=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}{\left[p(i,j)\right]}^2$

10)相关系数： $T_{10}=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}(i-{\mathrm{\σ}}_F)(j-{\mathrm{\σ}}_G)p(i,j)$

11)灰度熵： $T_{11}=-\left\{\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\right[\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}p(i,j)]{\log }_2\left[\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}p(i,j)\right]\}$

12)梯度熵： $T_{12}=-\left\{\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}\right[\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}p(i,j)]{\log }_2\left[\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}p(i,j)\right]\}$

13)混合熵： $T_{13}=-\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}p(i,j){\log }_{2}p(i,j)$

14)惯性： $T_{14}=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}{(i-j)}^{2}p(i,j)$

15)逆差矩： $T_{15}=\underset{i=1}{\overset{N_g}{\mathrm{\Σ}}}\underset{j=1}{\overset{N_s}{\mathrm{\Σ}}}\frac{1}{1+{(i-j)}^2}p(i,j)$

(c)对病变区域图像，使用二维离散正交S变换提取34个图像特征。病变区域图像通常为二维数字图像，以下以一幅二维数字图像为例进行进一步说明，用f(x，y)表示一幅二维数字图像，假设其大小为M×N。

首先，计算f(x，y)的二维离散傅立叶变换，表示为：

$F(u,v)=\underset{x=1}{\overset{N}{\mathrm{\Σ}}}\underset{y=1}{\overset{M}{\mathrm{\Σ}}}f(x,y)e^{-\left(\sqrt{-1}\right)2\mathrm{\π}(\frac{\mathrm{ux}}{N}+\frac{\mathrm{vy}}{M})}$

其中，u，v分别表示频域中的坐标变量，F(u，v)为二维数字图像的二维离散傅立叶频谱。二维离散傅立叶逆变换可以表示为：

$f(x,y)=\frac{1}{\mathrm{NM}}\underset{u=1}{\overset{N}{\mathrm{\Σ}}}\underset{v=1}{\overset{M}{\mathrm{\Σ}}}F(u,v)e^{-\left(\sqrt{-1}\right)2\mathrm{\π}(\frac{\mathrm{ux}}{N}+\frac{\mathrm{vy}}{M})}.$

其中，变换的结果f(x，y)是二维数字图像。对于大小为N×N的二维数字图像f′(x，y)，其二维离散正交S变换可以通过计算其的二维离散傅立叶变换来求取，步骤如下：

1)计算二维数字图像f′(x，y)的二维离散傅立叶变换，结果为二维离散傅立叶频谱F；

2)对F进行频率子带划分，并且对频率子带进行圆周移位，得到新的频谱。设定频率阶数p_x，p_y∈[2，log₂ N]，且p_x，p_y都是正整数。计算频率子带的频率范围为

然后对频率子带进行圆周移位，表示如下：

是对二维离散傅立叶频谱F进行频率子带划分，并且对频率子带进行圆周移位得到新的频谱。

3)在新的频谱中，对圆周移位后的每个频率子带进行二维离散傅立叶逆变换，得到二维离散正交S变换频谱，公式如下：

$D_{p_x,p_y}[x^{\′},y^{\′}]=\frac{1}{\sqrt{2^{p_x+p_y-2}}}\×\underset{m={-2}^{p_x-2}}{\overset{2^{p_x-2}-1}{\mathrm{\Σ}}}\underset{n=-2^{p_y-2}}{\overset{2^{p_y-2}-1}{\mathrm{\Σ}}}{F}_{p_x,p_y}[\frac{m}{N},\frac{m}{N}]e^{2\mathrm{\π}(\frac{{\mathrm{mx}}^{\′}}{2^{p_x-1}}+\frac{{\mathrm{ny}}^{\′}}{2^{p_y-1}})}.$

其中，x′，y′是二维离散正交S变换频谱的坐标变量，

是二维离散正交S变换频谱。以二维离散正交S变换频谱的中心作为原点，将二维离散正交S变换频谱变换到极坐标系，对极坐标系下半径相等的二维离散正交S变换频谱的分量进行求和，得到二维离散正交S变换谱向量；对二维离散正交S变换谱向量进行以自然对数e为底数的对数变换，然后进行归一化变换，归一化变换表达式如下：

$Z=\frac{Y-M_Y}{S_Y}$

其中，Y是对二维离散正交S变换谱向量进行以自然对数e为底数的对数变换后的结果，M_Y为Y的平均值，S_Y为Y的标准差，z为对Y进行归一化变换后的结果。归一化变换后，得到新的谱向量。新的谱向量的每个分量就是使用二维离散正交S变换提取的二维数字图像的图像特征。

通过步骤(a)、(b)和(c)，每个病变区域图像上可以提取图像特征有62个，其中包括使用灰度共生矩阵提取的图像特征13个，使用灰度-梯度共生矩阵提取的图像特征15个，使用二维离散正交S变换提取的图像特征34个。病变区域图像的图像特征的详见表1。对病变区域图像集中的每一个病变区域图像提取62个图像特征，并且将每一个病变区域图像的P53蛋白表达阳性或阴性标记进行量化，具体为：如果病变区域图像标记为P53蛋白表达阳性，则将该病变区域图像的类别标记为+1；如果病变区域图像标记为P53蛋白表达阴性，则将该病变区域图像的类别标记为-1。通过上述处理，所有病变区域图像的图像特征和类别标记结合，构成病变区域样本集。病变区域样本集可以表示为S＝{(f_i，l_i)|i＝1，2，…，n}，其中f_i是第i个病变区域图像的62个的图像特征，表示为f_i＝[f_i1，f_i2，…，f_ij，…，f_i62]，f_ij是第i个病变区域图像的第j个图像特征的特征值；l_i是第i个病变区域图像的类别标记；如果第i个病变区域图像标记为P53蛋白表达阳性，则l_i＝+1；如果第i个病变区域图像标记为P53蛋白表达阴性，则l_i＝-1；n是病变区域样本集的样本个数。

表1.病变区域图像的图像特征

(4)从所述病变区域样本集中任选n个样本作为训练样本集，利用所述训练样本集对分类器进行训练得到训练后的分类器，其中，n为正整数，且

m表示所述病变区域样本集的样本个数。

以下具体举例说明：在病变区域样本集中，T1加权序列上，P53蛋白表达阳性和阴性的病变区域样本个数均为66个；T1增强序列上，P53蛋白表达阳性和阴性的病变区域样本个数均为67个；FLAIR序列上，P53蛋白表达阳性和阴性的病变区域样本个数均为72个。那么从病变区域样本集中任选的样本构成训练样本集，训练样本集的样本个数优选如下：T1加权序列上，P53蛋白表达阳性和阴性的训练样本个数均为33个；T1增强序列上，P53蛋白表达阳性和阴性的训练样本个数均为34个；FLAIR序列上，P53蛋白表达阳性和阴性的训练样本个数均为36个。

选择支持向量机的参数，并且利用训练样本集训练支持向量机。支持向量机的核函数优选为径向基函数，惩罚因子C和核宽度σ两个参数采用网格搜索的方法进行选取。网格搜索的方法具体为，将惩罚因子C和核宽度σ分别取N个值和M个值，对N×M个C和σ的组合分别训练不同的支持向量机，采用交叉检验方法评估支持向量机的推广能力，选择最高分类准确率所对应的C和σ的组合，作为支持向量机的最优参数。其中，交叉检验方法，通常在样本个数较少的样本集进行分类时使用，以下以一个样本集为例进行进一步说明：交叉检验方法就是把样本集随机分成k个部分，选择一个部分作为验证样本，剩下的k-1个部分作为训练样本，重复使用训练样本训练支持向量机和对验证样本进行验证的过程，一直持续到每个部分都被当作一次验证样本；最后根据k个部分的每个样本的分类结果计算交叉检验方法得到的分类准确率，用于评估支持向量机的推广能力。需要注意的是，参数选择过程中仅使用训练样本集，本发明优选采用N＝5，M＝6，惩罚因子的取值范围优选为C∈{2^-1，2⁰，2¹，2²，2³}，核宽度的取值范围优选为σ∈{2^-4，2^-3，2^-2，2^-1，2⁰，2¹}，交叉验证方法中k的取值优选为5。

然后，根据最优参数设定支持向量机的参数，使用T1加权序列、T1增强序列或FLAIR序列中任一种序列的训练样本集对支持向量机进行训练，得到训练后的支持向量机。

如果对训练样本集进行特征优化，则可以有效降低特征集的复杂程度，减小分类器的运算负担，提高检测的效率。在对训练样本集进行特征优化后，重新选择支持向量机的参数，并重新训练支持向量机。

对训练样本集进行特征优化时，训练样本集可以表示为S_m＝{(f_i，l_i)|i＝1，2，…，n}，其中f_i是第i个训练样本的m维的特征向量，表示为f_i＝[f_i1，f_i2，…，f_ij，…，f_im]，f_ij是第i个训练样本的第j个特征的特征值；l_i是第i个训练样本的类别标记；对于P53蛋白表达阳性的训练样本，l_i＝1；对于P53蛋白表达阴性的训练样本，l_i＝…1；n是训练样本数量。特征优化的步骤如下：

1)将待优化的训练样本集表示为训练样本集S_p，S_p中的特征有p维，令p＝m；

2)从S_p中取出所有训练样本的第q维特征f_iq，i＝1，2，…，n，得到新的训练样本集

令q在1，2，…，p中逐个取值，得到多个新的训练样本集

使用k＝1的最邻近法对新的训练样本集

分别进行留一法分类，得到分类误差其中，最小分类误差所对应的新的训练样本集为p-1维训练样本集，即S_p-1。

3)判断训练样本集S_p-1的特征维数是否大于或等于2。如果是，则令p＝p-1，将步骤2)得到的p-1维训练样本集作为待优化的训练样本集，重复步骤2)；如果否，则进入步骤4)；

4)经过上述步骤，得到训练样本集S_m，S_m-1，…，S₂，S₁。使用支持向量机分别对训练样本集S_m，S_m-1，…，S₂，S₁进行留一法分类，得到分类准确率P_m，P_m-1，…，P₂，P₁。

5)选取最高的分类准确率对应的训练样本集S，则S包括的特征为最优特征集。

在特征优化的步骤中，留一法，即留一交叉检验法(leave-one-out cross validation strategy，LOOCV)，通常对样本个数较少的样本集进行分类时使用留一法，以下以一个样本集为例进行进一步说明。只使用样本集中的一个样本作为验证样本，而剩余样本作为训练样本，重复训练支持向量机和对验证样本进行验证的过程，一直持续到每个样本都被当作一次验证样本；根据每个样本的验证结果计算得到留一法的分类误差或者分类准确率。步骤4)中支持向量机的核函数优选为径向基函数，惩罚因子C和核宽度σ两个参数采用网格搜索的方法进行选取，惩罚因子的取值范围优选为C∈{2^-1，2⁰，2¹，2²，2³}，核宽度的取值范围优选为σ∈{2^-4，2^-3，2^-2，2^-1，2⁰，2¹}，网格搜索中交叉检验方法的k取值优选为5。

通过对训练样本集进行特征优化，特征维数大幅下降，特征的复杂程度得到有效地降低。例如，本发明中T1加权序列、T1增强序列、FLAIR序列上的训练样本集的特征维数为62。经过特征优化，T1加权序列上优化后的训练样本集的特征维数为4，T1增强序列上优化后的训练样本集的特征维数为10，FLAIR序列上优化后的训练样本集的特征维数为14，特征优化大幅降低了特征的复杂程度。

特征优化后，使用优化后的训练样本集对支持向量机进行训练。首先，选择支持向量机的参数。支持向量机的核函数优选为径向基函数，惩罚因子C和核宽度σ两个参数采用网格搜索的方法进行选取。网格搜索的方法具体为，将惩罚因子C和核宽度σ分别取N个值和M个值，对N×M个C和σ的组合分别训练不同的支持向量机，采用交叉检验方法评估支持向量机的推广能力，选择最高分类准确率所对应的C和σ的组合，作为支持向量机的最优参数。需要注意的是，参数选择过程中仅使用优化后的训练样本集，本发明优选采用N＝5，M＝6，惩罚因子的取值范围优选为C∈{2^-1，2⁰，2¹，2²，2³}，核宽度的取值范围优选为σ∈{2^-4，2^-3，2^-2，2^-1，2⁰，2¹}，交叉检验方法中k的取值优选为5。然后，根据最优参数，设定支持向量机的参数，使用T1加权序列、T1增强序列或FLAIR序列中任一种序列的优化后的训练样本集对支持向量机进行训练，得到优化训练后的支持向量机。

(5)将病变区域样本集中除所述训练样本集以外的其他样本作为验证样本集，使用训练后的分类器对所述验证样本集进行分类，获得各验证样本的脑瘤P53蛋白的表达状况。

以下具体举例说明：在病变区域样本集中，T1加权序列上，P53蛋白表达阳性和阴性的病变区域样本个数均为66个，如果从中任选P53蛋白表达阳性和阴性的病变区域样本均为33个作为训练样本集，那么验证样本集中的P53蛋白表达阳性和阴性的样本个数均为33个；T1增强序列上，P53蛋白表达阳性和阴性的病变区域样本个数均为67个，如果从中任选P53蛋白表达阳性和阴性的病变区域样本均为34个作为训练样本集，那么验证样本集中的P53蛋白表达阳性和阴性的样本个数均为33个；FLAIR序列上，P53蛋白表达阳性和阴性的病变区域样本个数均为72个，如果从中任选P53蛋白表达阳性和阴性的病变区域样本均为36个作为训练样本集，那么验证样本集中的P53蛋白表达阳性和阴性的样本个数均为36个。

如果未对训练样本集进行特征优化，则使用训练后的支持向量机对验证样本进行分类，可以获得各验证样本的脑胶质瘤P53蛋白表达状况。需要注意的是，若步骤(4)中训练支持向量机采用T1加权序列的训练样本，则验证样本必须取自T1加权序列；若步骤(4)中训练支持向量机采用T1增强序列的训练样本，则验证样本必须取自T1增强序列；若步骤(4)中训练支持向量机采用FLAIR序列的训练样本，则验证样本必须取自FLAIR序列。

如果在分类之前对训练样本集进行特征优化，则需要使用优化后的训练样本集对支持向量机进行训练，构成优化训练后的支持向量机。然后，根据优化后的训练样本集所保留的特征选取验证样本集的特征，构成新的验证样本集。最后，使用优化训练后的支持向量机对新的验证样本进行分类，获得每个新的验证样本的脑胶质瘤P53蛋白表达状况。需要注意的是，若步骤(4)中训练支持向量机采用T1加权序列的特征优化后的训练样本，则新的验证样本必须取自T1加权序列；若步骤(4)中训练支持向量机采用T1增强序列的特征优化后的训练样本，则新的验证样本必须取自T1增强序列；若步骤(4)中训练支持向量机采用FLAIR序列的特征优化后的训练样本，则新的验证样本必须取自FLAIR序列。

另外，如果未对训练样本集进行特征优化，则使用训练后的支持向量机对验证样本进行分类，仍然可以检测脑胶质瘤P53蛋白表达状况。由于未经特征优化的训练样本集给分类器带来了很大的计算负担。例如，T1加权序列上，验证样本集中的P53蛋白表达阳性和阴性的样本个数均为33个，那么使用训练后的支持向量机对验证样本进行分类，需要0.28秒，而使用优化训练后的支持向量机对新的验证样本进行分类，只需要0.14秒；T1增强序列上，验证样本集中的P53蛋白表达阳性和阴性的样本个数均为33个，那么使用训练后的支持向量机对验证样本进行分类，需要0.29秒，而使用优化训练后的支持向量机对新的验证样本进行分类，只需要0.14秒；FLAIR序列上，验证样本集中的P53蛋白表达阳性和阴性的样本个数均为36个，那么使用训练后的支持向量机对验证样本进行分类，需要0.30秒，而使用优化训练后的支持向量机对新的验证样本进行分类，只需要0.15秒。因此，为了提高方法的效率，本发明使用特征优化后的训练样本集对支持向量机进行训练，并使用优化训练后的支持向量机对新的验证样本进行分类。

最后，本发明可使用受试者工作特征(receiver operating characteristic，ROC)评价分类结果的准确性，使用ROC曲线下面积A_z来反映分类的有效性。当A_z取值等于0.5时，分类没有意义，相当于随机分类；当A_z取值等于1是，分类完善，对所有样本都能准确分类。A_z值越大，则分类效果越好。

经过计算得到，对来自T1加权序列的新的验证样本进行分类，A_z达到0.7808；对来自T1增强序列的新的验证样本进行分类，A_z达到0.8128；对来自FLAIR序列的新的验证样本进行分类，A_z达到0.6997。

从图2可以看出，本发明中，支持向量机对新的验证样本都能做到准确的分类，尤其在T1增强序列上具有比较高的识别率。可见，本发明能够准确的识别脑胶质瘤P53蛋白表达状况。

Claims (1)

1.一种脑瘤P53蛋白表达检测方法，其特征在于，包括包括以下步骤：

（1）采集脑瘤患者的磁共振图像，所述磁共振图像包括T1加权序列、T1增强序列、FLAIR序列中的任一种或任几种；

（2）在步骤（1）的任一种所述序列中截取所述磁共振图像中的病变区域图像，由所述病变区域图像构成病变区域图像集，并将病变区域图像标记为P53蛋白表达阳性或阴性；

（3）通过对所述病变区域图像集进行图像分析，提取所述病变区域图像中的图像特征，由所述图像特征构成病变区域样本集；

（4）从所述病变区域样本集中任选                                                

个样本作为训练样本集，利用所述训练样本集对分类器进行训练得到训练后的分类器，其中，

为正整数且

,

表示所述病变区域样本集的样本个数；

（5）将所述病变区域样本集中除所述训练样本集以外的其他样本作为验证样本集，使用所述训练后的分类器对所述验证样本集进行分类，获得各验证样本的脑瘤P53蛋白的表达状况。

Images