CN115100093A - 一种基于梯度滤波的医学图像融合方法 - Google Patents
一种基于梯度滤波的医学图像融合方法 Download PDFInfo
- Publication number
- CN115100093A CN115100093A CN202210903636.1A CN202210903636A CN115100093A CN 115100093 A CN115100093 A CN 115100093A CN 202210903636 A CN202210903636 A CN 202210903636A CN 115100093 A CN115100093 A CN 115100093A
- Authority
- CN
- China
- Prior art keywords
- image
- layer
- fusion
- gradient
- images
- 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.)
- Withdrawn
Links
- 238000001914 filtration Methods 0.000 title claims abstract description 24
- 238000007500 overflow downdraw method Methods 0.000 title claims abstract description 15
- 230000004927 fusion Effects 0.000 claims abstract description 58
- 238000000034 method Methods 0.000 claims abstract description 13
- 238000013528 artificial neural network Methods 0.000 claims abstract description 8
- 230000008878 coupling Effects 0.000 claims abstract description 4
- 238000010168 coupling process Methods 0.000 claims abstract description 4
- 238000005859 coupling reaction Methods 0.000 claims abstract description 4
- 230000002708 enhancing effect Effects 0.000 claims abstract description 4
- 210000002569 neuron Anatomy 0.000 claims description 30
- 238000000354 decomposition reaction Methods 0.000 claims description 16
- 239000011159 matrix material Substances 0.000 claims description 15
- 230000000694 effects Effects 0.000 claims description 11
- 238000013178 mathematical model Methods 0.000 claims description 9
- 238000004364 calculation method Methods 0.000 claims description 7
- 230000003044 adaptive effect Effects 0.000 claims description 5
- 238000012935 Averaging Methods 0.000 claims description 3
- 230000003213 activating effect Effects 0.000 claims description 3
- 210000000857 visual cortex Anatomy 0.000 claims description 3
- 238000010304 firing Methods 0.000 claims description 2
- 238000012886 linear function Methods 0.000 claims description 2
- 101150040772 CALY gene Proteins 0.000 claims 1
- 238000005457 optimization Methods 0.000 abstract description 2
- 238000002591 computed tomography Methods 0.000 description 5
- 230000008901 benefit Effects 0.000 description 3
- 238000013527 convolutional neural network Methods 0.000 description 3
- AOYNUTHNTBLRMT-SLPGGIOYSA-N 2-deoxy-2-fluoro-aldehydo-D-glucose Chemical compound OC[C@@H](O)[C@@H](O)[C@H](O)[C@@H](F)C=O AOYNUTHNTBLRMT-SLPGGIOYSA-N 0.000 description 2
- 108010043121 Green Fluorescent Proteins Proteins 0.000 description 2
- 102000004144 Green Fluorescent Proteins Human genes 0.000 description 2
- 210000004556 brain Anatomy 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 239000005090 green fluorescent protein Substances 0.000 description 2
- 238000003384 imaging method Methods 0.000 description 2
- 230000003902 lesion Effects 0.000 description 2
- 238000002595 magnetic resonance imaging Methods 0.000 description 2
- 238000000513 principal component analysis Methods 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 230000000007 visual effect Effects 0.000 description 2
- 101000941994 Homo sapiens Protein cereblon Proteins 0.000 description 1
- OAICVXFJPJFONN-UHFFFAOYSA-N Phosphorus Chemical compound [P] OAICVXFJPJFONN-UHFFFAOYSA-N 0.000 description 1
- 102100032783 Protein cereblon Human genes 0.000 description 1
- 230000000740 bleeding effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000012512 characterization method Methods 0.000 description 1
- 230000000052 comparative effect Effects 0.000 description 1
- 230000000295 complement effect Effects 0.000 description 1
- 238000009792 diffusion process Methods 0.000 description 1
- 230000009977 dual effect Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000010365 information processing Effects 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 230000004066 metabolic change Effects 0.000 description 1
- 230000004060 metabolic process Effects 0.000 description 1
- 230000008520 organization Effects 0.000 description 1
- 238000002600 positron emission tomography Methods 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 238000002603 single-photon emission computed tomography Methods 0.000 description 1
- 210000004872 soft tissue Anatomy 0.000 description 1
- 238000001931 thermography Methods 0.000 description 1
- 210000001519 tissue Anatomy 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/50—Image enhancement or restoration using two or more images, e.g. averaging or subtraction
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/20—Image enhancement or restoration using local operators
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- 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/20—Special algorithmic details
- G06T2207/20081—Training; Learning
-
- 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/20—Special algorithmic details
- G06T2207/20084—Artificial neural networks [ANN]
-
- 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/20—Special algorithmic details
- G06T2207/20212—Image combination
- G06T2207/20221—Image fusion; Image merging
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Image Processing (AREA)
Abstract
本发明公开了一种基于梯度滤波的医学图像融合方法，首先输入两张已配准的多模态灰度医学图像，然后进行增强，设置阈值函数区分边缘和纹理的细节或噪声，并分解为近似层和残差层；然后近似层图像采用一种能量属性的低频融合策略得到近似层图像，残差层采用一种参数自适应的脉冲耦合神经网络的高频融合策略得到残差层图像；最后将近似层和残差层的融合子图像重建获得最终的融合图像。本发明解决了现有技术中存在的融合图像因为没有考虑对比度和噪点导致目标的融合图像不清晰，而通过基于模糊梯度阈值函数和全局优化的梯度滤波器分解原图像，融合时考虑不同尺度特征来保存良好的对比度和详细信息，有效降低边缘模糊和噪声造成的影响。
Description
技术领域
本发明属于医学图像融合技术领域，具体涉及一种基于梯度滤波的医学图像融合方法。
背景技术
近年来，医学图像处理已经逐渐成为图像领域中十分重要的研究问题，广受人们关注。该领域随着技术的发展，出现了许多先进的成像设备，它们可以提供不同形式的医学图像。例如计算机断层扫描(CT)图像可以清晰地呈现密集的结构信息，但是组织表征受限；核磁共振成像(MRI)可以清晰的软组织等，但是缺乏体代谢的动态信息，其中MR-T1、MR-T2等又反应不同信息；正电子发射断层扫描(PET)和单光子发射计算机断层扫描(SPECT)图像作为功能图像提供流血和代谢变化信息。使断层扫描图像信息更精确氟脱氧葡萄糖(FDG)和绿色荧光蛋白(GFP)等热成像和视觉光成像技术已被开发出来。因此，医学图像融合是将提取到的各个图像的特征信息整合的有效方法。
传统在空间域上有主成分分析(PCA)等可以用于图像融合，而在变换域方法中大多数属于多尺度变换(MST)，MST一般包括三个步骤：图像分解、系数融合和图像重建。存在许多经典方法，例如拉普拉斯金字塔(LP)，小波变换(WT)、双树复小波变换(DTCWT)和非下采样轮廓波变换(NSCT)等。NSCT也有不足之处，手工制作的频率子带融合规则存在一些局限性。由于医学图像存在不同亮度信息，需要考虑目标对比度和亮度信息。对于融合规则的选取，卷积神经网络(CNN)由于其出色的性能已成为图像处理和医学图像融合中的流行算法，但是不能用于空间域，否则会丢失源图像中的重要信息。而脉冲耦合神经网络(PCNN)与人类视觉系统的信息处理机制密切相关，由于医学图像具有更复杂的边缘、轮廓和细节，使得基于PCNN的融合策略比传统的直接基于分解系数的融合策略具有更大的优势。
发明内容
本发明的目的是提供一种基于梯度滤波的医学图像融合方法，解决了现有技术中存在的融合图像因为没有考虑对比度和噪点导致目标的融合图像不清晰，而通过基于模糊梯度阈值函数和全局优化的梯度滤波器分解原图像，融合时考虑不同尺度特征来保存良好的对比度和详细信息，有效降低边缘模糊和噪声造成的影响。
本发明所采用的技术方案是，基于梯度滤波的医学图像融合方法，具体按照以下步骤实施：
步骤1、输入两张已配准的多模态灰度医学图像；
步骤2、对于步骤1中的两张医学图像进行增强，获得细节和边缘信息更丰富的图像；
步骤3、设置阈值函数区分边缘和纹理的细节或噪声，并根据梯度阈值设计一个惩罚系数，根据对步骤2中增强的图像梯度滤波分解成为近似层和残差层；
步骤4、根据步骤3中分解得到的近似层和残差层分别设置融合规则，近似层图像采用一种能量属性的低频融合策略得到近似层图像，残差层采用一种参数自适应的脉冲耦合神经网络的高频融合策略得到残差层图像；
步骤5、将步骤4中近似层和残差层的融合子图像重建获得最终的融合图像。
本发明的特点还在于，
步骤2具体按照以下步骤实施：
步骤2.1、首先对步骤1中的两张医学图像3×3像素为窗口内像素的局部差值信息D(i,j)：
其中D为窗口内像素的局部差值信息，(i,j)为图像像素横纵坐标，Z(x)表示每窗口内各点像素；
步骤2.2、计算图像的局部梯度信息：
其中g1(x)为图像的局部梯度信息，P(i,j)表示源图像的像素，θ和δ是g1的行和列，重新映射g1得：
Y2(x(θ,δ))＝ω2x(θ,δ)+ω2+1 (3)
其中Y2为重新映射g1的值，ω2表示增强强度参数，ω2越大表示图像整体增强水平越明显，x(θ,δ)表示Y2的自变量，
步骤2.3、得到加权系数g2：
增强的图像为：
其中H表示增强图像，⊙表示卷积操作，A(i,j)是步骤1中的两张医学图像。
步骤3具体按照以下步骤实施：
步骤3.1、构建边缘保留滤波器，优化以下函数：
其中，O*是需要优化的目标函数，I是输入图像，O是滤波后的输出图像，对于m×n的灰度图像，c表示像素坐标c＝(x,y)，
和
分别是x和y目标图像的方向，|·|2表示L2范数；
步骤3.2、设置一个阈值来区分边缘和纹理或噪声，惩罚系数设置为：
其中P为惩罚系数，λ是惩罚系数的最大值，g∈{x,y}表示x或者y方向，α是梯度阈值，n是阈值函数斜率，thresh定义为阈值函数：
步骤3.3、得到最后的目标函数为：
步骤3.4、为求解目标函数O*，将上式写为：
其中o和i分别是输出和输入图像的列向量形式，T是转置符号，Dx和Dy分别是在x和y方向上的Toeplitz矩阵，Px和Py是对角矩阵；
步骤3.5、使用加权最小二乘法WLS，对于变量的导数获得o并将导数设置为等于0，设置成为迭代形式并求解下面线性函数：
{E+Lt}ot+1＝i (12)其中E是单位矩阵，Lt是拉普拉斯矩阵，并且o0＝i，
t是迭代次数，
步骤3.6、根据上式设置滤波器为O＝GletF(I,t,α,λ,n)；
使用梯度滤波获得的近似层图像An表示如下：
An＝GletF(S,f,α,λ,k) (13)
其中An为近似层图像，S是步骤2中得到的增强图像，f是迭代次数，α是梯度阈值，λ是最大惩罚系数，k是阈值函数斜率，
步骤3.7、从增强图像中减去近似层图像以获得相应的残差层图像Rn：
Rn＝Sn-An (14)
其中n表示第n幅图像。
步骤4具体按照以下步骤实施：
步骤4.1、对于图像的残差层而言，使用参数自适应脉冲耦合神经网络PA-PCNN模型用于图像融合，该模型基于改进的简化PCNN模型，受生物视觉皮层模型的启发而产生，是由若干个神经元互联所构成的反馈型网络，其每一神经元由分支树、调制耦合器和脉冲产生器构成，分支树馈送的数学模型如下：
Fij[n]＝Rij (15)
步骤4.2、分支树链接输入的数学模型如下：
步骤4.3、调制耦合器的数学模型如下：
步骤4.4、PCNN的输出模块被称为脉冲产生器，脉冲产生器的变阈值函数如下：
步骤4.5、硬限幅函数的表达式如下：
其中，下标ij、kl为神经元的标号，Rij为输入图像的强度，Yij、Lij、Uij和θij分别为神经元的反馈输入、链接输入、内部活动项和动态阈值，Wijkl为根据八个相邻神经元的先前放电状态设置的链接权矩阵，VL、Vθ为幅度常数，αf、αθ为相应的衰减系数，β为链接强度，n为迭代次数，Yij为二值输出；
步骤4.6、使用PA-PCNN产生一次自适应计算获得新的参数值，具体的参数自适应计算公式如下：
其中σ(R)表示范围为[0,1]的输入图像R的标准偏差，R'和Rmax分别表示归一化的大津阈值和输入图像的最大强度；
步骤4.7、确定迭代次数n，结合公式(3)至(8)对PA-PCNN进行初始化，即将Uij[0]、Yij[0]和θij[0]值设为0，待初始化完成，第一次迭代中所有的非零强度神经元都将被激活，在此基础上依次完成n次迭代计算，最后统计n次迭代中每一个神经元的点火次数，每次点火时间描述为：
Tij[n]＝Tij[n-1]+Yij[n] (25)
其中Tij[n]表示n次点火时间，整个迭代完成后，神经元的总点火次数Tij[N]用于表征对应高频分量输入的活动水平，Tij[N]数值越大，则此神经元的活动水平越高，融合带的规则如下：
其中HF表征输入的融合系数，A和B分别为步骤3中分解得到的两幅残差层图像，HA,B与步骤3分解得到残差层图像的尺寸相同，新融合子带的融合系数选择以神经元点火总次数为依据，将两幅图的高频分量所对应同一位置像素点(i,j)的Tij[N]做比较，若TA,ij[N]＞TB,ij[N]，则新融合子带取HA(i,j)；反之，新融合子带取HB(i,j)，最终得到残差层融合图像FR。
步骤4.8、对于图像的近似层而言，采用一种能量属性的融合策略：
首先计算低频分量的固有属性值，定义固有属性值IPA,B的数学公式如下：
其中，IPA,B为低频分量的固有属性值，A,B为步骤三中对近似层分解得到的两幅图像，μ是低频分量的平均值；M是低频分量的中值，其次计算低频分量的能量属性函数，能量属性函数定义如下：
其中，EA,B为步骤3分解得到的两幅近似层图像的能量属性函数，a为调制参数，LA,B与步骤3分解得到近似层图像的尺寸相同，融合得到新的低频子带，融合后的低频子带通过能量属性函数和原低频分量加权平均得到，具体公式如下：
其中LF(i,j)即是新的低频融合子带系数，最终得到近似层融合图像FA。
步骤5具体按照以下步骤实施：
将步骤4中近似层和残差层的融合子图像重建获得最终的融合图像F：
F＝FA+FR (32)
其中F为最终融合图像，FA为近似层融合子图像，FR为残差层融合子图像。
本发明的有益效果是，基于梯度滤波的医学图像融合方法，作为多尺度变换分为三个部分，由于传统将图像分解为结构信息和细节信息，但是普通的平均加权滤波对局部亮度进行平均，导致亮度扩散和边缘模糊，因此考虑使用梯度滤波方法，分解成为存在亮度分布的图像和具有微小梯度的图像。由于医学图像本身图像特性，图像特征信息可能不太明显，准确高效地获取医学图像特征信息、改善图像纹理细节、提升对比度、减少噪声显得尤为重要，因此先对源图像进行的图像增强，使源图像的纹理和边缘信息更为明显。而对于梯度滤波器的使用考虑到了图像的像素强度分布和图像的总梯度，为了区分图像的总梯度边缘和纹理或噪声，引入一个阈值，并根据梯度阈值设计了一个惩罚系数。对于融合规则，PCNN模型可以处理像素强度，它还能够区分不同源图像的高频系数的绝对，所以PCNN在图像融合中的优势即相邻像素的活动水平信息可以相互传输，以实现鲁棒性更强的活动测量。
附图说明
图1是本发明基于梯度滤波的医学图像融合方法；
图2是残差层融合规则流程图；
图3是本发明PCNN简化模型方法；
图4是近似层融合规则流程图；
图5是本发明实例1中分解图像随梯度阈值α参数变化结果图；
图6是本发明实例1中分解图像随最大惩罚系数λ和阈值函数斜率n参数变化结果图；
图7是本发明实例1中其他方法与本方法的对比实验结果。
具体实施方式
下面结合附图和具体实施方式对本发明进行详细说明。
本发明基于梯度滤波的医学图像融合方法，流程图如图1所示，具体按照以下步骤实施：
步骤1、输入两张已配准的多模态灰度医学图像；步骤1中输入的多模态医学图像可以是CT、MRI等图像，MR图像也存在不同的图像类型，多用CT图像与MR图像融合，互补组织结构和细节。对于有病灶的图像，例如脑部图像，MRT2可很好观察病变。
步骤2、对于步骤1中的两张医学图像进行增强，获得细节和边缘信息更丰富的图像；
步骤2具体按照以下步骤实施：
步骤2.1、首先对步骤1中的两张医学图像3×3像素为窗口内像素的局部差值信息D(i,j)：
其中D为窗口内像素的局部差值信息，(i,j)为图像像素横纵坐标，Z(x)表示每窗口内各点像素；
步骤2.2、计算图像的局部梯度信息：
其中g1(x)为图像的局部梯度信息，P(i,j)表示源图像的像素，θ和δ是g1的行和列，重新映射g1得：
Y2(x(θ,δ))＝ω2x(θ,δ)+ω2+1 (3)
其中Y2为重新映射g1的值，ω2表示增强强度参数，ω2越大表示图像整体增强水平越明显，x(θ,δ)表示Y2的自变量，
步骤2.3、得到加权系数g2：
增强的图像为：
其中H表示增强图像，⊙表示卷积操作，A(i,j)是步骤1中的两张医学图像。
步骤3、设置阈值函数区分边缘和纹理的细节或噪声，并根据梯度阈值设计一个惩罚系数，根据对步骤2中增强的图像梯度滤波分解成为近似层和残差层；
步骤3具体按照以下步骤实施：
构建边缘保留滤波器，选择保留输入图像的亮度分布并消除小的梯度变化。对于输出图像，应使其像素强度分布与给定输入图像的像素强度分布相似，还应该限制滤波输出图像的总梯度，优化以下函数：
其中O*是需要优化的目标函数，I是输入图像，O是滤波后的输出图像，对于m×n的灰度图像，c表示像素坐标c＝(x,y)，
和
分别是x和y目标图像的方向，||·||2表示L2范数。
设置一个阈值来区分边缘和纹理或噪声，惩罚系数设置为：
其中P为惩罚系数，λ是惩罚系数的最大值，g∈{x,y}表示x或者y方向，α是梯度阈值，n是阈值函数斜率，thresh定义为阈值函数：
得到最后的目标函数为：
为求解目标函数O*可将上式写为：
其中o和i分别是输出和输入图像的列向量形式，T是转置符号，Dx和Dy分别是在x和y方向上的Toeplitz矩阵，Px和Py是对角矩阵。
使用加权最小二乘法(WLS)，对于变量的导数来获得o并将导数设置为等于0.设置成为迭代形式并求解下面线性函数：
{E+Lt}ot+1＝i (12)
其中E是单位矩阵，Lt是拉普拉斯矩阵，并且o0＝i，
t是迭代次数，根据上式设置滤波器为O＝GletF(I,t,α,λ,n).
梯度滤波将图像分解为具有亮度信息的近似层和梯度信息的残差层，使用梯度滤波获得的近似层图像An表示如下：
An＝GletF(S,f,α,λ,k) (13)
其中An为近似层图像，S是步骤2中得到的增强图像，f是迭代次数，α是梯度阈值，λ是最大惩罚系数，k是阈值函数斜率，从增强图像中减去近似层图像以获得相应的残差层图像Rn：
Rn＝Sn-An (14)
其中n表示第n幅图像。
步骤4、根据步骤3中分解得到的近似层和残差层分别设置融合规则，近似层图像采用一种能量属性的低频融合策略得到近似层图像，残差层采用一种参数自适应的脉冲耦合神经网络的高频融合策略得到残差层图像；
步骤4具体按照以下步骤实施：
步骤4.1，根据分解得到的图像设置融合规则，对于图像的残差层而言，融合规则如图2所示。存在图像更多的高频信息，使用参数自适应脉冲耦合神经网络PA-PCNN模型用于图像融合，该模型基于改进的简化PCNN模型，受生物视觉皮层模型的启发而产生，是由若干个神经元互联所构成的反馈型网络，其每一神经元由分支树、调制耦合器和脉冲产生器构成，如图3所示，分支树馈送的数学模型如下：
Fij[n]＝Rij (15)
分支树链接输入的数学模型如下：
调制耦合器的数学模型如下：
PCNN的输出模块被称为脉冲产生器，脉冲产生器的变阈值函数如下：
硬限幅函数的表达式如下：
其中，下标ij、kl为神经元的标号，Rij为输入图像的强度，Yij、Lij、Uij和θij分别为神经元的反馈输入、链接输入、内部活动项和动态阈值，Wijkl为根据八个相邻神经元的先前放电状态设置的链接权矩阵，VL、Vθ为幅度常数，αf、αθ为相应的衰减系数，β为链接强度，n为迭代次数，Yij为二值输出；
使用PA-PCNN产生一次自适应计算获得新的参数值，具体的参数自适应计算公式如下：
其中σ(R)表示范围为[0,1]的输入图像R的标准偏差，R'和Rmax分别表示归一化的大津阈值和输入图像的最大强度；
确定迭代次数n，结合公式(3)至(8)对PA-PCNN进行初始化，即将Uij[0]、Yij[0]和θij[0]值设为0，待初始化完成，第一次迭代中所有的非零强度神经元都将被激活，在此基础上依次完成n次迭代计算，最后统计n次迭代中每一个神经元的点火次数，每次点火时间可以描述为：
Tij[n]＝Tij[n-1]+Yij[n] (25)
其中Tij[n]表示n次点火时间，整个迭代完成后，神经元的总点火次数Tij[N]用于表征对应高频分量输入的活动水平，Tij[N]数值越大，则此神经元的活动水平越高，融合带的规则如下：
其中HF表征输入的融合系数，A和B分别为步骤3中分解得到的两幅残差层图像，HA,B与步骤3分解得到残差层图像的尺寸相同，新融合子带的融合系数选择以神经元点火总次数为依据，将两幅图的高频分量所对应同一位置像素点(i,j)的Tij[N]做比较，若TA,ij[N]＞TB,ij[N]，则新融合子带取HA(i,j)；反之，新融合子带取HB(i,j)，最终得到残差层融合图像FR。
步骤4.2、对于图像的近似层而言，更多存在图像的低频信息，即采用一种能量属性的融合策略，如图4所示：
首先计算低频分量的固有属性值，定义固有属性值IPA,B的数学公式如下：
其中，IPA,B为低频分量的固有属性值，A,B为步骤三中对近似层分解得到的两幅图像，μ是低频分量的平均值；M是低频分量的中值，其次计算低频分量的能量属性函数，能量属性函数定义如下：
其中，EA,B为步骤3分解得到的两幅近似层图像的能量属性函数，a为调制参数，LA,B与步骤3分解得到近似层图像的尺寸相同，融合得到新的低频子带，融合后的低频子带通过能量属性函数和原低频分量加权平均得到，具体公式如下：
其中LF(i,j)即是新的低频融合子带系数，最终得到近似层融合图像FA。
步骤5、将步骤4中近似层和残差层的融合子图像重建获得最终的融合图像。
步骤5具体按照以下步骤实施：
将步骤4中近似层和残差层的融合子图像重建获得最终的融合图像F：
F＝FA+FR (32)
其中F为最终融合图像，FA为近似层融合子图像，FR为残差层融合子图像。
实施例1
执行步骤1，输入已配准的医学CT和MR-T2脑图像。
执行步骤2，对两幅图像进行预处理。
执行步骤3，设置迭代次数f＝3、梯度阈值α＝0.1，最大惩罚系数λ＝60和阈值函数斜率n＝11。在一定范围内，如α增加，越来越多的渐变被删除，图像性能变得特别模糊，如α达到一定水平，图像模糊的程度不会改变，如图5所示。应删除低于梯度阈值的梯度，同时保留大于梯度阈值的梯度。参数λ表示小梯度的去除程度。在一定范围内，最大惩罚系数越大，去除小梯度越明显。参数n控制接近梯度阈值的惩罚系数的下降速度。较大的n是，惩罚系数的下降速度越快，表明纹理和边缘之间的边界是清晰的，如图6所示。
执行步骤4～步骤5，对于梯度滤波分解的图像分别进行融合，设置了三组对比实验，三种方法分别为CNN、NSCT-PC、NSCT-PCNN，具体结果如图7中ours所示。
Claims (5)
1.基于梯度滤波的医学图像融合方法，其特征在于，具体按照以下步骤实施：
步骤1、输入两张已配准的多模态灰度医学图像；
步骤2、对于步骤1中的两张医学图像进行增强，获得细节和边缘信息更丰富的图像；
步骤3、设置阈值函数区分边缘和纹理的细节或噪声，并根据梯度阈值设计一个惩罚系数，根据对步骤2中增强的图像梯度滤波分解成为近似层和残差层；
步骤4、根据步骤3中分解得到的近似层和残差层分别设置融合规则，近似层图像采用一种能量属性的低频融合策略得到近似层图像，残差层采用一种参数自适应的脉冲耦合神经网络的高频融合策略得到残差层图像；
步骤5、将步骤4中近似层和残差层的融合子图像重建获得最终的融合图像。
2.根据权利要求1所述的基于梯度滤波的医学图像融合方法，其特征在于，所述步骤2具体按照以下步骤实施：
步骤2.1、首先对步骤1中的两张医学图像3×3像素为窗口内像素的局部差值信息D(i,j)：
其中D为窗口内像素的局部差值信息，(i,j)为图像像素横纵坐标，Z(x)表示每窗口内各点像素；
步骤2.2、计算图像的局部梯度信息：
其中g1(x)为图像的局部梯度信息，P(i,j)表示源图像的像素，θ和δ是g1的行和列，重新映射g1得：
Y2(x(θ,δ))＝ω2x(θ,δ)+ω2+1 (3)
其中Y2为重新映射g1的值，ω2表示增强强度参数，ω2越大表示图像整体增强水平越明显，x(θ,δ)表示Y2的自变量，
步骤2.3、得到加权系数g2：
增强的图像为：
其中H表示增强图像，⊙表示卷积操作，A(i,j)是步骤1中的两张医学图像。
3.根据权利要求2所述的基于梯度滤波的医学图像融合方法，其特征在于，所述步骤3具体按照以下步骤实施：
步骤3.1、构建边缘保留滤波器，优化以下函数：
其中，O*是需要优化的目标函数，I是输入图像，O是滤波后的输出图像，对于m×n的灰度图像，c表示像素坐标c＝(x,y)，
和
分别是x和y目标图像的方向，||·||2表示L2范数；
步骤3.2、设置一个阈值来区分边缘和纹理或噪声，惩罚系数设置为：
其中P为惩罚系数，λ是惩罚系数的最大值，g∈{x,y}表示x或者y方向，α是梯度阈值，n是阈值函数斜率，thresh定义为阈值函数：
步骤3.3、得到最后的目标函数为：
步骤3.4、为求解目标函数O*，将上式写为：
其中o和i分别是输出和输入图像的列向量形式，T是转置符号，Dx和Dy分别是在x和y方向上的Toeplitz矩阵，Px和Py是对角矩阵；
步骤3.5、使用加权最小二乘法WLS，对于变量的导数获得o并将导数设置为等于0，设置成为迭代形式并求解下面线性函数：
{E+Lt}ot+1＝i (12)其中E是单位矩阵，Lt是拉普拉斯矩阵，并且o0＝i，
t是迭代次数，
步骤3.6、根据上式设置滤波器为O＝GletF(I,t,α,λ,n)；
使用梯度滤波获得的近似层图像An表示如下：
An＝GletF(S,f,α,λ,k) (13)
其中An为近似层图像，S是步骤2中得到的增强图像，f是迭代次数，α是梯度阈值，λ是最大惩罚系数，k是阈值函数斜率，
步骤3.7、从增强图像中减去近似层图像以获得相应的残差层图像Rn：
Rn＝Sn-An (14)
其中n表示第n幅图像。
4.根据权利要求3所述的基于梯度滤波的医学图像融合方法，其特征在于，所述步骤4具体按照以下步骤实施：
步骤4.1、对于图像的残差层而言，使用参数自适应脉冲耦合神经网络PA-PCNN模型用于图像融合，该模型基于改进的简化PCNN模型，受生物视觉皮层模型的启发而产生，是由若干个神经元互联所构成的反馈型网络，其每一神经元由分支树、调制耦合器和脉冲产生器构成，分支树馈送的数学模型如下：
Fij[n]＝Rij (15)
步骤4.2、分支树链接输入的数学模型如下：
步骤4.3、调制耦合器的数学模型如下：
步骤4.4、PCNN的输出模块被称为脉冲产生器，脉冲产生器的变阈值函数如下：
步骤4.5、硬限幅函数的表达式如下：
其中，下标ij、kl为神经元的标号，Rij为输入图像的强度，Yij、Lij、Uij和θij分别为神经元的反馈输入、链接输入、内部活动项和动态阈值，Wijkl为根据八个相邻神经元的先前放电状态设置的链接权矩阵，VL、Vθ为幅度常数，αf、αθ为相应的衰减系数，β为链接强度，n为迭代次数，Yij为二值输出；
步骤4.6、使用PA-PCNN产生一次自适应计算获得新的参数值，具体的参数自适应计算公式如下：
其中σ(R)表示范围为[0,1]的输入图像R的标准偏差，R'和Rmax分别表示归一化的大津阈值和输入图像的最大强度；
步骤4.7、确定迭代次数n，结合公式(3)至(8)对PA-PCNN进行初始化，即将Uij[0]、Yij[0]和θij[0]值设为0，待初始化完成，第一次迭代中所有的非零强度神经元都将被激活，在此基础上依次完成n次迭代计算，最后统计n次迭代中每一个神经元的点火次数，每次点火时间描述为：
Tij[n]＝Tij[n-1]+Yij[n] (25)
其中Tij[n]表示n次点火时间，整个迭代完成后，神经元的总点火次数Tij[N]用于表征对应高频分量输入的活动水平，Tij[N]数值越大，则此神经元的活动水平越高，融合带的规则如下：
其中HF表征输入的融合系数，A和B分别为步骤3中分解得到的两幅残差层图像，HA,B与步骤3分解得到残差层图像的尺寸相同，新融合子带的融合系数选择以神经元点火总次数为依据，将两幅图的高频分量所对应同一位置像素点(i,j)的Tij[N]做比较，若TA,ij[N]＞TB,ij[N]，则新融合子带取HA(i,j)；反之，新融合子带取HB(i,j)，最终得到残差层融合图像FR；
步骤4.8、对于图像的近似层而言，采用一种能量属性的融合策略：
首先计算低频分量的固有属性值，定义固有属性值IPA,B的数学公式如下：
其中，IPA,B为低频分量的固有属性值，A,B为步骤三中对近似层分解得到的两幅图像，μ是低频分量的平均值；M是低频分量的中值，其次计算低频分量的能量属性函数，能量属性函数定义如下：
其中，EA,B为步骤3分解得到的两幅近似层图像的能量属性函数，a为调制参数，LA,B与步骤3分解得到近似层图像的尺寸相同，融合得到新的低频子带，融合后的低频子带通过能量属性函数和原低频分量加权平均得到，具体公式如下：
其中LF(i,j)即是新的低频融合子带系数，最终得到近似层融合图像FA。
5.根据权利要求4所述的基于梯度滤波的医学图像融合方法，其特征在于，所述步骤5具体按照以下步骤实施：
将步骤4中近似层和残差层的融合子图像重建获得最终的融合图像F：
F＝FA+FR (32)
其中F为最终融合图像，FA为近似层融合子图像，FR为残差层融合子图像。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210903636.1A CN115100093A (zh) | 2022-07-28 | 2022-07-28 | 一种基于梯度滤波的医学图像融合方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210903636.1A CN115100093A (zh) | 2022-07-28 | 2022-07-28 | 一种基于梯度滤波的医学图像融合方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN115100093A true CN115100093A (zh) | 2022-09-23 |
Family
ID=83301545
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210903636.1A Withdrawn CN115100093A (zh) | 2022-07-28 | 2022-07-28 | 一种基于梯度滤波的医学图像融合方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN115100093A (zh) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116843596A (zh) * | 2023-08-28 | 2023-10-03 | 浙江大学杭州国际科创中心 | X射线光栅多模态图像自适应融合方法、系统及装置 |
| CN116883803A (zh) * | 2023-09-07 | 2023-10-13 | 南京诺源医疗器械有限公司 | 一种用于胶质瘤边缘获取的图像融合方法及系统 |
| CN117197014A (zh) * | 2023-09-12 | 2023-12-08 | 南京诺源医疗器械有限公司 | 一种可降噪的肺部医学图像融合方法、系统和电子设备 |
-
2022
- 2022-07-28 CN CN202210903636.1A patent/CN115100093A/zh not_active Withdrawn
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116843596A (zh) * | 2023-08-28 | 2023-10-03 | 浙江大学杭州国际科创中心 | X射线光栅多模态图像自适应融合方法、系统及装置 |
| CN116843596B (zh) * | 2023-08-28 | 2023-11-14 | 浙江大学杭州国际科创中心 | X射线光栅多模态图像自适应融合方法、系统及装置 |
| CN116883803A (zh) * | 2023-09-07 | 2023-10-13 | 南京诺源医疗器械有限公司 | 一种用于胶质瘤边缘获取的图像融合方法及系统 |
| CN116883803B (zh) * | 2023-09-07 | 2023-12-05 | 南京诺源医疗器械有限公司 | 一种用于胶质瘤边缘获取的图像融合方法及系统 |
| CN117197014A (zh) * | 2023-09-12 | 2023-12-08 | 南京诺源医疗器械有限公司 | 一种可降噪的肺部医学图像融合方法、系统和电子设备 |
| CN117197014B (zh) * | 2023-09-12 | 2024-02-20 | 南京诺源医疗器械有限公司 | 一种可降噪的肺部医学图像融合方法、系统和电子设备 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Liu et al. | A medical image fusion method based on convolutional neural networks | |
| WO2022047625A1 (zh) | 一种图像处理方法、系统和计算机存储介质 | |
| Yang et al. | Multimodal medical image fusion based on fuzzy discrimination with structural patch decomposition | |
| CN115100093A (zh) | 一种基于梯度滤波的医学图像融合方法 | |
| CN108537794B (zh) | 医学图像数据处理方法、装置和计算机可读存储介质 | |
| CN110827216A (zh) | 图像去噪的多生成器生成对抗网络学习方法 | |
| Grampurohit et al. | Brain tumor detection using deep learning models | |
| CN109934887B (zh) | 一种基于改进的脉冲耦合神经网络的医学图像融合方法 | |
| Padmavathi et al. | A novel medical image fusion by combining TV-L1 decomposed textures based on adaptive weighting scheme | |
| CN109754403A (zh) | 一种ct图像内的肿瘤自动分割方法及系统 | |
| Du et al. | Three-layer image representation by an enhanced illumination-based image fusion method | |
| Angoth et al. | A novel wavelet based image fusion for brain tumor detection | |
| CN110503614A (zh) | 一种基于稀疏字典学习的磁共振图像去噪方法 | |
| CN110660063A (zh) | 多图像融合的肿瘤三维位置精准定位系统 | |
| Wang et al. | Functional and anatomical image fusion based on gradient enhanced decomposition model | |
| KR102142934B1 (ko) | 가중 최소 자승 필터와 희소 표현을 이용한 영상 융합 장치 및 방법 | |
| Yang et al. | Edge-preserving image filtering based on soft clustering | |
| Ullah et al. | Fast local Laplacian filtering based enhanced medical image fusion using parameter-adaptive PCNN and local features-based fuzzy weighted matrices | |
| CN116645283A (zh) | 基于自监督感知损失多尺度卷积神经网络的低剂量ct图像去噪方法 | |
| CN115457359A (zh) | 基于自适应对抗生成网络的pet-mri图像融合方法 | |
| Choubey et al. | A hybrid filtering technique in medical image denoising: Blending of neural network and fuzzy inference | |
| CN116894783A (zh) | 基于时变约束的对抗生成网络模型的金属伪影去除方法 | |
| CN116342444A (zh) | 一种双通道多模态图像融合方法及融合成像终端机 | |
| Kaur et al. | Multimodal Medical Image Fusion based on Gray Wolf Optimization and Hilbert Transform | |
| CN114387257A (zh) | 肺部影像中肺叶区域的分割方法、系统、装置和介质 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| WW01 | Invention patent application withdrawn after publication | ||
| WW01 | Invention patent application withdrawn after publication |
Application publication date: 20220923 |