引　言

水源微生物广泛分布于自然与人工水体环境中，它们无法用肉眼直接观察，例如水中的原生动物和真菌的大小通常介于0.1～100 μm，只能在显微镜下进行观察。作为水质的关键评价指标之一，水源微生物的繁殖情况可以很好地反映水体的污染程度。因此，推进水源微生物分类识别的相关研究对于水源的生物安全性监测、水环境治理等具有非常重要的意义。

考虑到传统机器学习方法对图形处理单元的要求较低，能够以较低成本部署在设备中，且所需的样本量小，训练时间短，许多学者仍在利用机器学习方法来提高微生物检测的准确性。Rahmayuna等^[1]提出了一种利用支持向量机实现细菌种类的属级分类技术，对4种细菌的预测准确率达到了90.33%。Dhindsa等^[2]利用像素聚类和Kirsch滤波分割图像中的微生物，对不同分类器提取相应的特征，再比较分类回归树等多种分类算法，结果表明，分类回归树算法最优，对10种微生物的识别准确率达到了98%。本研究通过边缘检测来提取图像中的微生物，并提取多种图像特征并优化了部分特征参数，基于多特征融合的方法制作微生物数据集，采用粒子群优化算法的支持向量机（support vector machine optimized by particle swarm optimization，PSO-SVM）来训练识别模型，并将其与网格搜寻支持向量机（grid search support vector machine, GS-SVM）^[3]、误差反向传播算法（back propagation neural network, BP）和集成支持向量机（ensemble SVM）^[4]比较，建立了一种高效的微生物识别方法。

1 水样采集与微生物图像集制作

微生物采集于浦东威立雅自来水公司，由于夏季是微生物繁殖的高峰期，所以选择该时期进行采样工作，采样类型包括过滤水和出厂水。首先将滤袋放置在出水管道上，持续富集3 d，之后将采样带回实验室并转移到培养皿中，在Axio Vert.A1型倒置显微镜（德国卡尔·蔡司公司）下观察，照明光源为显微镜自带的明场光源，显微物镜放大倍率为5×。使用Canon EOS 80D相机（日本佳能株式会社）对微生物进行拍照，ISO速度为800，曝光时间为1/200 s，图像分辨率为4000×6000。最终得到8种微生物，共640张图像（每种80张），图1展示了这8种微生物的实例。采用每种微生物的50张图像用作训练集，30张图像用作测试集，为减小计算量，提高后续的处理速率，图像尺寸统一调整为600×900。用Sefexa软件来绘制评价标准图像。

2 水源微生物识别系统整体流程

提出了一种水源微生物自动识别系统，该系统的整体流程如图2所示。首先，输入制作完成的微生物图像集，通过自编写的图像分割算法对微生物图像执行图像增强、平滑滤波、边缘检测、形态学处理等操作，获得目标微生物的二值掩膜，并提取目标微生物所在区域，得到分割后的微生物图像。随后，对分割后的图像提取6种特征描述子，包括几何特征、内部结构直方图、傅里叶描述子、胡不变矩、灰度共生矩阵和旋转不变局部二值模式。由于6种特征的总维数不高，通过连接的方式进行特征融合。接下来，导入融合后的特征数据，通过粒子群算法优化SVM训练参数，建立优化后的微生物识别模型，最终由该模型输出微生物的识别结果。

3 微生物图像分割与特征提取 3.1 微生物图像分割

设计了两种微生物图像分割方法，分别是全自动图像分割方法与半自动图像分割方法。半自动分割方法是在全自动分割方法基础上，新增了可以手动画出感兴趣的分割区域的功能，并只在该区域内进行后续的图像分割，这样可以在复杂环境中更加准确地分割出目标微生物。针对不同的环境复杂程度可采用不同的分割方法，例如：在流动的水体中，微生物和杂质密度较低，适合全自动分割方法；在富集的水样中，目标微生物会被大量杂质覆盖，适合半自动分割方法。

对于全自动分割方法，首先输入原始微生物图像，再将原始图像转化为灰度图，之后通过点运算增强对比度，高斯滤波去噪，随后用Canny边缘检测算子^[5]来检测微生物边界，通过形态学处理闭合边缘图像并获得最终的分割结果。对于半自动分割方法，其流程如图3所示：第1步，读取原始图像；第2步，用户通过光标指针手动画出感兴趣的区域；第3步，将原始图像转化为灰度图，高斯平滑处理，点运算增强对比度；第4步，裁剪出所画的感兴趣区域，并求出其他区域的灰度均值，用来填充黑色背景区域；第5步，用Canny双阈值边缘检测提取微生物边缘轮廓，并使用形态学处理方法闭合边界；最后，填充微生物区域内部以获得分割后的微生物图像。

Canny边缘检测是最常用的边缘检测方法之一，它有着良好的抗噪声干扰性，对图像的弱边缘检测能力较好。基于梯度理论^[6]，假定f（x, y）代表微生物图像中坐标（x, y）处的灰度值。则该坐标处的梯度被定义为

$ \begin{array}{c}\nabla f=\left({G}_{x},{G}_{y}\right)=\left(\dfrac{{\partial} f}{{\text{∂}} x},\dfrac{{\partial} f}{{\text{∂}} y}\right) \end{array} $

(1)

相应的，梯度的幅值大小为

$ \begin{array}{c}M\left(\nabla f\right)=\sqrt{{G}_{x}^{2}+{G}_{y}^{2}}=\sqrt{{\left(\dfrac{{\partial} f}{{\partial} x}\right)}^{2}+{\left(\dfrac{{\partial} f}{{\partial} y}\right)}^{2}} \end{array} $

(2)

该坐标处的梯度方向为

$ \alpha(x, y)=\arctan \left(\dfrac{G_y}{G_x}\right) $

(3)

式中，G_x和G_y分别为横轴x方向和纵轴y方向的梯度分量。根据式（2）和（3），可以得到整幅图像所有像素点的梯度幅值和方向。随后Canny算子在0°、45°、90°和135° 4个梯度方向上对数据执行非极大值抑制，对于中心像素点（x, y），选定其周围的8邻域像素范围，得到该点的梯度方向所在直线与邻域相交的两个端点。由于实际图像是离散二维矩阵，端点不一定真实存在，需要进行插值。将点（x,y）与两个端点比较，若中心点梯度幅值是3个点中的最大值则保留，反之置0。随后设定适当的低阈值与和高阈值，若像素点梯度幅值小于低阈值则舍弃；若梯度幅值大于高阈值则保留；若梯度幅值介于两者之间，则求出该点的8邻域像素的梯度幅值，如果梯度幅值高于高阈值，则保留该点，否则舍弃。

3.2 图像分割评价

对于以上两种图像分割方法，采用相似度（similarity）、灵敏度（sensitivity）以及特异度（specificity）3个指标进行评价，评价时的分割图像与评价标准图像均为二值图像，1代表感兴趣像素，0代表背景像素。其相似度式定义分别为

$ \begin{array}{c}{V}_{\mathrm{sim}}=\dfrac{{T}_{{\rm{p}}}}{{T}_{{\rm{p}}}+{F}_{{\rm{p}}}+{F}_{{\rm{n}}}} \end{array} $

(4)

式中：V_sim表示分割图像与评价标准图像的相似度；T_p表示被正确标记的像素总数；F_p 表示被错误标记为感兴趣像素的背景像素总数；F_n 表示被错误标记为背景像素的感兴趣像素总数。

$ \begin{array}{c}{V}_{\mathrm{sen}}=\dfrac{{T}_{{\rm{p}}}}{{T}_{{\rm{p}}}+{F}_{{\rm{n}}}}\ \end{array} $

(5)

式中：V_sen表示灵敏度，即分割图像被正确标记的感兴趣像素总数与评价标准图像中的实际感兴趣像素总数之比。

$ \begin{array}{c}{V}_{\mathrm{spe}}=\dfrac{{T}_{{\rm{n}}}}{{T}_{{\rm{n}}}+{F}_{{\rm{p}}}} \end{array} $

(6)

式中：V_spe表示特异度，即分割图像被正确标记的背景像素总数与实际背景像素总数之比；T_n表示被正确标记的背景像素总数。

3.3 微生物特征提取

从分割后的微生物图像中提取几何特征、内部结构直方图特征、傅里叶描述子、胡不变矩、灰度共生矩阵和旋转不变局部二值模式，共6种图像特征。

几何特征（geometrical feature, GF）包括：（1）微生物的面积A；（2）微生物边界轮廓的周长P；（3）微生物的圆度（circularity）4πA/P²；（4）与微生物区域具有相同归一化二阶中心矩的椭圆长轴长度l₁和短轴长度l₂。

内部结构直方图（internal structure histogram, ISH）特征：如图4所示，在微生物边界轮廓上标记K个等距分布的采样点，将任意3个点组合，计算它们的内部结构角（例如图4所示的角θ），并将所有能得到的内切角按照一定的角度范围进行计数，由此可以得到ISH特征。

傅里叶描述子（Fourier descriptor, FD）：在微生物边界轮廓上标记M个等距分布的采样点，计算所有采样点对应的中心点坐标，并求出每个采样点到中心点的欧几里德距离，按顺序排列为{d₀, d₁, ···, d_m}，并通过离散傅里叶变换提取它们的频域分量作为傅里叶描述子

$ a_s=\dfrac{1}{M} \sum_{t=0}^{M-1} d_s \exp \left(-\dfrac{{\rm{j}} 2 \pi s t}{M}\right) $

(7)

式中：a_s表示第s个复系数，s = 0, 1, ···, M−1；j表示单位虚部；t表示第t个采样点到中心点的欧几里德距离。为了降低噪声、平移和旋转对系数a_s的影响，对其采取内部标准化

$ \begin{array}{c}{b}_{s}=\dfrac{{a}_{s}}{{a}_{0}} \end{array} $

(8)

由于微生物边界轮廓可表示为实函数，所以b_s有共轭对称性，最终得到一组M/2维特征向量（b₀, b₁, ···, $ {b}_{\frac{M}{2}-1} $ ）。

胡不变矩（Hu moments, Hu）：Hu是一种具有平移、旋转、尺度不变性的特征描述子^[7]。在离散条件下，Chen^[8]对其进行了改进，得到8个不变矩特征，使得提取的图像形状特征更加全面。前7个参数在以往的文献中已提及，不再赘述，第8个不变矩计算式为

$ \begin{split} \varPhi_8=& 2 \eta_{11}\left[\left(\eta_{30}+\eta_{12}\right)^2-\left(\eta_{03}+\eta_{21}\right)^2\right]-\\ & 2\left(\eta_{20}-\eta_{02}\right)\left(\eta_{30}+\eta_{12}\right)\left(\eta_{21}+\eta_{03}\right) \end{split} $

(9)

式中， $ {\eta }_{pq} $ 表示图像的p+q阶归一化中心矩。

灰度共生矩阵（gray level co-occurrence matrix, GLCM）：GLCM是一种纹理特征，它能反映图像中任意两像素点的灰度值的空间相关性^[9]。为了减小计算量，把原始图像的灰度级压缩并量化为16级。在0°、45°、90°和135° 4个方向上计算共生矩阵，对它们分别归一化处理，再对每个方向的共生矩阵计算6个统计量，即能量、熵、对比度、惯性矩、相关度以及逆差分矩，最终得到一组24维的特征向量。

旋转不变局部二值模式（rotation invariant local binary pattern, RI-LBP）：RI-LBP特征是原始LBP特征的一种改进型，具有计算简单、可适应不同尺度、灰度和旋转不变性的优点^[10]。通过将灰度图像分块，然后在每个块内的每一个像素点的圆形邻域内进行编码采样，之后对采样得到的LBP值进行循环移位操作，得到若干个不同的值，取其中最小的一个作为最终的LBP值，最后对每个块的LBP值进行直方图统计得到旋转不变LBP特征。

假设以灰度图的某一点（x_c, y_c）为中心，在其半径为R的圆形邻域内，采样点个数为Q，则RI-LBP的计算方法如下

$ L B P_{Q, R}\left(x_c, y_c\right)=\sum_{q=1}^Q s(I(q)-I(c)) 2^q $

(10)

$ \begin{split} R I- L B P_{Q, R}^{r o t}&=\min \{\operatorname{ROR}\left(L B P_{Q, R}, i\right) \\ i& =0, \cdots, P-1 \} \end{split} $

(11)

式中：ROR表示沿顺时针方向旋转LBP算子i次；I（q）和I（c）分别表示第q个采样点和中心点（x_c, y_c）的灰度值；s（ž）为门限函数。由于采样区域是圆形，采样点可能并非实际存在的点，要通过线性插值来完成采样点的估计。

3.4 特征参数优化与特征融合

调整3种特征的参数，分别是ISH的采样点个数和角度统计范围、FD的采样点个数、RI-LBP的邻域半径与采样点个数。随后用SVM训练不同参数条件下提取的所有微生物特征数据并测试，为了降低结果的偶然性，每组实验重复10次，每次实验的训练集和测试集均为随机取样，以10次实验的平均准确率作为评价标准，选择准确率最优的特征参数用于后续实验。准确率的定义为^[11]

$ V_{\mathrm{acc}}=\dfrac{T_{\mathrm{p}}+T_{\mathrm{n}}}{T_{\mathrm{p}}+F_{\mathrm{p}}+F_{\mathrm{n}}+T_{\mathrm{n}}} $

(12)

式中：T_p表示把实际为正样本的微生物分类为正；T_n表示把实际为负样本的微生物分类为负；F_p表示把实际为负样本的微生物分类为正；F_n表示把实际为正样本的微生物分类为负。

之后，提取参数优化后的特征，将所有6种特征通过连接的方式融合起来。

4 PSO-SVM水源微生物识别方法 4.1 粒子群优化算法

粒子群优化算法（PSO）是通过模拟鸟类和鱼类种群的觅食行为而发展起来的一种基于群体协作的全局随机搜索算法^[12]。其原理是通过随机初始化粒子，然后迭代更新粒子的速度和位置来搜寻最优解，每个粒子沿着个体最优与全局最优的方向移动，在多维空间的搜寻过程中，每个粒子更新各自的速度和位置的式分别为

$ \begin{split} v_{i d}^k=& w v_{i d}^{k-1}+c_1 {rand}_1\left(p_{i d}^*-x_{i d}^{k-1}\right)+\\ & c_2 {rand}_2\left(g_d^*-x_{i d}^{k-1}\right) \end{split} $

(13)

$ \begin{array}{c}{x}_{id}^{k}={x}_{id}^{k-1}+{v}_{id}^{k} \end{array} $

(14)

式中：id表示在第d维空间中的粒子i；k为迭代次数；w为惯性权重；c₁和c₂为两个加速度常量；rand₁与rand₂是在（0, 1）中均匀分布的两个随机数；p*与g*分别表示个体最优位置与全局最优位置。

4.2 支持向量机

支持向量机（SVM）是一种基于统计学习理论的监督学习算法。它在求解过程中加入正则项来优化结构风险，适合处理小样本问题，广泛应用于分类和回归分析等领域^[13]。因此，对于本研究中有限的微生物图像数据，选择SVM来进行分类识别，由于径向基函数（radial basis function，RBF）具有普适性^[14]，使用它作为核函数。SVM的决策函数可表示为^[15-16]

$ f(x)=\operatorname{sign}\left(\sum_{i=1}^N \alpha_i y_i K\left(x, x_i\right)\right)+b $

(15)

式中：α_i为拉格朗日乘子；N为用于模型训练的样本数；b为常数偏置项；x_i为训练样本的特征数据；y_i为训练样本的类标签；K为RBF核函数

$K\left(x, x_i\right)=\exp \left(-\gamma\left\|x-x_i\right\|^2\right) $

(16)

式中，γ（γ > 0）为核参数。上述优化问题可以用拉格朗日泛函表述为

$\begin{split} & \max _\alpha\left(\sum_{i=1}^N \alpha_i-\frac{1}{2} \sum_{i=1}^N \sum_{j=1}^N \alpha_i \alpha_j y_i y_j K\left({\boldsymbol{x}}_i, {\boldsymbol{x}}_j\right)\right), \\ & \alpha=(\alpha _1, \cdots, \alpha _i, \cdots, \alpha _N) ; \sum_{i=1}^N \alpha_i y_i=0 ; 0 {\text{≤}} \alpha_i {\text{≤}} C \end{split}$

(17)

式中：y_j表示序号为j的单独样本的标签；x_i和x_j分别是对应于α_i和α_j的支持向量；C是惩罚因子。

4.3 粒子群优化算法的支持向量机（PSO-SVM）

在此项研究中，使用PSO算法优化支持向量机的核参数γ以及惩罚因子C，建立PSO-SVM识别模型。其流程如图5所示，具体步骤为：（1）输入微生物的特征数据集；（2）PSO算法参数设定，包括粒子群的粒子总数、最大迭代次数、惯性权重、局部和全局学习因子、待优化参数C和γ的限制范围；（3）初始化粒子群，随机初始化粒子的速度和位置，迭代次数置零；（4）循环执行步骤6到步骤9，每循环一次则迭代次数加1；（5）利用PSO优化的SVM参数C和γ建立微生物图像特征的SVM分类模型，并评估测试集的识别准确率；（6）计算每个粒子的适应度；（7）计算粒子的个体和全局优化适应度，并根据式（13）和式（14）更新每个粒子的速度和位置；（8）判断循环是否终止，如果达到两代之间的最小误差或者最大迭代次数，则终止循环，反之则返回步骤5；（9）输出最优的参数C和γ；（10）建立优化的SVM模型；（11）输出微生物分类识别结果。

5 结果与讨论 5.1 特征参数优化

为了获得客观的评估结果，使用评价标准图像进行特征参数的优化。表1展示了在不同采样点数量以及不同角度统计范围下提取的ISH特征的识别准确率，结果表明，采样点为100个，角度分区范围为{0～15°, ···, 165°～180°}时的准确率最高，达到了79.41%，此时的特征维数为12维。表2展示了在不同采样点数量下提取的FD特征的识别准确率，结果表明，当采样点数量为150个时，准确率最高，达到88.50%，此时的特征维数为75维。表3展示了不同邻域半径和采样点个数条件下RI-LBP特征的识别准确率，当邻域半径为1，采样点为8个时，准确率最高，达到72.16%，此时特征维数为256维。

5.2 图像分割评价

图6展示了两种分割方法的实际分割效果以及评价标准图，第1行是原始图像，第2行是评价标准图像，第3行是全自动分割图像，第4行是半自动分割图像。原始图像中目标微生物的周围有不同程度的杂质干扰，通过半自动分割方法可以将目标微生物非常精准地分割出来，分割后的图像形态与评价标准图像非常接近，而全自动分割方法无法去除与目标微生物紧挨着的大体积杂质，分割效果较差。

由式（4）—（6）可计算得到全自动分割图像与半自动分割图像相较于评价标准图像的相似度、灵敏度和特异度。评价结果如表4所示：半自动分割图像与评价标准图的相似度为90.94%，比全自动分割方法提升了4.96%；半自动分割图像的灵敏度为97.24%，比全自动分割方法提升了0.70%；半自动分割图像的特异度为99.89%，比全自动分割方法提升了0.27%。因此，半自动分割方法可以更好地分割微生物图像。

为了进一步评估分割性能，将以上两种不同分割方法分割后的微生物图像制成特征数据集，利用SVM训练和测试，求出单类别特征的识别准确率，训练过程采用5折交叉验证^[17]，实验重复10次，取平均值。结果如图7所示，除RI-LBP特征以外，半自动分割方法的其他所有特征的识别准确率均高于全自动方法，其中GF的识别结果最优，其准确率为90.84%。与评价标准图相比，半自动分割图像的GF、Hu、GLCM 3种特征的识别准确率差异均在1%以内，ISH的准确率差异约3%，RI-LBP的准确率差异约4%。此外，对于半自动分割图像，其所有6种特征的平均准确率为77.35%，与评价标准图的相应值相比仅低了2.62%；全自动分割方法的所有特征平均准确率与评价标准图的相应值比较，低了7.35%。因此，对拍摄于富集后的水体微生物图像，适用于半自动分割方法。因为大量的微生物与杂质堆积在一起，分割难度非常大，而通过手动选择感兴趣区域，可以有效降低图像噪声的影响，实现对目标微生物的精准分割。对于拍摄于流动水体的微生物图像，则适合全自动图像分割方法，因为该方法无需人工干预，能实时处理图像，且微生物和杂质的密度较低，大多数微生物能被正确地分割。

5.3 PSO-SVM微生物识别模型

评价标准图提取的微生物的6种特征数据，经过特征融合后作为特征数据集。利用PSO算法优化SVM的核参数γ以及惩罚因子C，得到PSO-SVM微生物识别模型。PSO算法的参数设置如下：粒子的种群大小为24；加速常数C₁与C₂均为1.5；惯性权重为0.6；最大迭代次数为100。在PSO优化过程中，使用SVM模型的识别准确率作为粒子适应度。

利用PSO优化SVM参数过程中的适应度曲线如图8所示，适应度随着粒子的迭代次数增加而逐渐收敛，在第72代时达到最大值，由此获得最优的SVM参数C = 21.1652，γ = 0.025504。随后利用训练完成的PSO-SVM模型对8种微生物的测试集进行分类识别，采用5折交叉验证，测试10次，取平均值作为最终测试结果，并将其与GS-SVM、BP、Ensemble SVM 3种识别方法的测试结果比较。

结果如图9所示，PSO-SVM识别方法对于此8种微生物具有最好的分类性能，平均识别率为98.60%，其中：月形腔轮虫的识别率为99.19%；溞的识别率为98.13%；猛水蚤的识别率为98.77%；钩状狭甲轮虫的识别率为97.17%；未知微生物1的识别率为97.92%；未知微生物2的识别率为97.60%；红虫和颤藻的识别率最高，均为100%。GS-SVM模型的平均识别率为96.28%，比PSO-SVM低了2.32%；Ensemble SVM的识别率为91.57%，比PSO-SVM低了7.03%；BP神经网络的平均识别率为81.33%，在4种算法中最低，这是因为微生物样本数量较少所致。

实验结果表明，PSO-SVM识别方法的识别性能优于其他3种分类算法，更适用于水源微生物识别。利用全自动分割图像提取的特征数据建立的PSO-SVM模型平均识别率为92.25%，相较于评价标准图低了6.35%，这是由于本文中的原始微生物图像均采集于富集后的水样，其微生物与杂质比较密集，微生物图像的全自动分割效果较差所致。利用半自动分割图像提取的特征数据建立的PSO-SVM模型，其平均识别率为97.08%，仅比评价标准图的识别率低1.52%。

6 结　论

提出了一种基于多特征融合和PSO-SVM的水源微生物自动识别方法。考虑到水源环境的不同复杂程度，设计了全自动和半自动分割两种微生物图像分割方法，分别适用于流动水体中拍摄的微生物图像以及水体富集后拍摄的微生物图像。针对8种微生物优化了ISH、FD、RI-LBP的特征参数，通过连接的方式将所有特征融合起来，最后利用该特征数据建立了PSO算法优化的SVM模型，实现对8种水源微生物的分类和识别。实验结果表明，该模型能够准确地识别其中的每一种微生物，平均识别准确率达到97.08%。因此，该技术有望用于水源微生物的高效自动检测，以提高居民用水的安全保障和水源环境监测效率。