第6章医学图像配准(6.3-6.6)

3x3反向模板
第一次扫描： 用正向模板，沿图象从左向右，自上向下遍历，在模板 每一位臵，模板元素与图象上对应象素值相加，得到一组5个临时和。 取其最小值作为该象素之新灰度值，模板空的元素什么都不作。 第二次扫描： 模板从右向左，自底向上，运算方法相同，两次扫描 后结果即为图象的距离变换灰度图。 此3x3模板产生的距离变换灰度值为笛卡儿距离的3倍。经规范化，与 真实笛卡儿距离最大偏差8%。
典型优化算法： • 单纯形法（Simplex） 特点是容易实现，能够给出良好特性。该方法在M维空间中移动具有M+1 个点的单纯形，M是参数的个数。对单纯形做反射、扩展或收缩去捕获最 小值。 • Powell法
6.5 基于最大互信息的 多模医学图像配准
互信息是信息理论中的一个基本概念， 通常用于描述两个系统间的统计相关性，或 者是在一个系统中包含的另一个系统的信息 的多少。 在多模医学图像配准问题中，虽然两幅 图像来源于不同的成像设备，但是它们基于 共同的人体解剖信息，所以当两幅图像的空 间位臵完全一致时，其中一幅图像表达的关 于另一幅图像的信息，也就是对应像素灰度 的互信息应为最大。
式中S是相似性测度，T是空间变换。 配准的过程可归结为寻求最佳空间变换：
T arg max S (T )
T

由于空间变换包含多个参数，这是一个多参数最 优化问题， 一般由迭代过程实现。
T ' T T
增量T可以用计算梯度的方法计算
dS (T ) T dT
6.5.2 以互信息为相似性测度
1 n f (c, R, t ) yi (cRxi t ) n i 1
2
两幅图像的配准可以通过函数f对参数c，R和t的最小化求解。

xy
1 n ( yi y )(xi x ) T n i 1
令Σ的奇异值分解（SVD）为 UDVT，这里U和V是归一正交矩阵，D是对角矩阵。 则
如何选取部分采样点也有两种方法。一个是 随机采样法；另一个是从原有数据格点中有规 律的每隔几点选取一点。采样子集技术固然能 几十倍、甚至几百倍的加快计算速度。但是本 来用直方图近似概率分布的误差会因过少的数 据点，使误差进一步扩大。
实验表明，对于CT图像在XY平面上的分辨率为 512*512，在Z轴方向采样28-34层， MR图像在XY 平面上的分辨率为256*256，在Z轴方向采样20-26 层的图像，采样因子fx= fy=4, fz=2会取得满意的 配准结果。
6.3.3 基于极值点的刚体配准
对同一个病人用相同成像模式但不同位臵的两幅图像进行配准， 可以从这两幅图像上分别提取极值线及极值点, 再设法寻找这两组极值 点间的对应关系。最后，计算能使这些匹配线段重叠的最佳刚体变换。 一种典型的技术是采用几何散列（Hashing）表与Hough变换相结合的 方法。
倒角匹配算法的关键部分是距离变换（Distance Transformation）。距 离变换的类型决定了使用该变换后代价函数描述两幅图像特征间距离的准 确程度。 距离变换是将二值图转换为灰度图像的一种变换。每个象素的灰度值 是该象素与最近背景象素间的距离。常用的方法为倒角（Chamfer）算法。 倒角算法通过对图像两次扫描，实现距离变换。 1. 一维距离变换 对图像的一列像素做距离变换。设二值图像像素0：背景， 1：物体。
用联合概率分布 p AB (a, b) 和完全独立时的概率分布 p A (a) p B (b) 间的广义距离来估计互信息：
I ( A, B) p AB (a, b) log
a ,b
p AB (a, b) p A (a) p B (b)
对于离散的数字图像，联合概率分布
p AB (a, b)
前瞻性和回顾性医学图像配准
医学图像配准方法可分为前瞻性和回顾 性两种。 如果在成像前，在病人的颅骨附近固定 标志点或定位框架，成像后通过对齐标志点 使两幅图像配准，称为前瞻性配准算法。
前瞻性配准方法一般说来配准精度较高， 但它是有创伤的，而且操作比较复杂。
前瞻性和回顾性医学图像配准
如果成像时不采取任何措施，在成像 后应用某种算法寻找两幅图像的空间对应 关系，则称为回顾性配准。 例如特征点（面）法、不变矩配准法、 相关配准法等。 特征点（面）法原理简单，应用面广， 但是需要较多的人工介入，配准精度受特 征点（面）提取精度的限制。
平均代价函数:
1 C mean (T ) N
平方根代价函数:
i
F (T , r )
i 1
N
C RMS (T )
1 N 2 F (T , ri ) N i 1
上述方程是倒角匹配算法的核心。关于F，并不直接使用图像的灰度 值，而是先对输入图像分割，然后再对分割的结果做距离变换得到的。 在距离变换中，像素值给出到分割特征的最短距离。所以代价函数 是 两幅分割的图像特征间平均距离的精确估算。对此函数的最小化可以 实现两幅图像的配准。
互信息只是一个相似性测度，如何准 确计算互信息和实现不同模式图像间互 信息的最大化，尚存在许多重要的技术 问题需要解决。 实施最大互信息配准法的一些重要技 术问题,包括 不增加新数据点的格点采样子集 不产生分数灰度值的PV插值技术 出界点策略 在搜索策略上采用无须计算梯度的 Powell算法
• • • •
6.4.2 图像分割与代价函数的优化
采用何种分割方法取决于实际应用。倒角匹配算法只要求将图像F分 割成二值图像作为距离变换的输入。还要从分割得到的图像G产生一组轮 廓点。轮廓点可用简单的轮廓跟踪法或对二值图像扫描搜索非零元素（线 条图中点的顺序无关紧要）。倒角匹配算法对所用的分割算法并无特殊的 约束。从效率因素考虑总是希望线条中的点数目相对少些为好。 倒角匹配算法对分割质量具高度强健性，这意味可以使用“低质量”自动 分割技术。即算法可以是全自动的，即使图像质量差些也无妨。
6.3.4 仅依赖基准点位臵的刚体配准
在许多情况下可以通过人为地在目标表面加记号，或者用交互式方法或自 动识别图像中明显的解剖标志点、几何图形明显转折点或图像灰度阶跃变化点 等，统称基准点。如果这类基准点数目足够多的话，可以用来引导图像配准。 其特点是仅依靠基准点的空间坐标而无须计算那些复杂的微分几何特性。
6.5.1 配准原理
对于在不同时间或/和不同条件下获取的 两幅图像A(X) 和B(X) 配准， 就是要定义一个相似性测度并寻找一个空 间变换关系，使得经过该空间变换后两幅 图像间的相似性达到最大。即使图像A上的 每一个点在图像B上都有唯一的点与之相对 应。并且这两点应对应同一解剖位臵。
S (T ) S ( A( X ), B(T ( X )))
6.3.1 极值线与极值点
极值线 (Extremal Lines)
极值函数e 是 曲率k在t方向的方向导数。
e1 k1 t1 e2 k 2 t 2
局部绝对值最大的曲率点是极值函数e的零交叉点。 极值线定义为e1 或e2 的零交叉点轨迹。脊线是正最大曲率极大 及负最大曲率极小的点轨迹。
第一次：左->右扫描，背景点保留0， 物体点为到最近0的距离。 第二次：右->左扫描，物体点取值=min(第一次结果，物体点到最近0的距 离) 得到的结果是用灰度表示的反映该像素与最近边界距离的图像。
2. 二维倒角变换
将模板与图像卷积实现二值图像到灰度图像的转换。
4
3
3
0
4
0 4 3 3 4
3x3正向模板
6.3 基于基准点的配准方法
三维医学图像的配准算法可以基于从图像提取的特征，或基 于图像强度的相似性最优化。首先从三维医学图像提取一个表面， 然后在该表面上提取基准线，可能的话再从这些基准线上选取基 准点。 脊线是表面上的轨迹（或轨迹的集合），该处的主曲率在相 应的主方向达到局部极大。脊线的提取要通过计算图像强度的三 阶导数完成。涉及微分几何的知识。
互信息配准的关键技术问题
1. 采样子集 (Sub-sampling) 关于待配准的两幅图像（可以是二维 或三维的）的定义域是整幅图像大小。对 于两幅图像的联合直方图及每幅图像的独 立直方图的计算应使用图像的全部数据点 灰度值。
理论上，也可以使用任一幅图像的子集或超集 （Super-sampling）。但是，对于大多数断层 扫描医学图像来说，3D体积数据集所包含的数 百万、乃至上千万的数据点极大地增加了计算 的负担，就目前的计算机速度来说，会使3D体 积图像的配准无法实时实现和临床应用。因此， 实用中都是只采用部分数据点计算图像的直方 图，称做采样子集。
极值点(Extremal Points)
f I,
极值点定义为三个隐含数表面的交点，这些表面是
e1 0 , e2 0
其中，f 是 等值面 。这些极值点对刚体变换具有不变性。
6.3.2 极值点的自动提取方法
(a)
(b)
(c)
第1步：提取等值表面：在小立方体元的内部提取等值表面。该等值表面（f=I）与 小立方体元的三个边相交（图a） 。 第2步：提取脊线：用小立方体元角点 e1值线性内插计算得到等值表面三角形顶点 的e1值。如果这些e1值符号相同，则在此小立方体元内没有 的极值线。否则，在e1 符号相反的小立方体元的两条边上寻找插值后e1 =0的点{ P1,P2}。联接 P1,P2的线段 就是近似的极值线（图b）。 第3步：计算极值点的位臵：由于P1,P2 两点都在小立方体元的表面上，可以用表面 内e2 值做双线性内插得到P1,P2 点上的 e2值。如果这两个e2 值符号相同，则在此小 立方体元内没有极值点存在。否则，可以沿线段{P1,P2 }对 e2值做内插。插值后该 线段上 e2 =0的点P就是所寻求的极值点（图c）。
R USV
T
t y cR x
c
1

2 x
tr ( DS )
6.4 倒角匹配图像配准法
倒角匹配（Chamfer Matching）是一种对分割特征进行匹配的经典图 像配准方法。后来在从卫星图像中寻找湖泊的轮廓方面得到有效的应 用。该方法的基本思想是将线条图与图像匹配。
6.4.1 代价函数与距离变换
a ,b
其中，a A, b B 如果
H ( A | B) 表示已知系统B时A的条件熵，那么
H ( A) 与 H ( A | B)
的差值，就代表了在系统B中包含的A的信息，即互信息。 因此两个系统间的互信息可以用下式描述：
I ( A, B) H ( A) H ( B) H ( A, B) H ( A) H ( A | B) H ( B) H ( B | A)
距离变换示例：
二值的哑铃形状图像经过距离变换转化成灰度图像。 数值的大小代表该像素到物体边界的距离。左图为原始 图象，右图为使用3x3的模板得到的距离变换图象
距离变换示例： 二值图像转换为灰度图像
左图给出带有噪声的二值图像的 距离变换。 原始图像是包含较多噪声干扰的 局部重叠的两个圆盘状对象。 图（b）和（c）分别是使用3x3 模板和5x5模板得到的距离变换 结果。由于原图噪声过多，所以 直接进行距离变换的效果不理想， 对原图做闭运算处理，起到了平 滑和去噪的作用，得到可以接受 的距离变换图象(d)。
互信息是信息理论中的一个基本概念，通 常用于描述两个系统间的统计相关性， 或者是在一个系统中包含的另一个系统的信 息的多少，可以用熵来描述。熵表达的是一 个系统的复杂性或者是不确定性。 系统A的熵定义为：
H ( A) p A (a) log p A (a)
a
两个系统的联合熵为：
H ( A, B) p AB (a, b) log p AB (a, b)
可以用归一化的联合直方图表示：
h(i, j ) p AB (i, j ) i, j h(i, j)
边缘概率分布 p A (a) 表示为： p A (i)
p
j
AB
(i, j )
边缘概率分布 PB (b) 表示为： p B ( j )
p
i
AB
(百度文库, j )
p AB (i, j ) I ( A, B) p AB (i, j ) log p A (i) p B ( j ) i, j