加权稀疏表示

加权稀疏表示（共7篇）

加权稀疏表示 篇1

摘要：针对低分辨率、低质量人脸图像重建问题,提出了一种新的基于稀疏表示的人脸超分辨率算法。在训练阶段,人脸的位置特征被用于保持人脸块的全局信息,人脸块间的几何结构被用于保持高低分辨率超完备冗余字典的流形结构,从而提高字典的表达能力;在重建阶段,K近邻加权稀疏表示被用于消除稀疏编码噪声,以提高高分辨率人脸图像重建系数的精度。实验结果表明,提出的方法取得了较好的主客观质量。

关键词：图约束字典,加权稀疏表示,人脸超分辨率

视频监控具有实时、可视、可记录及信息量大等特点,在安防、交通、生产管理等各应用领域发挥重大作用。监控视频中的人脸图像是辨识关键人物最直接的线索。但是在实际监控应用中,由于摄像头和人脸距离通常较远,视频监控系统的带宽和存储资源有限等原因,往往导致监控视频中人脸图像的分辨率较低,难以达到辨识要求。人脸超分辨率技术也叫做人脸幻象技术( Face Hallucination)[1],它能够在不改变硬件环境的情况下,从一幅或多幅低分辨率输入人脸图像中,重建出一幅高分辨率的人脸图像,达到改善人脸图像的清晰度的目的。该技术在安防监控、计算机视觉等领域中具有重要的应用。

自Baker和Ce Liu等人提出“人脸幻象”的概念之后,人脸超分辨率的研究就受到了广泛的关注,出现了一系列的人脸超分辨率算法,如基于子空间的方法[2,3,4]、基于流形的方法[5,6]和基于稀疏表示的方法[7,8,9]等。本研究中主要关注基于稀疏表示的人脸超分辨率方法。

稀疏表示( Sparse Representation) 被认为是一种新型、有效而且鲁棒的特征表示方式,成功地应用于一系列实际问题中,是目前最热门的研究课题之一[10]。Yang等首次将压缩感知的思想应用到超分辨率领域,提出了基于稀疏表示的自然图像超分辨率算法[7]。该方法充分利用稀疏表示符合人眼感知特性的优势,自动选取近邻个数,避免了人为指定近邻块个数会导致合成信息不足或过拟合的弊端,取得了较好的超分辨率效果。在Yang的基础上,Chang等提出基于稀疏表示和双字典学习的人脸素描图像合成方法[8],Jung等提出了基于稀疏表示和原始样本块字典的人脸超分辨率算法[9]。

尽管基于稀疏表示的人脸超分辨率方法比其他方法取得了更好的效果,但是,已有的基于稀疏表示的方法较少考虑样本数据内在的几何结构,这影响了冗余字典的表达能力。最近研究表明,保持数据内部的局部拓扑结构能有效改善稀疏编码的效果[11,12,13]。因此,本文提出了一种基于图约束字典和加权稀疏表示的人脸超分辨率算法,算法包括训练和重建两个阶段。在训练阶段,三方面措施被用于提高冗余字典的表达能力: 一是利用人脸的位置信息对样本进行聚类,保持人脸块的全局信息; 二是利用图约束稀疏编码进行字典的学习,保持人脸块的局部信息; 三是利用联合训练高低分辨率冗余字典,保持高低分辨率字典对在内部流形结构上的一致性。在重建阶段,加权稀疏编码均值被用于消除超分辨率重建过程中的稀疏编码噪声,以提高高分辨率人脸图像重建系数的精度。实验结果证明了该方法的有效性。

1 图约束稀疏编码

1. 1 图的构建

谱图理论[14]和流形学习理论[13]的研究表明,散布点的局部几何结构可以通过近邻图进行有效地建模。这里简要介绍近邻图的构建过程。

假设X = [x1,x2,…,xN]∈ RM ×Q是一个数据集,M表示数据维数,Q表示数据集中元素的个数。基于数据集X构建的图表示为G = ( V,W) ,其中,V = Q表示该图中顶点的个数,每个数据点作为图的一个顶点。每两个顶点之间用无向边wi,j连接,wi,j表示顶点i和顶点j之间的相似性,则数据集X的权重矩阵表示为W = [wi,j]Q ×Q。计算顶点之间权重的方法主要有三种,分别是: 0 - 1 权重法、热核权重法和点积权重法[12]。由于热核权重法更适合于图像处理,这里采用热核权重法计算顶点之间的相似性权重。欧几里德距离空间的热核权重表示为

式中: h为常数,h >0,其值取决于期望的稀疏编码系数权重的分布幅度,期望的稀疏编码系数权重的分布幅度大,则h取值大; ci为归一化。权重矩阵W表示为

式中:Nk( xi) 表示数据点xi的K近邻数据集合。权重矩阵W包含了整个数据空间的几何结构信息。基于以上图模型,将每个人脸图像块作为一个数据点,从而可以实现训练样本库的几何图的构建。

1. 2 图约束稀疏编码

假设X = [x1,x2,…,xN]∈ RM ×Q表示原始数据矩阵,D = [d1,d2,…,dL]∈ RM ×L表示字典矩阵,S = [s1,s2,…,sN]∈ RL ×Q表示数据集X在字典D上的稀疏编码系数矩阵。对数据集X进行稀疏编码的目的在于找到一个字典矩阵D和一个稀疏系数矩阵S,使得二者的乘积能够尽可能地逼近原始数据矩阵。稀疏编码的目标函数表示为

式中: 第一项为逼近误差项; 第二项为稀疏性约束项; λ 为正则化参数,用于平衡逼近误差和稀疏性之间的权重。

为了实现原始数据集内部的几何图结构向稀疏编码系数空间的映射,通常假设在欧几里德距离空间近邻的两个原始数据点xi和xj,在稀疏编码系数空间的表示系数si和sj也是近邻的。该假设得到了大量实践的验证[11,12,13]。因此,K近邻图约束可以表示为

式中: I为单元矩阵; L = ( I - W) ( I - W)T。将K近邻图约束结合到原始的稀疏编码框架中,新的优化问题可以表示为

式中: 常数 β 用于平衡K近邻图约束项在优化过程中的比重。式( 5) 同时拥有两个变量D和S,是一个非凸优化问题,无法正常求解。对上式通常采用迭代方法求解,即先固定一个变量,求解另一个变量。

2 基于图约束字典和加权稀疏表示的人脸超分辨率算法

2. 1 图约束位置块字典的学习

由于人脸较之自然图像更具有规则性。人脸图像分块的位置,如眼睛、鼻子、嘴等,包含了人脸图像重建的先验信息。因此,在构建学习字典的过程中,首先根据位置信息,对人脸图像分块进行聚类。同一位置上的分块作为一个训练集。

假设IH={IH,q}Qq=1=[IH,1,IH,2,…,IH,Q]∈RM×Q和IL={IL,q}Qq=1=[IL,1,IL,2,…,IL,Q]∈RN×Q分别表示对应的高低分辨率训练图像集矩阵,矩阵的列数Q表示训练图像数量,矩阵的行数M、N分别表示高低分辨率图像向量的维数,其中M=s2N,s表示下采样的倍数。训练集中的每幅图像分割成P个小块。根据人脸图像的位置信息,所有的训练块可以分为P个集合。高分辨率训练块集合表示为{I1H,q}Qq=1,{I2H,q}Qq=1,…,{IPH,q}Qq=1,低分辨率训练块集合表示为{I1L,q}Qq=1,{I2L,q}Qq=1,…,{IPL,q}Qq=1。由于人脸图像事先进行了对齐操作,相同位置块集合能提供更多的重建信息。

对于给定的位置P,对应的高、低分辨率块集合IHp={ IpH,q| 1 ≤ q ≤ Q} 和ILp= { IpL,q| 1 ≤ q ≤ Q} 被用于训练高、低分辨率字典DHp和DLp。基于图约束位置块的高、低分辨率字典分别可以通过式( 6) 、式( 7) 得到

为了保持高、低分辨率字典DHp和DLp中高、低分辨率图像块之间的映射关系,这里采用联合字典学习方法进行字典对的学习。通过联合字典学习,只需一次稀疏分解,可以同时训练得到高、低分辨率两个字典,并且两个字典共享相同的稀疏编码系数。图约束联合字典对学习的目标函数如下

式中: 是联合高低分辨率块之后得到的列向量; 是高低分辨率联合字典; W1和W2分别是高低分辨率图像块列向量的维数; α = [α1,α2,…,αi]是稀疏编码系数矩阵,每一列是一个图像块的稀疏表示向量。

因此,高、低分辨率训练样本集可以表示为IpH=DpHS和IpL=DpLS。采用图约束位置块联合字典对训练,保证了SH和SL的一致性。

2. 2 重建系数的优化求解

在图像超分辨率场景下,通常假设高低分辨率图像块拥有同样的稀疏编码系数。因此,低分辨率图像块的稀疏编码系数可以被映射到高分辨率图像块字典上,从而生成高分辨率图像块,即

然而,实际超分辨率应用中,高、低分辨率图像块的稀疏编码系数并不完全一致,之间存在稀疏编码噪声σS= SH- SL。通过抑制稀疏编码噪声,可以提高重建过程中,高、低分辨率编码系数的一致性,即可提高超分辨率图像重建的精确度。因此,在求解输入低分辨率图像块的稀疏编码系数时,可以加入稀疏编码噪声约束以提高重建系数求解的精确度。优化重建系数求解的目标函数为

式中: γ 是正则化常量; lp范数用于表示SH和SL之间的距离。由于SH是未知的,因此,稀疏编码噪声无法直接计算。文献[9]提出了利用SH的稀疏编码均值E[SHN( k)]表示SH的思路。假设稀疏编码噪声近似于零均值随机变量,那么高、低分辨率图像块的稀疏编码均值非常接近,即E[SHN( p)]≈ E[SLN( p)]。式( 10) 可表示为

这里采用加权的K近邻块的稀疏编码均值来表示E[SLN( p)],距离越远的近邻块,权重越小,反之,权值越大。输入图像块的K近邻稀疏编码均值,采用以下公式计算获得

式中: Sp,k是第k个近邻块的稀疏编码系数; ωp,k是第k个近邻块的稀疏编码系数的权重; Np表示图像块p的K近邻块组成的集合,k ∈ Np。

2. 3 基于图约束字典和加权稀疏表示的人脸超分辨率算法

本算法主要包括字典学习和高分辨率图像重建两个阶段。前者是离线完成,得到一系列高低分辨率人脸图像块冗余字典对; 后者是在线完成,通过将重建高分辨率图像块进行交叠融合,得到最终的高分辨率人脸图像。整个算法过程表示如下。

1) 离线训练阶段

( 1) 输入训练集IH,IL,正则化参数 λ,β。

( 2) 构建相似矩阵W,矩阵中的元素定义如下

( 3) 计算矩阵

( 4) 联合训练过完备字典对

(5)重复(2)~(4),获取所有位置的字典对序列。

(6)输出高、低分辨率字典对序列

2) 在线重建阶段

( 1) 输入高、低分辨率字典对序列

测试图像IT,正则化参数λ,γ。

(2)For each LR patch IpTfrom IT

计算输入图像块ITp的K近邻稀疏编码均值; 计算输入图像块ITp的稀疏编码系数; 生成高分辨率图像块I*pT=DpHSp。

( 3) End for。

( 4) 交叠生成的高分辨率图像块,得到最终高分辨图像块IT*。

( 5) 输出高分辨率图像IT*。

3 实验结果及分析

采用MATLAB 7. 8 ( R2009a) 作为仿真实验平台在PC机上( 双核3. 20 GHz CPU,4 Gbyte RAM ,Windows7 操作系统) 实现了提出的算法。实验中,采用AR[15]人脸数据库作为训练和测试图像。AR数据库有4 000 多张正面人脸图像组成,分别来自126 个人,其中,男性70 人、女性56 人,每个人获取26 张不同表情、不同照度和不同遮挡程度的图像。这些图片分两个时间段分别获取。本实验中,选取AR数据库中同一时间获取的100 人( 其中,男女各50 人) 不同表情、不同照度情况下的正面图像用于实验,每人选取7 张图像。所有的人脸图像均被裁剪为120 × 160 像素大小,选取12 个特征点进行仿射变换对齐。部分对齐和归一化处理后的AR人脸图像如图1 所示。

将对齐后的人脸图像,进行降质处理,得到与高分辨率图像对应的低分辨率人脸图像。降质过程如下

式中: X,y分别表示高、低分辨率人脸图像; B表示模糊核为8 ×8 的平均模糊操作; U表示4 倍下采样操作; n表示均方差为10 的高斯加性噪声。在高、低分辨率人脸图像对中,随机选取90 个人的630 对人脸图像作为训练样本,剩下的10 个人的70 对人脸图像作为测试图像。图像分块的大小为8 ×8 像素,相邻块交叠32 像素。

为了验证算法的有效性,将本算法重建图像的主客观质量,与双三次插值算法、Jung提出的算法[9]进行比较。Jung方法中的正则化参数设置为0. 05。本文提出方法中 λ = 10,β = 0. 03,γ = 0. 05,K = 5 。

图2 为主观质量比较结果。从图中可以看出,双三次插值算法重建的结果非常平滑,但是不够清晰,也难以辨识。这主要是由于双三次插值算法在放大重建低分辨率图像过程中,并没有增加任何新的信息。Jung提出的方法都是基于块位置的人脸超分辨率方法,较之双三次插值方法,结果图像的清晰度有了明显的改进。但是,Jung提出的方法在重建图像的轮廓边缘都有明显的鬼影效应,重建人脸图像有或轻或重的眼镜假象。较之参考算法,本文提出的基于图约束稀疏编码的人脸超分辨率方法明显地改善了重建结果图像的清晰度,也减少了重建带来的人工效应和鬼影效果。本文方法取得了比参考算法更好的主观质量。

表1 为客观质量比较结果。较之参照算法,本文提出的方法重建的结果图像在客观质量方面具有较高的峰值信噪比( PSNR) 和结构相似性( SSIM) ( PSNR和SSIM值均是越大越好) ,这表明本文方法重建的结果图像更接近于原始的高分辨率图像。

总之,无论从人眼视觉效果,还是客观评价指标,均表明本文提出的方法可以更好地对人脸图像进行超分辨率处理,获得更好的图像重建质量。

4 小结

本文提出了一种基于图约束字典和加权稀疏编码的人脸超分辨率算法。该算法在字典训练阶段引入图约束正则项,利用加权稀疏表示提高重建系数精确度,从而改进了已有算法性能。主客观仿真实验结果证明了该方法的有效性。后续研究中,自适应正则化参数的选择将有待进一步完善。

基于异质局部特征的图像稀疏表示 篇2

关键词：异质局部特征,稀疏学习,视觉词典,基于内容的图像检索

图像特征被广泛应用于图像处理与分析中作为图像的有效表达方式, 其典型应用之一就是基于内容的图像检索 (Content-Based Image Retrieval, CBIR) 。给定查询图像, CBIR直接从图像库中查找与之视觉特征相似的图像, 这种“图找图”式的依据视觉特征相似度给出图像检索结果的方法克服了传统的基于关键字的图像检索技术“字找图”式的不足, 成为当前智能信息检索领域的研究热点[1]。

近年来, 图像局部特征因其良好的可重复性、可分辨性和鲁棒性得了广泛重视和飞速发展, 很多图像局部特征被相继提出 (如SIFT[2], HOG[3]等) , 并用于CBIR等视觉任务中, 弥补了图像全局统计特征 (如颜色、形状、纹理等) 的不足。然而, 图像上所提取的局部特征因其数量较多、维数较高也给大规模图像检索任务带来了新挑战。这种“高维局部特征集”表示图像的方式, 由于进行相似性度量的时间和空间复杂度高而难以适应大规模数据库环境下的图像检索任务, 而图像的全局表示形式 (Holistic Representation) 用一个向量表示一幅图像, 在这方面则有其天然优势, 因为任何两幅图像的相似性直接可用向量之间的距离函数来度量。

为了使图像的表示形式既能像局部特征一样有描述图像细节信息的能力, 又能像全局特征一样简洁明了, 本文提出利用“稀疏学习”的思想, 从训练图像的特征数据中建立超完备视觉词典, 采用局部稀疏编码 (Local Sparse Coding) 和最大值合并 (Max Pooling) 将图像“高维局部特征集”转化为更高维的稀疏特征向量[4], 然后通过直接计算向量相似性得到图像相似性, 并应用于CBIR系统中。另外, 本文不是使用单一局部特征, 而是选取了信息互补的不同局部特征构成异质局部特征, 从而能从多角度描述图像的内容, 在CBIR系统中能得到比单一局部特征更好地检索结果。

1 相关工作

如何从原始“像素级”表示的图像中提取有更强表示能力的图像特征一直是图像分析任务首先要解决的问题, 也一直是计算机视觉和模式识别领域的研究热点。

为了适应大规模图像数据库环境下的应用, 可借助学习机制将提取的图像底层局部特征的集合通过“多对一”映射 (编码) 成为一种全局表达方式, 即用一个向量来表示一幅图像, 以便在图像检索任务中使用通用的相似性度量方法来比较图像的相似性, 迅速返回检索结果。不同于通常所说的颜色、形状、纹理等全局统计特征, 这是一种构造在局部特征基础上的全局特征, 它仍能保留接近于图像底层的细节信息, 可看作是较高一层的特征表示, 这就向着“语义层”特征表示又前进了一步。目前这类特征中比较典型的例子是Bo VW (Bag of Visual Words) [5], S.Lazebnik等人则使用SPM (Spatial Pyramids Matching) [6]方法对其进行了改进, 在一定程度上加入了局部特征的空间分布信息。

另外, 压缩传感和稀疏表达理论近年来在信号处理、模式识别和计算机视觉领域中掀起新一轮热潮, 在人脸识别、场景分类等诸多应用中都取得了较好效果。其中一个核心的概念就是稀疏编码, 最早源自Barlow等人对生物视觉系统研究而提出的有效编码假设[7]。Olshausen和Field则进一步提出了著名的稀疏编码模型[8], 该模型通过基向量 (或基函数) 线性相加表示输入图像, 在最小均方差意义下使重构图像尽可能地与原图像相似, 同时要求表示系数尽量稀疏化。在此基础上, 很多研究者在稀疏编码模型的理论和应用方面做了大量的工作, 取得了丰硕成果, 也提出了许多改进的稀疏编码模型。

对于图像数据, 这些稀疏编码模型大多是从自然图像中“随机”选取若干图像块 (按像素灰度值排列成多维向量) 构成一个训练集合, 加以训练学习后得到基向量和对应图像的稀疏表示 (编码) 。随机采样的子图像块作为样本会给学习过程带来不稳定性 (比如引入背景或非目标区域噪声、对图像尺度、方向、视觉、亮度变化敏感等) , 从而学习到的基向量不一定具有代表性, 并可能存在大量噪声信息冗余。针对以上这些问题, 考虑到目前流行的图像局部特征 (如SIFT等) 本身就是对图像中感兴趣区域 (图像块) 的一种优于“像素级”的多维向量表示形式, 因此, 直接以图像局部特征作为训练样本, 并采用稀疏编码模型学习基向量和图像的稀疏表示, 是特征学习方法研究的新趋势。研究者们最近几年来在这方面做了一些尝试, 并在图像重构、图像分类等应用中取得了较好效果[9,10,11]。

但其训练数据 (局部特征) 大多是在图像上密集采样的结果, 一般都只使用单一的视觉特征。密集采样得到的特征数比基于兴趣点检测得到的特征数要多得多, 而且极易引入背景和非目标区域噪声, 另外单一视觉局部特征一般是精心设计的, 是对图像块某一属性的描述, 如果还能联合其他信息互补的局部特征, 如基于兴趣点检测的特征加上基于局部纹理或形状描述的特征, 则描述能力会更强。

本文正是以此为突破点, 运用稀疏学习的思想, 将SIFT (Scale Invariant Feature Transform) [2]、LBP (Local Binary Patterns) [12]和HOG (Histograms of Oriented Gradients) [3]等3种信息互补的图像局部特征视为异质局部特征进行融合, 最终以一个高维稀疏向量的全局表示形式描述图像多角度的视觉内容, 并将其应用于CBIR任务中。

2 图像异质局部特征的稀疏学习

2.1 局部特征的稀疏学习

图像的局部特征可以看作是对图像某一采样区域特性的向量描述。例如, SIFT特征是基于“团点”检测的, 对图像缩放、旋转、光照变化甚至遮挡和裁剪等均保持着较好的不变性;LBP特征反映了图像上像素点与其近邻像素点灰度值的大小关系, 描述了图像的局部纹理特性;HOG特征则描述了图像内容的局部形状或边缘特性。

图像局部特征稀疏学习的目的, 是利用学习机制将图像底层局部特征的集合通过“多对一”映射 (编码) 成为图像的全局稀疏表示形式, 以便在图像检索任务中使用通用的相似性度量方法来比较图像的相似性。图像局部特征的稀疏学习过程如图1所示。

一方面, 选取图像库中的部分图像作为训练图像, 提取底层局部特征, 通过聚类方法得到初始的超完备视觉词典, 然后利用初始化的视觉词典和训练图像的局部特征, 交替使用词典学习方法和稀疏分解算法, 通过不断的训练学习得到优化的超完备视觉词典和图像库中图像的稀疏特征 (即全局稀疏表示) ;另一方面, 应用系统的输入图像 (如CBIR系统的查询图像) 的局部特征被提取, 并利用训练好的词典对其进行稀疏学习, 得到输入图像的稀疏特征。随后, 这些稀疏特征可被用于各种具体计算机视觉应用中。

记X=[x1, x2, …, xn] (x1∈Ra×1) 为输入矩阵 (每列是一个输入向量) , 表示在d维空间中的一组包含n个局部特征向量的集合, B=[b1, b2, …, bk] (b1∈Ra×1) 为基矩阵 (每列是一个基向量) , 表示由K个基向量构成的视觉词典, S=[s1, s2, …, sn] (s1∈RK×1) 为系数矩阵 (每列是一个系数向量) , 表示利用视觉词典进行稀疏分解 (局部稀疏编码) 得到输入矩阵X的稀疏编码矩阵, 则以上稀疏学习的过程可以表示成下面的优化问题

式中:‖x1-Bsi‖2表示重构误差;si是稀疏性的惩罚函数;λ为规则化参数, 用于权衡重构误差和稀疏约束。该优化问题在S保持不变时是关于B的凸优化问题, 在B保持不变时是关于S的凸优化问题。一般通过交替固定B和S之一的同时优化另一个的方法来优化上述目标函数。

对于学习基矩阵B (即学习视觉词典) , 此时固定S, 该优化问题等价于平方约束最小二乘问题

对于学习系数矩阵S (即局部稀疏编码, 学习局部特征的稀疏编码矩阵) , 此时固定B, 该优化问题等价于L1规则化最小二乘问题

为了将图像用一个向量表示, 对学习到的局部特征的稀疏编码矩阵, 还要进行一个合并操作, 一般采用最大值合并 (Max Pooling) 方法[9]

式中:scj是sc (最终的高维稀疏向量) 的第j行元素;Sj是S的第j行第i列的矩阵元素;n是局部特征向量的数目。最大值合并相当于在对应基向量位置的最强响应, 许多图像分类任务已证实该方法行之有效[4], 故本文也采用最大值合并方法来合并各个稀疏编码, 从而得到整幅图像的稀疏表示。

2.2 异质局部特征的稀疏学习

不同的局部特征, 其设计思路不同, 对图像底层细节信息描述的角度也就不同。这种信息互补的特征组合可视为异质局部特征 (Heterogeneous Local Features) 。本文从众多的图像局部特征中, 选择了如前所述的SIFT (128维) 、LBP (采用P=8, R=1统一模式LBP, 58维) 和HOG (36维) 来构成异质局部特征加以研究。

为了融合图像的异质局部特征, 采用如图2所示的稀疏学习方法。

从图像数据库中选择一部分图像作为训练图像, 分别提取SIFT, LBP, HOG特征组成3个训练特征集, 分别得到3个超完备视觉词典B_sift, B_lbp, B_hog。对于训练好的每一个视觉词典, 利用其对图像的局部特征矩阵进行局部稀疏编码和最大值合并, 分别得到稀疏特征sc_sift, sc_lbp, sc_hog, 最后按照一定的权重进行首尾相连并进行归一化就能得到一个信息融合后的稀疏特征———单位向量sc_slh, 即为该图像的最终全局稀疏表示形式。

这样, 每幅图像仅用一个包含图像多角度局部信息的高维稀疏向量描述。图像相似性可直接用向量相似性来衡量。用这个稀疏特征向量来描述图像的特征, 相比单一特征对图像进行了更全面的描述, 又具备全局特征的形式, 因而这样的特征既具备了较强的图像局部信息描述能力, 又能够适应大规模数据库检索的需求。

3 应用实例:CBIR

将按以上稀疏学习方法得到的图像稀疏表示应用于基于内容的图像检索 (CBIR) 任务中。

3.1 图像库

采用标准图像库Zu Bud[13]。该库包含201栋建筑物的1 005幅图像, 每个建筑物各有5幅图像, 原始分辨率为640×480 (本文实验中将其缩小到320×240以减少数据量) , 是在不同季节和天气条件下从不同视角由两个不同相机拍摄的, 还特别拍摄一些被树木遮挡的图像。采用该库来做图像检索实验能够方便地评估图像特征的性能, 如尺度不变性、方向 (旋转) 不变性、视角不变性、光照不变性以及抗干扰能力等。

3.2 视觉词典学习

为减少计算量, 对每建筑物各取1幅图作为训练图像集, 对这201幅图像提取SIFT, LBP, HOG特征分别组成训练特征集, 通过K-Means聚类得初始化视觉词典, 并按式 (2) 进行词典学习, 分别得到3个具有K (K=1 000) 个基向量的超完备视觉词典。

3.3 图像的稀疏表示

对于学习好的每个视觉词典, 利用其对全部1 005幅图像的局部特征矩阵按式 (3) 和式 (4) 进行局部稀疏编码和最大值合并, 先将每幅图像用一个K维的稀疏向量进行表示;然后, 通过加权联接的方式融合3个稀疏向量, 并进行归一化, 从而形成图像的最终全局稀疏表示, 即3K维的稀疏单位向量。

3.4 图像检索过程和性能评价

为了便于统计结果和评价检索性能, 本文取每幅库图像作为查询图像, 这样图像检索过程简化为, 用两个稀疏单位向量的内积 (夹角余弦) 的来衡量两幅图像的相似度 (内积越大越相似) , 并按相似度从大到小返回指定数量的图像作为检索结果。根据Zu Bud库特点, 指定返回结果图像数T=5, 即等于实际相关图像数, 故本文实验中单次检索的查准率与查全率相同。这里采用平均查准率 (Average Precision, AP) 作为性能评价标准, 即

式中:ni是第幅查询图像检索出的相关图像数目, N=1 005。

3.5 实验结果及分析

表1为按6种不同的加权系数进行稀疏特征融合时图像检索实验的平均查准率。

图3和图4分别是实验中某幅查询图像利用单一SIFT稀疏特征及按0.5∶0.3∶0.2的权值进行异质特征融合后的5-近邻 (5-NN) 检索结果。

由此可见, 在本文设定的实验条件下, 相比单一局部特征, 综合利用异质局部特征进行图像检索, 能够得到更高的查准率, 异质局部特征对图像局部信息具有更全面的描述与区分能力。另外, 每幅图像均由一个高维的稀疏向量来表示, 因而只需要存储该向量中非零系数的值和索引, 且图像间的相似性直接用稀疏向量的距离函数来度量, 明显降低了直接用“局部特征集”表示图像时度量图像相似性的时空复杂度。

4 结论

本文提出了一种将图像的异质局部特征通过稀疏学习映射为图像全局稀疏表示形式的方法, 并将之应用于基于内容的图像检索任务中。文中选取了SIFT, LBP, HOG这3种典型的图像局部特征形成图像异质局部特征, 它们分别描述了图像的兴趣点特性、局部纹理特性和局部形状特性, 加权融合后对图像视觉内容形成了多角度、更全面的描述。

稀疏表示保持的鉴别特征选择算法 篇3

特征选择[1]用于从高维特征空间中选择特征子集,并保持特征子集的原始物理特性,根据使用类别标签与否,特征选择算法可分为非监督和监督两种,本文主要研究监督特征选择算法。经典的监督特征选择算法包括Relief F[2],Fisher Score[3]以及多簇 特征选择(Multi-Cluster Feature Selection,MCFS)[4]等,它们通过特征和类别标签之间的相关性来度量特征的重要性,但是大多数传统特征选择算法对每个特征的度量是独立进行的[3,5],并且将特征逐个添加至所选特征子空间,这种选择方式的局限性在于特征之间的相关性被忽略[4]。最近,l2,1范数正则化优化已经应用到特征选择算法,此类算法通过对特征选择矩阵进行l2,1范数最小化约束来选择特征[6,7]。

与此同时,稀疏表示作为一种基于部分数据的表示,已经吸引了越来越多的关注,并已广泛应用于模式识别和机器学习领域[8]。稀疏表示方法假设一个超完备字典中样本的稀疏线性组合可以重构一个给定的样本,例如Wright等提出的 基于稀疏 表示的分 类方法[9](Sparse Representation-based Classification,SRC),该方法的优化问题惩罚线性组合系数的l1范数,SRC尝试使用所有训练样本的稀疏线性组合来表示一个给定的测试样本,并且认为稀疏非零表示系数集中在测试样本的同类训练样本上。受到SRC的启发,很多基于稀疏表示的特征抽取算法出现,例如文献[10-11]提出的稀疏表示分类器引导的监督特征抽取算法,该算法旨在减少类内重构残差,并与此同时增加类间重构残差,但二者在目标函数的形式上有所不同,文献[10]采用比值方式文献[11]采用差值方式。与特征选择算法不同,特征抽取将原始特征进行转换从而实现数据降维,特征的原始物理特性发生变化。回顾经典的监督特征选择算法,却不存在与SRC直接关联的,本文提出了一种稀疏表示保持的鉴别特征选择(SRPFS)算法,旨在寻找一种线性映射使得在所选特征子空间中,样本的稀疏类内重构残差足够小并且稀疏类间重构残差足够大,并用于优化提出的l2,1范数最小化的目标函数。

1 基于稀疏表示的分类方法

假设样本集为X = [x1,x2,⋯,xn]∈ Rm × n,类别数为c ,类别标签向量z = [l1,l2,⋯,ln]T,其中ln表示X中第n个样本即xn的所属类,给定一个测试样本y ∈ Rm,然后使

用以训练样本为基础向量的超完备字典表示y ,如下:

假设式(1)的系统是欠定的( m< n ),通过求解如下的优化问题可以得到最稀疏解:

然而,公式(2)中的L0优化问题是NP难题而且非常耗时,最近的研究理论[12,13]表明式(2)也可以通过寻找以下优化问题的解决办法进行求解:

该优化问题可以在多项式时间内通过标准线性规

划算法来解决[14],或者使用一种更高效的算法[15],然后利用式(3)中求解的稀疏表示系数向量 δ 对y进行分类,令 φk:RM→ RM(k = 1,2,⋯,c) 表示一种能够从 δ 中选择出与第k类有关的稀疏表示系数的函数[9],即 φk(δ) ,然

后计算y及其第k类原型之间的残差:

如果rl′(y) = minkrk(y) ,SRC将y分到第l′ 类。

2 稀疏表示保持的鉴别特征选择

2.1问题描述

X中除第i个样本xi后,xi所属类为li剩余n - 1个样

本记为Xi= [x1,⋯,xi - 1,xi + 1,⋯,xn]∈ Rm × {n - 1},通过解决L1

优化问题 得到xi的关于Xi的稀疏表 示系数向 量

δiw定义为xi的稀疏类内表示系数向量,该向量的非零元素与li类训练样本相关;δib定义为xi的稀疏平均类间表示系数向量,该向量的非零元素与c类中除li类的剩余c - 1类训练样本相关,即c - 1类的稀疏平均类间表示系数向量;因此,xi的稀疏类内和类间重构分别表示为:

这里的目标是寻找一种特征选择矩阵U ∈ Rm × m进而选择出m′( m′ < m )个特征,U满足的条件为:元素只有‘0’或‘1’;每行或每列中‘1’的数目不超过1;只有m′ 行或列的元素为‘1’。

通过使用U可以使用特征选择后的类内以及类间训练样本对xi进行重构,即稀疏类内重构UTXiδiw以及稀疏类间重构UTXiδib,稀疏类内类间重构残差采用F-范数进行度量,表示如下:

基于SRC决策规则,希望在所选特征子空间中样本xi尽可能接近其稀疏类内重构并同时尽可能远离其稀疏类间重构,考虑所有样本,SRPFS的目标函数定义如下:

式中:β 是一个权衡参数;1m∈ Rm是一个元素为1的向量,然而式(10)是NP难题,因此将关于U的二元约束放松到l2,1范数最小化约束[6,7],此时目标函数可以重写为:

式中:α 是一个权衡参数;.2,1表示l2,1范数,l2,1正则化

项控制U的大小并同时保证U的行稀疏性(行元素接近于0),使SRPFS为数据表示选择出最具鉴别性的特征。

2.2优化

式(11)的向量形式表示如下:

式中 : Sw= [δ1w′,δ2w′,⋯,δnw′]∈ Rn × n, Sb= [δ1b′,δb2′,⋯,δbn′]∈Rn*n,δiw′以及 δib′分别定义为:

令:

对L(U) 关于U求导,可以得到下式:

t = t + 1 ;

∂L(U)

直到收敛准则满足;输出:U 。

= 2αPU - 2XSw(XT- STwXTU) + 2βXSb(XT- SbTXTU)

∂U

(16)

2.3 L(U) 的凹性研究

式中P是一个对角矩阵,第r (r = 1,2,⋯,m) 个轴元素为Prr= 1 (2U(r,:)2)。为了求解式(12)中的U ,对L(U)关于U求导,然而很难用理论证明L(U) 是凹函数,将∂L(U) ∂U置为0,得到关于U的更新规则:

∂L(U)

对式(16)中的

关于U求导,得到下式:

∂U

∂2L(U)

= 2αP + 2XSwSTwXT- 2βXSbSbTXT

(18)

∂U2

根据凹函数的性质,式(15)中的L(U) 是凹的,当且

U ← (αP + XSwSTwXT- βXSbSbTXT)-1(XSwXT- βXSbXT)

∂2L(U)

仅当式(18)中的

是正定的,令:

∂U2

(17)

G = 2XSwSTwXT- 2βXSbSbTXT

为了得到最优U ,重复上述过程直到收敛标准满足,即|L(Ut)- L(Ut + 1)|< ε ,算法1给出了优化的详细过程描述。

(19)

因此式(18)可以重写为:

∂2L(U)

= 2αP + G

(20)

∂U2

算法1:SRPFS算法

2αP是正定的因为它是一个轴元素为正数的对角矩阵,根据正定矩阵的定义,如果G是正定的很容易证明2αP + G是正定的,然而很难直接证明G的正定性,事实上通过在实验中对参数 β 进行控制来保证G的正定性,β 的取值在实验部分给出。在假设2αP + G是正定的前提下,通过下面的定理证明目标函数在算法1中的迭代过程中的收敛性:

输入:训练数据集X = [x1,x2,⋯,xn]∈ Rm × n,类别标

签向量z = [l1,l2,⋯,ln]T

,权衡参数 α ,β ;

将特征选 择矩阵U初始化为 单位矩阵 即

U0∈{0,1}m*m,迭代次数:t = 1 ;

|

n

通过式(5)计算稀疏表示系数向量 δi|

;

|

i = 1

通过式(13)和式(14)计算Sw以及Sb;重复:更新对角矩阵Pt,即:

定理1:式(12)中的目标函数值在算法1中的迭代过程中单调减小。

1

é

ùúúúú

ê2Ut - 1(1,:)2

证明:很容易证明式(12)就是解决以下的问题:

êêêêë

Pt=

⋯

+ αTr[UTPU]

mUin X - UTXSw2

- β X - UTXSb2

1

F

F

2Ut - 1(m,:)2ûú

(21)

通过式(17)更新Ut;

相应地,在第t次迭代中有:

2

2

- βX

- (Ut + 1)TXSb

+ αTr[UTPU] ⇒X

- (Ut + 1)TXSw

Ut + 1= mUin X - UTXSw2F

- β X - UTXSb2

+

F



F

F2

2

αTréë(Ut + 1)TPtUt + 1ùûX

- (Ut)TXSw

- βX- (Ut)TXSb

+ αTré(Ut)TPtUtù

F



F

ë

û

即:

utr+ 122

utr222utr2

2

2

2

2

X

- (Ut + 1)TXSw

- βX

- (Ut + 1)TXSb

- βX- (Ut)TXSb

X- (Ut)TXSw

+ α∑

+ α∑

⇒

F



F

F

æ

F



2utr2

r

r

utr+ 122

ö÷÷ø

2

2

+ α∑utr+ 12- αçç∑utr+ 12-∑r r r

X

- (Ut + 1)TXSw

- (Ut + 1)TXSb

- βX



F



F

2utr2

è

utr222utr2

æ

ö÷÷ø

2

2

- βX- (Ut)TXSb

+ α∑utr2- αçç∑utr2-∑r r r

X

- (Ut)TXSw

F

F



è

根据[7]中的引理,对于任意非零向量u ,ut∈ Rm,下面的不等式成立:

因此,有以下不等式成立:

即:

utr22

+ α∑utr+ 12X

+ αUt + 12,1X- (Ut)TXSw

它表明式(12)中的目标函数值在算法1中的迭代过程中单调减小。

3 实验

在本节中,通过实验验证算法SRPFS的性能,首先将SRPFS与经典的监督特征选择算法进行比较,然后分析SRPFS的收敛性。

3.1实验设置

4个公共数据集:Wine[16],Breast Cancer Wisconsin(Diagnostic)[16],Connectionist Bench (Sonar,Mines vs.Rocks)[16]以及COIL20[17],Wine,Breast Cancer Wisconsin和Connectionist Bench来自标准UCI库;来自哥伦比亚图像数据库的COIL20包含20个对象,数据集的描述在表1中给出。

将SRPFS与All Features,Fisher Score,MCFS,以及Relief F进行比较,实验中为保证式(20)中G的正定性,β 在4个数据集上分别设置为10-3,10-5,10-3,10-2,使用快速 迭代收缩 阈值算法(Fast Iterative Shrinkage and Thresholding Algorithm ,FISTA)[16]求解式(5),FISTA中的规范 化参数设 置为1,α 的调整范 围为{1,10-1,10-2},对于MCFS以及Relief F邻居样本数设置为5,由于Connectionist Bench和COIL20的特征数大于50,相应的所 选特征数 分别设为 {1,2,⋯,30} ,{1,2,⋯,512} ,即最大值取数据集维度的50%。

3.2分类识别率比较

对于每个数据集,随机选择每类样本集的5种方法在4个数据集上的平均最高识别率(±std)的比较,如表2所示。选择的样本中80%做训练集,剩余样本做测试集,为了证明不同算法的可靠性,将训练集以及测试集的选择过 程重复10次 ,All Features,Fisher Score,MCFS,Relief F以及SRPFS在4个数据集上的平均最高识别率及标准差在表2中给出,可以看出所有的特征选择算法优于All Features,因此,特征选择算法有助于提高识别率,由于SRPFS中保持了样本之间的稀疏相关性,SRPFS从识别率和稳定性两方面的性能明显优于其他方法。

3.3收敛性

在本节中,通过实验证明所提出的迭代算法单调减小目标函数值,直到收敛,图1展示了式(12)中的目标函数值在4个数据集上的收敛曲线图,可以看出目标函数在数次迭代后收敛。

4 结语

在本文中,提出了一种新的监督特征选择算法,称为稀疏表示保持的鉴别特征选择(SRPFS),其目的是选择鉴别性特征子集,使得在所选特征空间中,样本的稀疏类内重构残差和稀疏类间重构残差的差值最小化。通过实验验证SRPFS的性能并与其他4种方法即All Features,Fisher Score,MCFS,以及Relief F在4个公共数据集上进行比较,实验表明SRPFS在识别率以及稳定性方面明显优于其他方法。在未来,考虑将SRPFS的思想应用到非监督特征选择算法研究中,由于不使用样本的类别标签这将是一个更大的挑战。

摘要：稀疏表示作为一种基于部分数据的表示,已经吸引了越来越多的关注,并广泛应用于模式识别和机器学习领域。提出一种新的算法,称为稀疏表示保持的鉴别特征选择(SRPFS),其目的是选择鉴别性特征子集,使得在所选特征子空间中,样本的稀疏类内重构残差和稀疏类间重构残差的差值最小化。与传统算法选择特征的独立性方式不同,该算法以批处理方式选择最具鉴别性的特征,并用于优化提出的l2,1范数最小化的目标函数。在标准UCI数据集和哥伦比亚图像数据库的实验结果表明,该算法在识别性能和稳定性方面优于其他经典特征选择算法。

加权稀疏表示 篇4

为了克服上述缺点,一般采用阵元非均匀排布的稀疏阵列。将天线阵列以较少的阵元数进行稀疏布置,能够有效减弱阵元间的互耦效应,增大阵列的孔径从而提高角度分辨率以及减少系统成本和降低软硬件复杂度。阵列的稀疏布置往往使得波束方向图的旁瓣电平抬高,因此,稀疏阵列综合的目的是在满足期望的波束辐射特性要求下,求解使得阵列最大稀疏化的阵元位置和阵元幅度的最佳分布。然而,非均匀间隔的稀疏阵列综合属于多变量的非线性优化问题[3],处理起来非常困难。传统的智能优化算法如遗传算法[4]、模拟退火法[5]、粒子群算 法[6]等常应用于求解稀疏阵列的综合问题。它们都是先设计使得阵列合成波束峰值旁瓣电平最小化的阵列综合目标函数,然后对该目标函数进行搜索处理,寻找满足条件的最佳的阵元分布位置及其激励幅度。由于这些传统智能优化算法都是基于随机性的自然算法,在优化过程中需要较多的搜索处理运算,导致此类算法的收敛速度慢,并容易陷入局部最优解等问题。

最近,稀疏信号重构理论开始在阵列信号处理中得到了应用,为阵列综合提供了一种新的思路。 由于稀疏阵列的阵元在空间上是离散和稀疏分布的,因此稀疏阵列综合在本质上等效于稀疏信号重构问题。文献[7—9]利用基于迭代加权范数最小化的稀疏信号重构方法来求解阵列的综合问题, 在每次迭代中利用凸优化软件计算出用于下次迭代的阵列加权向量,当满足迭代终止条件时,由阵列加权向量的非零值确定阵列的阵元位置及其激励,此类方法利用稀疏阵列天线的稀疏物理特性,能以较少的迭代次数获得稀疏程度更高的天线阵列。针对凸优化软件的使用会导致阵列综合的耗时增加和通用性变差的问题,文献[10]提出一种新的基于迭代加权范数的阵列综合方法,该方法利用复数求导结合启发式近似方法获得每次迭代过程中阵列加权向量的闭式解,无需使用软件工具从而实现稀疏阵列的快速综合。与范数方法相比,范数方法能以更少的测量数据来准确重构稀疏信号[11]。因此,文献[12]在MIMO ( multiple-inputmultiple-output) 雷达中利用迭代加权范数最小化方法来准确估计空间稀疏分布目标的角度,距离和多普勒信息,从而提高了MIMO雷达的成像质量。本文将迭代加权范数[12]应用于求解阵列综合问题中,建立了迭代加权范数最小化的阵列综合模型,并利用拉格朗日乘数法[12]求解出用于下次迭代的阵列加权向量闭式解,当达到迭代终止条件时,由当前的阵列加权向量的非零值来获得阵元位置和阵元幅度的最佳分布。仿真结果表明,相比文献[10]提出的基于迭代加权范数的阵列综合方法,本文方法能以更少的迭代次数和阵元数来满足阵列波束辐射特性,从而有效提高了阵列综合的性能。同时,讨论了在满足相同辐射特性的情况下利用本文方法进行稀疏阵列综合时所需的迭代次数和阵列稀疏程度与q取值的变化关系。

1阵列模型

在阵列综合前,假设所有N个各向同性初始阵元都均匀分布一条直线上,各阵元位置分别为x1, x2,…,xN,则其远场方向图可以表示为

式( 1) 中,xn= ( n - 1) d表示第n个阵元相距第一个阵元的距离,其中d为阵元间距; k = 2π/λ( λ 为波长) 表示空间波数; θ 表示阵列的扫描角( 偏离阵列法线方向的角度) ,- π/2 < θ ≤π/2; wn表示第n个阵元的激励,如果wn= 0 ,则等效于该位置上无天线阵元存在。阵列综合的目的是需要寻求一组非零元素 最少的阵 列加权向 量w= [w1,w2,…,wN]T,使得综合出的阵列方向图尽可能接近期望的阵列波束辐射特性,其中[·]Τ表示矢量转置。

在观测角度区间 - π/2 < θ ≤ π/2内等角度间隔选取P个观测点,即为 θ1,θ2,…,θP,则由这P个角度采样点所对应的阵列流形矩阵( 导向矢量阵) 可表示为

式( 2) 中,对应角度采样点 θl( l = 1,2,…,L) 的阵列的方向向量表示为

在这P个角度观测点上的波束方向图采样值用F = [F( θ1) F( θ2) …F( θP) ]T表示,则F可表示为

2稀疏阵列综合方法

为了使得综合后的阵列在P个角度观测点上的方向图F = [F( θ1) F( θ2) … F( θP) ]T尽可能接近期望 的阵列波 束方向图Fd= [Fd( θ1) Fd( θ2) … Fd( θP) ]T,同时希望阵列加权向量w尽可能稀疏,即w非零元素尽可能少。则采用迭代加权 q( 0 < q < 1) 范数最小化[12]来表示该阵列综合模型,即为

式中,{ zi> 0}Ni= 1表示一组加权系数; {·}( l-1)表示第 ( l - 1) 次迭代; ‖·‖2‖2范数; ε 表示最大能够接受的误差。在式( 5) 中,如果q等于1,则式( 5) 和式( 6) 所描述的阵列综合模型等效于文献[10]中基于迭代加权范数最小化的阵列综合模型。

由于式( 6) 属于有约束的阵列综合模型,不便求导计算。因此,通过构造拉格朗日乘数法[11]将式 ( 6) 约束综合模型转化为无约束综合问题,

式中,μ ∈ R+表示正则化参数。将式( 7) 表示的无约束综合目标函数f( w,z) 进行展开,

式( 8) 中,(·)H表示共轭转置。

为了得到阵列加权向量w的表达式,在式( 8) 中对wH进行求偏导数,可得

式( 9) 中,Q的值可由式( 10) 来表示

式( 10) 中,diag( ·) 表示对角化操作。令 ∂f( w, z) / ∂wH= 0 ,得到如下非线性方程

假设已知在第l次迭代中的加权系数 { z(il)}Ni= 1和已经估计获得阵列加权向量 { w(il)}Ni= 1,则通过式 ( 11) 可以计算出下次迭代中的阵列加权向量[12],

式( 12) 中,I为P阶单位矩阵,加权系数zi和正则化参数 μ 分别由式( 13) 和式( 14) 所决定

式中,参数 δ > 0是一个非常小的常量,避免当加权向量w( l)中出现零元素值时导致算法运行终止。当满足迭代终止条件时,由当前迭代中计算得到的阵列加权向量非零值确定阵列的阵元位置及其激励。

本文提出的基于迭代加权范数最小化的稀疏阵列综合方法具体实施流程如下:

( 1) 初始化: 根据初始化阵元数N和阵元间距d,确定阵列流形矩阵A; 在第一次迭代中即 ( l = 1 ) 时,加权系数z( l)和w( l)均为全1向量。

( 2) 利用式( 13) 和式( 14) 计算更新加权系数z(il)和正则化参数 μ 。

( 3) 更新阵列加权向量w( l+1):

( 4) 重复步骤( 2) 和( 3) ,直到迭代次数l达到lmax或者满足下面的终止条件:

式( 16) 中,ξ 为设定的最小误差值。

( 5) 确定稀疏阵列的阵元位置和激励: 由达到终止条件时加权向量w估计值中非零元素所在的位置确定为综合后稀疏阵列的阵元位置,而该非零元素值的大小即为该阵元的激励幅度值。

3仿真实验

为了验证本文提出的基于迭代加权范数最小化阵列综合方法的有效性,本文设计了一些仿真对比实验。仿真参数设置如下: 初始阵元间距d = 0. 25 λ ,各阵元在0到175λ 范围内以d为间距均匀地分布,则共有701个初始天线单元,角度采样点数P = 49 ,常量 δ = 1. 0 × 10-6,最小误差值 ξ = 1. 0 × 10-4。满足迭代终止条件时,将阵列加权向量w中小于0. 000 1倍max( w) 的元素全部设置为零,表示对应的位置不放置天线阵元。设置期望对称波束方向图的主瓣宽度|sinθ| ≤ 0. 05 ,即主瓣宽度限制在[- 2. 86°,2. 86°]内,旁瓣电平限制在 - 25 d B。由文献[10]可知,基于迭代加权 1范数的稀疏阵列综合方法其实是本文方法在q = 1时的特例。因此,为了方便,在下面仿真分析中q = 1对应文献[10]提出的基于迭代加权 1范数的稀疏阵列综合方法。

图1分别给出了由本文方法 ( q = 0. 5) 和基于迭代加权范数的稀疏阵列综合方法[10]重构得到的阵列阵元位置分布及其对应激励值。图2为两种方法综合后的波束方向图与期望波束方向图。由图1和2可知,在达到期望波束辐射特性的情况下,本文方法综合后得到的稀疏阵列的阵元数为43,而基于迭代加权范数的稀疏阵列综合方法需要50个阵元数,因此本文方法节省了7个阵元,而且本文方法仅需17次迭代运算就可完成阵列综合,比基于迭代加权范数的稀疏阵列综合方法少花了8次迭代运算。从综合后阵列的稀疏程度以及所需迭代运算次数上来看,本文方法优于基于迭代加权范数的稀疏阵列综合方法。

图3和图4分别在满足相同辐射特性的情况下利用本文方法进行稀疏阵列综合时所需的迭代次数以及阵列稀疏程度与q取值的变化关系,其中q = 1对应基于迭代加权范数的稀疏阵列综合方法。由图3和图4可知,与基于迭代加权范数的阵列综合方法相比,本文方法能以更少的迭代次数来综合出稀疏程度更高的稀疏阵列来达到期望的波束辐射特性。兼顾阵列综合的稀疏程度和所需迭代运算次数,本文方法在q = 0. 5时阵列综合的性能最为理想。

4结束语

本文将迭代加权范数最小化方法应用于稀疏阵列综合中,并利用拉格朗日乘数法求解该阵列综合数学模型,以获得用于下次迭代的阵列加权向量闭式解,反复迭代直至满足给定的终止条件时,由迭代结束时的阵列加权向量的非零值来获得阵元位置和阵元幅度的最佳分布。与基于迭代加权范数的稀疏阵列综合方法相比,本文方法能以更少的迭代次数重构出稀疏程度更高的阵列来满足期望的波束辐射特性,而且无需使用优化工具来求解阵列综合问题,提高了阵列综合的实时性和通用性。

摘要：针对非均匀稀疏阵列综合问题,提出一种利用迭代加权lq(0<q<1)范数最小化的阵列综合方法。该方法利用稀疏阵列天线的稀疏物理特性,将稀疏阵列综合问题转化为一系列迭代加权lq(0<q<1)范数最小化的稀疏重构问题,并在每次迭代中求解出用于下次迭代的阵列加权向量闭式解,由满足迭代终止条件时的阵列加权向量的非零值来确定阵列的阵元位置及其激励幅度。仿真结果表明,与基于迭代加权l1范数的阵列综合方法相比,该方法在满足辐射特性前提下能以更少的迭代次数来综合出稀疏程度更高的稀疏阵列。

加权稀疏表示 篇5

差分图像在许多领域得到了广泛的应用,比如:视频压缩,生物医学诊断,天文学,遥感等[1,2,3,4]。本文提出一种新的采用压缩测量值得到差分图像的新方法。差分图像是由目标场景在相邻时间点的图像相减得到的,从而能够得到目标场景随时间的变换。与传统方法需要对目标场景进行成像不同,该方法通过采样矩阵得到压缩测量值,然后利用图像时空互相关性或者差分图像的稀疏性,采用线性或者是非线性的优化方法得到差分图像的估计值。压缩测量值的获取过程是将高维的目标场景数据投影到一个低维空间的过程,实现了数据的实时压缩,从而降低了设备的成本。采用该方法的另外一个优点就是能够降低系统的噪声敏感性,在信噪比较低的时候能够得到较好的差分图像估计结果。

该方法主要包含两方面的内容:首先设计最优的采样矩阵以得到压缩测量数据[5,6,7];采用闭环的迭代方法得到差分图像的估计值,包括基于l2和l1的方法。基于l2范数的方法能够由压缩测量值直接估计出差分图像而不需要首先重建目标场景。这种方法主要利用了目标场景子在连续时间点空间和时间上的相关性。采用基于l1范数的方法主要是能够利用差分图像的稀疏性,从而可以将问题转化为带有不同约束条件的线性逆问题进行求解。

1 信号模型

差分图像就是目标场景在连续时间点图像相减所构成的图像,如图1所示。图1(a)和(b)分别为包含运动目标的场景在连续时间点的图像,(c)即为差分图像。

假设xk和xk+1分别为目标场景在时间点tk和tk+1时的图像,则可以定义第k幅差分图像为:

广义的差分图像定义为目标场景在时间点tk和tk+L所成图像的差别:

与传统方法相比,压缩测量能够降低采样数据的维数,但是压缩测量方法是一种离散数据处理方法:将成像场景离散化为分辨率为δr的单元,然后将离散数据构成一个N×1的向量。压缩测量值的获得需要采用一个M×N的随机投影矩阵将高位数据投影到低维数据空间上,使得压缩测量值包含了恢复原图像所需要的全部信息。

假设x1和x2为两个连续时间点t1和t2的目标场景。在两个时间点目标场景均被离散化为N×1的向量。假设Φ1和Φ2为相应的M×N(M<N)的随机采样矩阵。则压缩测量数据可以表示为:

其中,n1和n2表示方差为δ2的0均值高斯噪声。

2 l2范数优化方法

采用基于l2范数优化方法就是由得到的压缩测量值y1和y2采用优化的方法得到差分图像Δx1,使得Δx1与真实差分图像差值的l2范数最小。

差分图像的求解可以直接由压缩测量值y1和y2得到。采用这种方法不仅能够利用像素间的空间相关性,而且能够利用不同时间点的时间相关性,因而能够得到更加准确的差分图像。假设所估计的差分图像为:

其中,W1和W2为联合最优化估计函数。我们称这种方法为直接差分图像估计方法(DDIE)。

在DDIE方法中我们假设场景在最初采样的时刻是确知的,因此在时间点t1我们可以假设采样矩阵Φ1为单位矩阵,而从t2及以后的时刻采用压缩的方法得到测量值。则式(3)和(4)可以表示为:

采用贝叶斯的估计方法,则目标函数为:

求解式(8)相对于W1和W2的微分,并令其等于0,则可以得到:

矩阵Rβ反映了由于噪声的存在所造成的相关信息的丢失。如果噪声为0,则Rβ为0。然而实际应用中不可避免的存在噪声,从而会造成相关信息的丢失,而Rβ正好反映了相关信息丢失的程度。在存在噪声的情况下,寻找最优的采样矩阵是不切实际的,因此通常采用一些次优的采样矩阵,如由主分量构成的采样矩阵(PC)、由差分图像的主分量构成的采样矩阵(DPC)以及加权的主分量构成的采样矩阵(WPC,WDPC)等[8,9,10]。

3 l1范数优化方法

采用基于l2范数的估计方法的优点在于其能够提供一种闭环的线性估计方法使得目标场景的MSE最小化,然而基于l2范数的估计方法并没有利用差分图像的稀疏性特点。除此之外在l2范数的估计方法中采用的时空相关矩阵虽然能够得到较好的估计结果,但是假设图像静止,这是不切实际的。为此可以采用l1范数的估计方法。采用基于l1范数的估计方法能够充分利用差分图像的稀疏性特点。除此之外,还能有效地采用基于l1范数的全变差方法有效地估计出差分图像的边缘特征。

基于11范数的差分图像估计方法可以看做式(11)所示的线性逆问题:

其中,D为线性变换。线性逆问题的目标就是由含噪测量值y估计出s。式(11)是典型的基于l1范数的信号重建问题,然而式(3)和式(4)所示的差分图像重建问题与式(11)不同。为了采用基于l1范数的方法,可以将式(3)和式(4)改写为:

其中:

假设x在稀疏矩阵Ψ上可以稀疏表示:

则式(12)可以表示为:

式(15)可以通过式(16)求解:

其中,R(s)为正则化分量。

为了充分利用差分图像的稀疏性,可以定义式(16)所示的正则化分量:

如果采用基于全变差的方法求解式(16),则正则化分量可以定义为:

其中,Riso和Rniso分别为各向同性分量和异性分量。和为水平方向和垂直方向的一阶差分算子。

差分图像的估计同样可以采用基于过完备字典的稀疏表示方法。假设稀疏表示字典为Ψ,且Ψ为对称双正交小波变换矩阵,正则化分量为。该方法可以通过基追踪方法(BPDN)求解[11]。

本文给出部分代码如下:

4 仿真分析

本节通过仿真来验证本文所提出的基于l2范数和l1的差分图像估计方法。将郊区的一幅视频图像作为仿真器的输入图像,如图1所示。然后采用PanasonicPV-GS500摄像机对目标场景成像。我们采用传统方法所成图像作为原始图像,然后通过仿真的方式得到压缩采样的光学成像系统的原因是在增加采用不同采样矩阵的灵活性。这种灵活性有利于评估本文所采用方法的性能。同时,我们需要考虑到计算的难易程度。为了降低计算量,不是直接对源图像进行处理,而是将原图像分解为:小块图像进行处理。仿真中所采用的压采样矩阵为PCs,DPCs,WPCs,WDPCs,高斯随机采样矩阵(GPR),以及由计算得到的最优采样矩阵。

图2所示为采用l2范数的差分图像估计方法。原图像大小为N=480×720,分块图像大小为8×8,压缩采样的数据维数为M=5,信噪比为10dB。比较了采样矩阵分别为WDPC和最优采样矩阵时的差分图像估计结果。从图2中可以看出当采样矩阵为最优采样矩阵时估计得到的差分图像与真实的差分图像更接近。采用WPDC作为采样矩阵依然能够得到很好的估计结果。同时需要主要到,WPDC矩阵构造简单,计算效率高;而最优采样矩阵计算非常困难。图中采用RMSE作为衡量图像估计效果的标准。

图3比较了采用l2范数方法时,压缩测量值对各种不同的采样矩阵的估计方法的结果。从图中可以看到,压缩测量数据维数的增加能够显著提高估计的结果。然而从图4中可以看出,并不是采样数据维数的增加就会使得估计结果无限提高。从图4中可以看出,采样数据维数的增加起初会使得估计结果提高;然而当采样数据维数增加到一定程度时反而会使得估计结果变差。这种现象是强调光子数量限制的结果。当M较小时,采样数据维数的增加能够提供更多的信息。然而优于光子的总数量是一定的,随着采样数据维数的增加,每一采样样本所包含的额外信息量减少,从而使得当采样数据维数增加时,估计结果反而降低。当采样矩阵为PC,DPC,WPC时都会出现这样的情况,然而当采样矩阵为WPDC和最优采样矩阵时,估计结果并不会变差。这是因为对于给定的SNR,WPDC和最优采样矩阵是最优的。只有当测量样本所提供的信息能够提高估计效果时才会被考虑,否则就会被忽略掉。对比图4(a)和(b)可以看出,当M相同的时候,分块越大,所得结果的质量越高。当分块大小为16×16时,M=1所能达到的性能和分块大小为8×8时M=5所能达到的性能接近。

从图3中能够看出不同的采样矩阵的有效性。从图中可以看出,采用根据噪声的统计特性对采样矩阵进行加权确实能够有效的提高估计的结果。这种提高在图3(b)中表现得更加明显。采用最优采样矩阵能够进一步提高估计的效果。

采用l2范数方法的优点在于其能够提供一种线性估计方法,使得能够快算简单的估计出差分图像。其缺点在于没有充分利用差分图像的稀疏特性。接下来验证一下采用l1范数方法的差分图像估计结果。所采用的采样矩阵与l2范数方法一致。

图5所示为采用不同的基于l1范数方法的差分图像估计结果。所采用的分块图像大小为8×8。从图中可以看出,三种方法都能够有效地估计出差分图像,其中TV方法的估计效果最优。采用各向同性和各向异性的正则化参数所得到的估计结果相似。

图6所示为采用TV方法时,采样数据维数不同对重建结果的比较。分块大小为8×8,采样数据维数分别为M=1,5,10,20。当M=1,5时,所有的采样矩阵的得到的估计结果相似。事实上,当SNR较低时,所有的采样矩阵都能得到比传统成像方法更好的成像结果。当M=10,20时,不同的采样矩阵会得到不同的估计结果。采样矩阵为PC和DPC时,估计结果类似。这是因为基于l1范数的方法不是直接估计差分图像,而是联合估计出两个时间点的图像。同样,采样加权的方法能够提高估计的效果,但是采样矩阵为WPC和WPDC时。

图7所示为采样矩阵为WPDC时,基于l1范数方法的估计结果比较。从图中可以看出,基于TV的方法能够得到最好的估计结果。但是,所得到的估计结果的提高有限。通过使差分图像的梯度最小化,基于TV的方法能够捕捉到差分图像间的强度变化。然而,差分图像本身是稀疏的,因此利用稀疏性限制同样能够得到较好的估计结果。和其它两种方法相比,BPDN方法的估计结果最差。这主要是因为BPDN方法没有直接利用差分图像的稀疏特性,而是计算两个不同时间点时场景的联合稀疏表示。

5 结语

文中提出了一种基于l1范数和l2范数的利用压缩测量样本的差分图像估计方法,并且对两种方法进行了质化和量化的分析。两种方法都有各自的优点。基于l2范数的方法能够提供一种闭环的线性估计方法,使得计算简单、方便。基于l1范数的方法充分利用了差分图像的稀疏性特点。在l1范数方法中,分析了三种不同方法的重建结果,其中基于TV方法的估计结果最优。对于不同的采样矩阵,估计结果会不同。其中最优采样矩阵的估计结果最优,但是计算量大,不易于构造。采用根据噪声的统计特性进行采样矩阵的加权能够有效地提高估计的结果。

摘要：差分图像能够显示目标场景随时间的变化。采用一种基于图像特征的方法,由目标场景的一组压缩测量数据得到差分图像。该方法主要包含两方面的内容:首先设计最优的采样矩阵以得到压缩测量数据;采用闭环的迭代方法得到差分图像的估计值,包括基于l2和l1的方法。采用基于l2范数的方法能够由压缩测量值直接估计出差分图像而不需要首先重建目标场景。这种方法主要利用了目标场景在连续时间点空间和时间上的相关性。基于l1范数的方法主要采用一种改进的全变差方法和基追踪降噪方法。仿真表明了该方法的有效性。

关键词：压缩感知,差分图像,全变差,基追踪

参考文献

[1]Xue D W,Lu Z M.Difference image watermarking based reversible image authentication with tampering localization capability[J].Int.J. Comput.Sci.Eng.Syst.,2008,2(3):213-218.

[2]Lee S K,Suh Y H,Ho Y S.Reversible image authentication based on watermarking[C]// Proc.IEEE Int.Conf.Multimedia Expo.,2006, 1321 -1324.

[3]Hahn D,Daum V,Hornegger J,Bautz W,et al,Difference imaging of inter- and intra-ictal spect images for the localization of seizure onset in epilepsy[C]//Proc.IEEE Nucl.Sci.Symp.Conf.Rec.,2007, 4331 -4335.

[4]Neifeld M A,Shankar P.Feature-specific imaging[J].Appl.Opt., 2003,42(17):3379 -3389.

[5]Candes E J,Wakin M B.An introduction to compressive sampling[J]. IEEE Transactions on Signal Processing,2008,25(2):21 -30.

[6]Babacan S D,Morula R.Bayesian compressive sensing using Laplace priors [J].IEEE Transactions on Image Processing,2010,19(1);53 -63.

[7]Mishali M,Eldar Y C.Xampling:Compressed Sensing of Analog Signals [M]."Compressed Sensing:Theory and Applications",Cambridge University Press,2011:78.

[8]Carvajalino J M D.Sapiro G.Learning to sense sparse signals:Simultaneous sensing matrix and sparsifying dictionary optimization[J].IEEE Trans.Image Process.,2009,18(7):1395 - 1408.

[9]Elad M.Yavneh I.A weighted average of sparse representations is better than the sparsest one alone[J].IEEE Trans.Inf.Theory,2009,18 (7):1395-1408.

[10]Hyder M M.Mahata K.An improved smoothed approximation algorithm for sparse representation[J].IEEE Trans.Signal Process.,2010,58 (4):2194-2205.

基于鲁棒稀疏表示的人脸识别算法 篇6

关键词：人脸识别,LBP算子,稀疏表示,RSRC

0 引言

LBP(局部二值模式)算子是一种有效的图像纹理描述算子,后来Ahonen[1]等将其成功地应用到人脸识别领域,并得到了显著的成果。该方法首先是将人脸分成若干个区域,然后对每个区域提取LBP算子,最后把若干个区域的LBP算子连接到一起,形成整个图像的LBP算子,这也就是我们所说的LBP特征。

SR(稀疏表示)最近是计算机视觉领域比较热门的研究话题,Ma Y等[2]将其成功用到了人脸识别领域。简单理解为,测试人脸图像可以用它所在类的训练样本图像线性表示,用l0范数或者l1范数最小化问题来表示编码系数矢量的稀疏性,传统的稀疏表示模型用公式表示:

其中,ε>0是一个常数,y是测试人脸图像,D是训练样本构成的字典矩阵,β是编码系数矢量。该模型存在的问题是对编码系数矢量β施加l1范数约束是否能够有效地表示人脸图像的稀疏性;对于该问题研究者已经做了很多工作。例如:Liu等人[3]对编码系数β施加一个非负约束;Gao等人[4]在稀疏编码里加入了系数的拉普拉斯余项;在现实的人脸识别问题中由于噪音、表情的干扰,l1范数项或l2范数项有可能不能鲁棒地和有效得衡量重构人脸图像的准确性。

本文提出了一种基于LBP特征的鲁棒SR人脸识别算法,首先提取训练样本和测试样本的LBP特征,在原有稀疏表示分类器(SRC)里添加一个权值函数W来去除人脸图像中的噪声或异常值,最后用鲁棒稀疏表示分类器(RSRC)来对测试样本进行分类。

1 LBP算子

LBP算子最先是由Ojala等首先提出的,是一种用来描述图像纹理的方法。原始的LBP算子是基于以下思想:相邻的图像像素之间具有相似性,因此我们可以首先选定一个3×3的邻域,用该邻域的周围的像素值和邻域的中心像素值相比较,如果该点像素值大于邻域中心像素值,那么我们记为1;反之,如果该点像素值小于邻域中心像素值那我们记为0。依次计算完所有点的邻域周围像素值,然后我们从左上角按顺时针排序得到一个8位的二进制数,将此二进制数转化为十进制数,即为该邻域的LBP算子值。

2 RSRC分类器

定义训练样本矩阵D=[r1,r2,…,rn],其中,行矢量r1是D的第i行元素,定义编码误差e=y-Dβ=[e1,e2,…,en],则e中的元素为ei=yi-riβ。假设e1,e2,…,en是独立同分布的,并且概率密度为fθ(ei),其中θ表示分布参数集,不考虑编码系数矢量β,编码误差e的最大似然估计表示为 ,两边取对数, ,其中ρθ(ei)=-lnfθ (ei)。编码系数矢量的最大似然估计,即鲁棒稀疏表示(RSRC),可以用公式表示为: 。RSRC模型实际上是一个稀疏约束最大似然估计问题,换句话说它是一个更一般的稀疏编码模型,传统的稀疏编码模型式(1)当编码误差服从高斯或者拉格朗日分布时的一种特殊情况,为了解决式(2)我们可以对施加稀疏约束,一个关键的问题在于如何决定概率分布函数f θ(ei),简单地假设f θ(ei)服从高斯或者朗格朗日分布过于简单并且也不够有效,这里我们把式(4)的最小化问题转化为一个迭代赋权值的稀疏编码(IRSRC)问题: ,其中,W是诱导权值矩阵,

Wi,j是一个对角矩阵它的元素由下面的表达式给出: 考虑到权值矩阵W是一个对角矩阵,它的元素是给测试图像像素赋权值的,我们把异常值像素(如:遮挡、腐蚀部分的像素)赋予较低的权重,我们采用SVM[5]的权值函数 其中μ和σ是正的标量,参数μ控制从1到0减少的速率, μ控制分界点的位置。综合式(7)、式(6)和ρθ(0)=0 ,得到 传统的SRC模型由于迭代次数很难估计,因此导致识别率较低。我们通过下式很好地解决了这个问题: ,其中,ω 是一个很小的正标量,当等式成立时迭代结束。最后我们根据文献[2]相同的方法来分类测试图像所属类别。

计算稀疏编码系数β的算法如下:

输入:用单位l2范数标准化测试图像y,训练样本矩阵D(D的每一列具有单位l2范数),

初始化为yini,yini=md是所有训练样本均值图像。

输出:β

3 实验验证及结果分析

AT&T人脸库:由剑桥大学AT&T实验室创建,包含40人共400张面部图像,每张图片尺寸是92×112,部分志愿者的图像包括了姿态、表情和面部饰物的变化。该人脸库在人脸识别研究中经常被人们采用,图1为AT&T人脸库里部分人脸图像。

我们采取与参考文献[6]一样的方法取每个对象的5张图像用于训练,剩余的5张图像用于测试,训练样本集和测试样本集均由200张图像构成。首先将训练样本和测试样本图像划分为3×3的block,然后分别对每个区域选取LBP8,2u2模式下提取LBP算子,再将这些区域的LBP算子组合成整幅图像的LBP特征,再利用RSRC分类器对测试样本进行分类,将实验结果与其他几种识别算法做了比较,如表1所示。为了比较RSRC分类器和其他分类器的分类效果,将人脸图像划分为不同的block。

从表1可以看出,几种不同方法在AT&T人脸库上均取得比较好的识别结果,我们的方法最后实验测得的识别率是97.5%,明显优于其他几类识别算法。可以得出我们的方法在人脸识别方面具有较强的识别能力的结论。

4 小结

加权稀疏表示 篇7

电阻抗成像(Electrical Impedance Tomography,EIT)技术的原理为通过安装在物体表面的电极对物体有规律地施加激励,若物体内部存在阻抗变化,则会引起表面测量电极电位的变化,基于测量电极电位变化,结合相应重建算法可以获得物体内部阻抗变化图像。EIT技术应用领域广泛,涉及肺功能成像,脑部功能成像,乳腺癌检测,腹部脏器功能成像等多个领域[1]。其中,颅脑电阻抗技术可以检测出脑功能性活动引起的阻抗变化,也可以检测出颅脑疾病如脑出血等引起的阻抗变化,具有诱人的应用前景[2]。

第四军医大学课题组长期致力于颅脑EIT研究,在国际上率先开展了颅脑EIT的临床实验研究,并取得了一系列的研究成果[3,4]。课题组研究发现相比实验室环境,临床实验中颅脑EIT面临更为严重的噪声与干扰。随机噪声、电极位移、接触阻抗的变化等多种因素均会影响到EIT的重建图像质量,故在临床实验背景下对重建算法的性能提出了更高的要求。在过去的几年中,围绕着提高重建图像质量和稳定性,国内外多个研究小组对重建算法开展了深入研究,在传统二阶范数加权算法(L2算法)的基础上提出了一系列优化算法。随着稀疏重建理论的发展,稀疏加权算法(L1算法)被应用到EIT,相比二阶范数加权算法,稀疏加权算法的抗噪声能力强,可以避免对重建图像边界造成过多模糊[5,6,7]。2009年Adler小组提出一种名为GREIT的优化算法并将其应用到肺部阻抗成像[8]。

基于简单仿真和部分实验研究,这些优化算法相对于传统算法,对EIT图像质量均有改善。然而,在颅脑EIT成像背景下针对这两种优化算法开展的算法比较工作尚未有文献报道。对此,本文基于三维颅脑仿真模型,模拟了EIT监测脑出血过程。我们在颅脑仿真模型不同位置设置了与血液电导率相同的球形目标,在添加一定的高斯噪声后,使用传统二阶范数加权算法、稀疏加权算法与GREIT算法进行图像重建,并基于图像噪声、位置误差、形变误差3项指标对重建结果进行了量化评价。根据重建图像与评价指标,我们总结了3种算法的不同特点并针对哪种算法更适合颅脑EIT进行了分析,为以后颅脑电阻抗算法的优选与优化打下了基础。

1 材料与方法

1.1 EIT算法原理

EIT算法可以分为正问题和逆问题求解。正问题指的是已知场域内部阻抗分布和边界条件,求解场域电势分布,其常用的求解方法为有限元方法[9]。逆问题又称作图像重建,则是给定边界测量数据和边界条件,求解场域内部阻抗分布。在EIT中,当场域内部阻抗变化较小时,边界测量数据与场域阻抗变化分布如下线性方程:

其中ρ表示场域电导率变化分布,v为边界测量数据,S是敏感系数矩阵。

由于EIT技术中测量数远远超过未知数,因此逆问题求解存在病态性,不能直接对(1)式求逆,采取的解决方法是引入正则项来构造如下的目标方程,然后求满足下式的最优解:

为阻抗变化分布的最优估计,λ为正则化参数,Φ(ρ)是正则化项。

1.2 常用的电阻抗重建算法

随着EIT算法的发展,形成了多种重建算法,不同算法具有不同的成像特点。本文选取了传统二阶范数加权算法、稀疏加权算法和GREIT算法3种常用的重建算法用于算法比较。

1.2.1 传统二阶范数加权算法

按照二阶范数正则化加权方法[10],电阻抗目标函数有如下形式:

对于这种求二阶范数加权最优化问题,常用的阻尼最小二乘解[11]的形式为:

W为正则化矩阵,可以控制解的一些特性。

1.2.2 稀疏加权算法

按照稀疏重建理论,其目标函数如下:

由于一阶范数的存在,不能直接对上式进行求导,针对处理不可导项产生了一系列的稀疏加权算法,本文选取的是成像较好、收敛速度较快的Split Bregman算法(SB算法)[12]。

根据Split Bregman算法[13],电导率分布可以按照以下3个步骤迭代求解:

其中(6)为二阶范数优化问题通过求导可得:

(7)为L1正则问题,可通过软阈值算子策略进行求解

其中Shrink为软阈值算子

1.2.3 GREIT算法

EIT线性算法一般满足如下形式:

其中x为所求电导率分布,y为边界电位变化,R为重构矩阵。为了寻取最优的重构矩阵R,GREIT算法首先定义一系列的训练目标t(i),然后寻求使下式误差最小的R作为优化的重构矩阵:

对于每一个训练目标,均可以生成与之对应的边界电位y和理想的重构值分布x。w是加权矩阵,代表训练目标集合t中每一个目标的权重。

对上式进行求导:

由上式,重构矩阵R为:

已知训练目标集合t与理想重构分布x、边界电位y满足如下关系:

D为训练目标集合与理想重构分布之间的映射矩阵,J为敏感系数矩阵,将上述两个式子带入

Σt*为训练目标集合的协方差,Σn为噪声的协方差,则可得重构矩阵R

由(19)可知,相比传统算法,GREIT算法的核心在于映射矩阵D。D的作用为使训练目标集合t可以映射到符合EIT特点且均匀性较好的理想重构分布x。

1.3 三维真实颅脑仿真模型构建

基于Mimics软件对颅脑CT进行分割和三维表面重建,然后使用Solid Works对三维表面模型实体转换,将实体化的模型导入Comsol4.4生成可用于颅脑EIT分析的有限元模型[14,15],见图1。构建出的模型结构完整,包括头皮层、颅骨、脑脊液、脑实质以及脑室,见图2。

分别对不同组织赋予不同的电导率[16],基于该模型,可以较为真实的模拟颅脑疾病及颅脑电阻抗技术应用过程。

1.4 仿真实验设计

1.4.1 真实颅脑边界重构模板

选取电极所在平面对三维模型进行切割,导出切割面。手动选取中心点与边界点,然后对颅脑切面图进行自适应剖分,生成重构模板,见图3。

1.4.2 图像评价指标

我们参考了Adler等在2009年提出的算法评价指标体系并使用图像噪声(Image Noise,IN)、位置误差(Location Error,LE)、形变误差(Shape Error,SE)3项指标对重建图像进行评价。首先,定义重建图像中单元像素值大于最大像素值四分之一的集合为感兴趣区域集合,记作ΩROI,感兴趣区域以外的区域为Ωother。

为重建图像非感兴趣区域灰度变化的平均值,为整幅图像灰度变化的平均值。

dPre为参考图像中预设目标距离中心点的距离,dRe为重建图像重建目标距离中心点的距离,l为颅脑长轴长度,LE越小位置误差越小。

本文形变误差的计算方法为重建图像感兴趣区域单元面积AΩROI相对于参考图像预设目标在电极平面的投影面积APre的变化。

2 结果

首先使用简单圆域模型,在理想条件下进行仿真。在Comsol中构建圆域模型,圆域半径为0.2,背景电导率为1,见图4(a)。在圆域1/2处设置电导率为1.2,半径为0.01的圆形目标,生成边界电位数据后分别使用L2加权算法、L1-SB算法、GREIT算法进行图像重建,见图4。

由图4简单圆域仿真结果可知,理想情况下,3种重建算法均能获得质量较好的重建图像,可以对目标准确定位。但由重建图像对应的扩散函数来看,3种算法重建特征存在较大差异。稀疏重建算法扩散函数在背景处几乎为0,而在目标边界处变化极为剧烈,对应的重建图像背景伪影少,重建目标边界明确;GREIT算法在目标内部扩散函数变化平缓,对应的重建目标较为均匀。

然后基于三维真实颅脑模型,进行仿真。分别在脑实质不同位置设置体积为5 cm3,电导率为0.7Ω/m的球形目标模拟脑出血,并在生成的重建数据中添加0.1%的高斯噪声,使用L2算法、稀疏加权算法、GREIT算法进行图像重建,见图5。

由图5重建结果可知,相比理想重构,基于颅脑三维仿真模型并添加高斯噪声生成的数据经重建后图像出现了扰动与伪影,尤其是在中心处的目标发生了较大的形变。对于3种算法的重建结果,传统L2加权算法获取的重建图像质量低于另外两种优化的算法;L1-SB算法对应的图像伪影最少,目标形变也较小;GRETI重建目标均匀性好。

与重建图像结果对应,在图像噪声、位置误差和形状误差三项指标上,稀疏加权算法的性能要优于其他两种算法,见图6。

3 讨论

本文在颅脑EIT应用背景下,基于三维颅脑仿真模型比较了传统二阶加权算法、稀疏加权算法与GREIT算法的性能。根据仿真实验结果可知,相比于传统的二阶范数加权算法,两种优化的算法对颅脑电阻抗图像重建均有一定的改善作用,但是基于不同的数学模型,两者的重建特性不同。稀疏重建算法的突出特点为抗噪性能好,可以保留边界,抑制由重建带来的形状误差,但是相比其他两种算法稀疏重建算法的需要调节的参数较多,在实际使用中增加了调参的负担[17];而GREIT算法在加权矩阵D的作用下其解分布比较均匀,噪声相对较少。对于颅脑EIT成像,相比均匀性,抑制伪影与减少目标形变对于临床应用更为重要,故我们认为稀疏加权重建算法更适用于颅脑EIT。

另外,我们的对比研究中存在几点需要进一步完善的地方。首先,本文所有的参数选择均为手动选择,基于该方法选取的参数可能不是最优。但只要保证3种算法进行重建时采用的参数相同,算法对比工作就不会受到本质影响。其次,本文在仿真模型上验证了算法的性能,下一步我们将开展物理模型实验、活体实验来验证算法性能。