第5次课 医学图像恢复技术

1、医医 学学 图图 像像 处处 理理Medical Image Processing陈家益陈家益()广东医学院 信息工程学院 计算机科学教研室主要内容：主要内容：5.4 医学图像失真的校正医学图像失真的校正5.2 医学图像中的噪声及滤除医学图像中的噪声及滤除5.3 系统退化模型系统退化模型5.1 图像恢复的基本概念图像恢复的基本概念第五章 医学图像恢复技术 在成像过程中受到各种非完美系统因素影响和人体本身的限制，医学图像会存在诸多伪影和失真，亟待进行图像的恢复。对图像丢失的信息知道的越多，越容易恢复到清晰的图像。 图像恢复技术与图像增强技术的区别： 目的都是改善图像，只是改善的方法途径不一样。

2、图像增强不考虑图像质量下降的原因，只是增强感兴趣的区域，是主观的过程。 图像恢复要充分考虑图像质量下降的原因，过程是客观的。（1）图像退化的概念 图像在形成、传输和记录的过程中，由于受到各种因素的影响，产生模糊、失真或出现噪声，致使图像质量下降，无法完全反映目标的真实信息。（1）图像退化的概念图像退化原因： 成像设备与物体的相对运动（聚焦不准）、传感器特性的非线性、感光胶卷的非线性和胶片颗粒噪声等。典型的退化类型： 噪声污染、边缘模糊、运动模糊。（2）图像恢复 根据图像退化的原因，建立退化过程的数学模型，并采用相反的过程（逆过程）进行处理，以改善退化图像的质量并尽可能恢复原始图像的过程。（2）

3、图像恢复 根据图像退化因素的不同，采用不同的恢复方法。 根据方法所在域区分，分为空域、频域或小波域恢复； 根据方法所用恢复模型的约束条件区分，分为无约束恢复和有约束恢复； 根据是否需要外来干预，分为自动恢复和交互式恢复。（3）图像退化/恢复过程的模型 图像恢复可以表示为建立退化模型并推导恢复模型的过程。图像退化/恢复过程的模型表示为：退化函数h(x,y)代表的是线性、移不变成像系统，n(x,y)是噪声。( , )( , )*( , )( , )g x yh x yf x yn x yf(x,y)退化函数退化函数H H+g(x,y)n(x,y)恢复滤波器恢复滤波器f(x,y)的的最优估计最优估计

4、退化过程退化过程恢复过程恢复过程 噪声直接影响图像密度分辨率和空间分辨率，掩盖了图像某些重要的信息，给图像分析识别和诊断带来不利影响。滤除噪声是医学图像恢复技术一个重要的分支。（1）噪声模型 噪声分为多种类型，其统计特性可由概率密度函数（pdf）描述。滤除噪声关键是确定图像中噪声的类型和相应的统计参数。（1）噪声模型高斯噪声 高斯噪声的概率密度服从高斯分布（正态分布）。设变量z服从高斯分布，高斯噪声的PDF为：Z为图像的灰度值，u为 z 的均值或期望值， 是z的标准差， 是z的方差。22() /21( )2z up ze 2 （1）噪声模型高斯噪声+（1）噪声模型高斯噪声 g=imnoise(

5、f,gaussian,m,var) 将均值m,方差为var的高斯噪声加到图像f上，默认值为均值是0，方差是0.01的噪声。（1）噪声模型高斯噪声=02=0.0025=02=0.01=0=0.0025=0.2=0.0025（1）噪声模型伽马噪声 伽马噪声常见于光子成像的探测过程，服从伽马分布。 其中a0，b是正整数，当参数确定时，伽马分布的概率密度函数的均值和方差分别为：和 。 10(z)(1)!00bbaza zezpbz ba 22ba （1）噪声模型伽马噪声 1(z)(1;5)(1)!bbaza zpeabb （1）噪声模型指数噪声 当伽马噪声的概率密度函数中b=1就成为指数分布，指数分布

6、主要出现在激光成像系统中。指数噪声服从指数分布。 指数概率密度的均值和方差分别为： 和 0(z)00azaezpz 1a 221a （1）噪声模型指数噪声 a(z)(a1)azpae （1）噪声模型指数噪声f=imread(gray.jpg);m,n=size(f); a=0.1;noise= (-1/a)*log(1-rand(m,n);g=double(f)+noise;imshow(uint8(g); （1）噪声模型瑞利噪声 瑞利分布用来描述光子的散射特性或平坦衰减信号。瑞利噪声服从瑞利分布。 瑞利分布的概率密度的均值和方差分别为： 和 2(z a) /2(z a)(z)0bezapbz

7、a / 4ab2(4)4b （1）噪声模型瑞利噪声 2(z a) /2(z)(z a)(a1; b5)bpeb（1）噪声模型瑞利噪声 g=raylrnd(B,M,N) 生成大小为M*N的服从B为参数的瑞利分布的随机数,此处B应为瑞利分布的均值u(x)=B*sqrt（pi/2）中的B。f=imread(gray.jpg);g=raylrnd(20,size(f);rayl=double(f)+g;imshow(uint8(rayl); （1）噪声模型脉冲噪声 脉冲噪声是由非连续的、持续时间段和幅度大的不规则脉冲或噪声尖峰组成的分布，脉冲噪声可以是正的，也可以是负的。脉冲噪声通常由电磁干扰和通信系

8、统传输故障等引起。脉冲噪声服从脉冲分布。 (z)0elseabPzapPzb （1）噪声模型脉冲噪声 如果ba，则灰度b显示为亮点；反则b显示为一个暗点。如果Pa或Pb为0，称为单极脉冲噪声；如果Pa或Pb都不为0，称为双极脉冲噪声，又称椒盐噪声。 (z)0elseabPzapPzb （1）噪声模型脉冲噪声 （1）噪声模型脉冲噪声 g=imnoise(f,salt&pepper,d) 用椒盐噪声污染图像f，其中d是噪声密度（即包括噪声值的图像区域的百分比）。因此，大约有d*numel（f）个像素受到影响。默认的噪声密度为0.05。 （1）噪声模型均匀噪声均匀噪声服从均匀分布。均匀分布概

9、率密度函数的均值和方差分别是 和 1(z)0azbpbaelse 2ab 22()12ba （1）噪声模型均匀噪声 （1）噪声模型均匀噪声f=imread(gray.jpg);a=10;b=50;noise=a+(b-a)*rand(m,n);g=uint8(double(f)+noise);imshow(g); （2）噪声模型的参数估计 参数的估计就是求均值 u 和方差 值（如高斯噪声），噪声密度d（如椒盐噪声），参数a、b的值（如指数、伽马噪声）。 有些参数可以从成像传感器的性能指标中得到；如果不知道这些参数，但对成像系统的噪声特性有先验知识或者可以获得噪声图像，考虑到不同类型噪声的均值、

10、方差与参数之间的关系，可以根据模型由均值和方差得到相应的模型参数。 2 （2）噪声模型的参数估计 但如果无法获得成像系统的先验知识，可以用下面的方法推断或估计噪声的类型和参数。 利用带噪图像中的均匀区域（背景近似为常数）的直方图对噪声的类型和参数进行估计。假设S为一均匀区域，利用下面两个式子对该区域的均值和方差进行估计： 22=()=( - )()iiiiiizSzSz p zzp z（2）噪声模型的参数估计 直方图可以看做是灰度PDF分布的近似或估计。均匀区域的灰度相近，因此区域内的灰度变化主要是由噪声引起，通过直方图的形状可以得到与噪声分布相近的类型。 如果直方图的形状近似于高斯分布，通过

11、 和 得到均值和方差就是所需要的参数；如果是其他的噪声类型，则可以通过均值和方差计算出相应的参数a和b的值。对于脉冲噪声，最好选择图像灰度中等的均匀区域。22=()=( - )()iiiiiizSzSz p zzp z00.20.40.60.81051015（2）噪声模型的参数估计f=imread(noise.jpg);B=roipoly(f);g=uint8(zeros(size(f); %或 g=f; g(:)=0;g(B)=f(B);imshow(g); %或imshow(f.*uint8(B);p=imhist(f(B); p=p./sum(B(:); z=0:255; z=z./25

12、5; u=sum(z.*p); v=sum(z-u).2).*p); 22=()=( - )()iiiiiizSzSz p zzp z（3）噪声图像的恢复 当图像唯一存在的退化是噪声时，退化的图像就成为原图像与噪声的相加。噪声在数学上主要表述为加性噪声和乘性噪声。对于乘性噪声，经对数处理后，可以转化为加性噪声。 ( , )( , )*( , )( , )g x yh x yf x yn x y( , )( , )( , )g x yf x yn x y（3）噪声图像的恢复1）平滑滤波器滤除噪声 常用的滤波方法包括空域滤波和频域滤波。 空域滤波主要采用各种图像平滑模板对图像进行处理，以滤除噪声，

13、包括均值滤波和中值滤波。 频域滤波先对图像进行傅里叶变换到频率域，然后用适当的频率带通滤波器进行处理，再变换回空间域的图像。 （3）噪声图像的恢复1）平滑滤波器滤除噪声 均值滤波适合处理高斯或均匀分布的随机噪声。但均值滤波对所有的像素点都同等对待，在对噪声平均的同时，也将图像的边界点和细节平均了。因此在滤除噪声的同时，会使图像变模糊。滤波器越大，滤除噪声效果越好，但是图像越模糊。（3）噪声图像的恢复1）平滑滤波器滤除噪声 中值滤波属于顺序统计滤波器，根据滤波器区域内像素大小的排列顺序来确定滤波器的输出。因此适合处理单极或双极脉冲噪声。 脉冲噪声是叠加在图像上的非线性噪声，呈现出较高的灰度值（盐

14、噪声）和/或较低的灰度值（胡椒噪声），因此选取邻域内位于中间的灰度值，去除过高或过低的影响，可有效滤除脉冲噪声。中值滤波引起的模糊很小，能更好地保留图像的细节信息。（3）噪声图像的恢复2）均值滤波器：算术平均：几何平均：调和均值：反调和均值：( , )1( , )( , )xys tSg x yf s tmn 1( , )( , )( , )xymns tSg x yf s t ( , )( , )1( , )xys tSmng x yf s t 1( , )( , )( , )( , )( , )xyxyQs tSQs tSf s tg x yf s t （3）噪声图像的恢复2）均值滤波器：

15、算术平均：平滑图像的局部变化，在减少噪声的同 时，带来模糊。几何平均：平滑程度与算术均值滤波器相当，但丢 失的细节更少。调和均值：对“盐”噪声滤波效果好，但不适合 “胡椒”噪声，适于处理类似高斯的噪 声。反调和均值：Q为正时，适于处理“胡椒”噪声； Q为负时，适于处理“盐”噪声；（3）噪声图像的恢复2）均值滤波器：算术平均：h=fspecial(average,m n);g=imfilter(f,h,replicate);（3）噪声图像的恢复2）均值滤波器：几何平均：f=im2double(f);g=exp(imfilter(log(f),ones(m,n),replicate).(1/m/n);g=im2uint8(g);含高斯噪声的退化图像含高斯噪声的退化图像算术均值滤波后图像算术均值滤波后图像几何均值滤波后图像几何均值滤波后图像（3）噪声图像的恢复2）均值滤波器：调和均值：f=im2double(f);g=m*n./imfilter(1./(f+eps),ones(m,n),replicate);g=im2uint8(g);（3）噪声图像的恢复2）均值滤波器：反调和均值：f=im2double(f); Q=1.5;g1=imfilter(f.(Q+1),ones(m,n),replicate);g2=imfilter(f.Q,ones(m,