三维图像的海森矩阵

标签： 三维图像 海森矩阵 二阶偏导数 高斯函数

---

海森矩阵（Hessian matrix）

---

雅可比矩阵

在向量分析中，雅可比矩阵是一阶偏导数以一定方式排列成的矩阵, 其行列式称为雅可比行列式。

海森矩阵

数学中，海森矩阵(Hessian matrix)是一个自变量为向量的实值函数的二阶偏导数组成的方块矩阵（假设其二阶偏导都存在）。
$$H(f)=\left[\begin{array}{cccc}\frac{{\partial }^{2}f}{\partial x_1^2}& \frac{{\partial }^{2}f}{\partial x_1{x}_2}& \cdots & \frac{{\partial }^{2}f}{\partial x_1{x}_n}\\ \frac{{\partial }^{2}f}{\partial x_2{x}_1}& \frac{{\partial }^{2}f}{\partial x_2^2}& \cdots & \frac{{\partial }^{2}f}{\partial x_2{x}_n}\\ ⋮& ⋮& \ddots & ⋮\\ \frac{{\partial }^{2}f}{\partial x_n{x}_1}& \frac{{\partial }^{2}f}{\partial x_n{x}_2}& \cdots & \frac{{\partial }^{2}f}{\partial x_n^2}\end{array}\right]$$

---

高斯求导

---

前言

通过上述公式可知，求海森矩阵的过程实际上就是求二阶偏导的过程。卷积中有一个重要的性质：卷积的微分特性—两个函数相卷积后的导数等于其中一个函数的导数与另一个函数的卷积
$$\frac{d}{dt}[f_1(t)\ast f_2(t)]=\frac{d{f}_1(t)}{dt}\ast f_2(t)=f_1(t)\ast \frac{d{f}_2(t)}{dt}$$
证明如下：
$$\frac{d}{dt}[f_1(t)\ast f_2(t)]=\frac{d}{dt}{\int }_{-\infty }^{\infty }{f}_1(\tau )\ast f_2(t-\tau )d\tau ={\int }_{-\infty }^{\infty }{f}_1(\tau )\ast \frac{d}{dt}{f}_2(t-\tau )d\tau =f_1(t)\ast \frac{d{f}_2(t)}{dt}$$
同理可证：
$$\frac{d}{dt}[f_1(t)\ast f_2(t)]=\frac{d}{dt}{\int }_{-\infty }^{\infty }{f}_1(\tau )\ast f_2(t-\tau )d\tau ={\int }_{-\infty }^{\infty }\frac{d}{dt}{f}_1(\tau )\ast f_2(t-\tau )d\tau =\frac{d{f}_1(t)}{dt}\ast f_2(t)$$
故上述微分特性得证。
推广易得：
$$\frac{{\partial }^2}{\partial t^2}[f_1(t)\ast f_2(t)]=\frac{{\partial }^2{f}_1(t)}{\partial t^2}\ast f_2(t)=f_1(t)\ast \frac{{\partial }^2{f}_2(t)}{\partial t^2}$$

推导

根据上式，令$f_1(t)$为高斯函数$G(x,y,z)=(\frac{1}{\sqrt{2\pi }\sigma }{)}^3{e}^{-\frac{x^2+y^2+z^2}{2{\sigma }^2}}$，$f_2(t)$为三维图片$I(x,y,z)$，则有:
$$\frac{{\partial }^2}{\partial x^2}[G(x,y,z)\ast I(x,y,z)]=\frac{{\partial }^{2}G(x,y,z)}{\partial x^2}\ast I(x,y,z)=G(x,y,z)\ast \frac{{\partial }^{2}I(x,y,z)}{\partial x^2}$$
可知，只要求得了$\frac{{\partial }^{2}G(x,y,z)}{\partial x^2}$,那么通过上式就可以得到$\frac{{\partial }^{2}I(x,y,z)}{\partial x^2}$为图像$I$在$(x,y,z)$点的对$x$的二阶偏导。其他方向的二阶偏导同理可求。
到这里，求图像的二阶偏导转换成了求高斯函数的二阶偏导。

---

跳过。。跳过。。一系列的求导过程((٩(//̀Д/́/)۶))

---

得到以下高斯函数的二阶偏导：
$$\frac{{\partial }^{2}G(x,y,z)}{\partial x^2}=\frac{1}{(\sqrt{2\pi }\sigma {)}^3}\frac{x^2-{\sigma }^2}{{\sigma }^4}{e}^{-\frac{x^2+y^2+z^2}{2{\sigma }^2}}=\frac{x^2-{\sigma }^2}{(\sqrt{2\pi }{)}^3{\sigma }^7}{e}^{-\frac{x^2+y^2+z^2}{2{\sigma }^2}}$$
$$\frac{{\partial }^{2}G(x,y,z)}{\partial x\partial y}=\frac{1}{(\sqrt{2\pi }\sigma {)}^3}\frac{xy}{{\sigma }^4}{e}^{-\frac{x^2+y^2+z^2}{2{\sigma }^2}}=\frac{xy}{(\sqrt{2\pi }{)}^3{\sigma }^7}{e}^{-\frac{x^2+y^2+z^2}{2{\sigma }^2}}$$
同理易得高斯函数$G(x,y,z)$对$y^2,z^2,yx,xz,zx,yz,zy$等方向的偏导。
可以发现：
$$\frac{{\partial }^{2}G(x,y,z)}{\partial x\partial y}=\frac{{\partial }^{2}G(x,y,z)}{\partial y\partial x}$$
$$\frac{{\partial }^{2}G(x,y,z)}{\partial x\partial z}=\frac{{\partial }^{2}G(x,y,z)}{\partial z\partial x}$$
$$\frac{{\partial }^{2}G(x,y,z)}{\partial y\partial z}=\frac{{\partial }^{2}G(x,y,z)}{\partial z\partial y}$$
至此，高斯函数所有的二阶偏导已经求得，然后利用matlab中的convn函数进行三维空间内的卷积(参数选择same保证结果和图像一致)，这也意味着黑森矩阵已经可以通过上述过程得到。

代码实现

%%  求高斯函数的二阶偏导数
%%  num为高斯核的大小
%%  sigma为高斯函数的方差
function [gau_xx,gau_yy,gau_zz,gau_xy,gau_xz,gau_yz]=gaus_creation_3D(num,sigma)

gau_xx=[];gau_yy=[];gau_zz=[];%初始化矩阵

gau_xy=[];gau_xz=[];gau_yz=[];%初始化矩阵

for i=1:1:2*num+1
    for j=1:1:2*num+1
        for k=1:1:2*num+1
            
            x=i-num-1;y=j-num-1;z=k-num-1;
            
            gau_xx(i,j,k)=1/power(sqrt(2*pi),3)*(-(sigma^2-x^2)/sigma^7)*exp(-(x^2+y^2+z^2)/2/sigma^2);
            gau_yy(i,j,k)=1/power(sqrt(2*pi),3)*(-(sigma^2-y^2)/sigma^7)*exp(-(x^2+y^2+z^2)/2/sigma^2);
            gau_zz(i,j,k)=1/power(sqrt(2*pi),3)*(-(sigma^2-z^2)/sigma^7)*exp(-(x^2+y^2+z^2)/2/sigma^2);
            
            gau_xy(i,j,k)=1/power(sqrt(2*pi),3)*(x*y/sigma^7)*exp(-(x^2+y^2+z^2)/2/sigma^2);
            gau_xz(i,j,k)=1/power(sqrt(2*pi),3)*(x*z/sigma^7)*exp(-(x^2+y^2+z^2)/2/sigma^2);
            gau_yz(i,j,k)=1/power(sqrt(2*pi),3)*(z*y/sigma^7)*exp(-(x^2+y^2+z^2)/2/sigma^2);
        end
    end
end
end

思考

关于最终的结果

从$\frac{{\partial }^2}{\partial x^2}[G(x,y,z)\ast I(x,y,z)]=\frac{{\partial }^{2}G(x,y,z)}{\partial x^2}\ast I(x,y,z)=G(x,y,z)\ast \frac{{\partial }^{2}I(x,y,z)}{\partial x^2}$式中可知最后的结果其实是图像的二阶偏导和高斯函数的卷积，并不只是单纯的图像二阶偏导。高斯函数在图像处理中常用于去除高斯噪声，它具有良好的低通滤波效果，一般在检测边缘之前常用高斯卷积来移除图像一些细节以及噪声。所以事实上，这里的卷积不会影响图像整体的结构，而且一定程度上对图像进行了去噪使图像质量更好(当然不可否认的是损失了一些图像细节)，如果是利用海森矩阵进行三维图像的线性结构或是面结构的检测，那么一定程度的去噪以及平滑处理将可能得到更好的结果。

---

以上。(づ●─●)づ