干旱与半干旱地区的风沙跃移不仅驱动了沙漠和戈壁的演化^［1］，且对所在区域的基础设施、生产、交通及生命健康造成严重威胁^［2-4］，因而探索风沙颗粒的运动规律对风沙灾害的防治有着重要的意义。

为了研究风沙颗粒运动规律，前人提出了基于不同相关策略的算法，如松弛算法^［5-6］、Voronoi图法^［7-9］及卡尔曼滤波的系列算法等^［10-12］。在这些算法中，卡尔曼滤波系列算法准确率和召回率优于以往算法^［10-12］，但仍然存在着在中等沙粒浓度下（~100颗粒/帧）算法召回率低而影响追踪准确率的问题^［10-12］。除了上述算法外，人工关联的方法具有准确率高而不易做批处理的局限性 ^［13］。在卡尔曼滤波系列算法中^［10-12］，基于MATLAB平台的灰度阈值算法^［10］、YOLOv8算法^［11］及改进的YOLOv5算法^［12］等方案被用来识别沙粒，但这些方案在中等浓度下都存在召回率低的问题^［14］。同样地，用于风沙颗粒识别的掩码相关算法可能也存在类似的问题^［15］。YOLOv8算法^［11］及改进的YOLOv5沙粒算法^［12］虽然能显著地提升沙粒的识别率，但因为需要借助并行计算来实施而影响其可用性。

最近，高斯混合模型^［16］和支持向量机^［17］等轻量化模型被用于计算机视觉领域的目标追踪^［18-20］。考虑到风沙颗粒识别的复杂性和这两个模型的不确定性，本文构建了高斯混合模型初筛-支持向量机分类-沙粒统计特征再确认的新算法，以发展计算成本低廉且适合中低浓度环境的沙粒识别算法。

用于新算法验证的风沙图像来自一个细沙床面的风洞输沙实验，拍摄于陕西师范大学风沙动力学实验室的开口式风洞。所用沙粒的粒径在150~200 µm，拍摄时风洞的自由风速为6.8 m·s^-1，对应的摩阻风速为0.72 m·s^-1。以绿色激光片光源（约1×10^-3 m厚度）作为照明，使用配备微型镜头（AF Micro-Nikkor，60 mm f/2.8D）的数字高速摄像机（PHANTOM V711，由AMTEK 公司制造）进行拍摄。数字高速摄影机布置于风洞试验段外距离风洞中心轴面约270 mm处。设置拍摄频率为3 000 Hz，共获得57 652帧分辨率为512 像素（3.9×10^-2 m）×318像素（2.9×10^-2 m）的目标图像。更多的实验细节见文献［10］。

和以前研究相似^{［10-11， 21］}，这里分别选择极低浓度E（~4个颗粒/帧）、低浓度L1（~26个颗粒/帧）、低浓度L2（~30个颗粒/帧）、中浓度M1（~103个颗粒/帧）、中浓度M2（~108个颗粒/帧）和中浓度M3（~134个颗粒/帧）等6段视频进行新算法的验证。除E视频为1 000帧图像外，其余视频段均为300帧图像。

新算法包括高斯混合模型、支持向量机和沙粒特征的统计模型等方案。高斯混合模型通过假定被识别目标的灰度分布符合多模态高斯模式来识别目标^［16］，其敏感性高但易产生错误的识别，因而被用来对风沙颗粒进行初步识别。作为一个轻量化的机器学习模型^［17］，支持向量机可以在个人电脑上运行，且和传统模型如高斯混合模型相比较，具有解决非线性问题如沙粒识别的独特优势，故采用支持向量机对高斯混合模型的识别效果进行再分类以提高沙粒识别的可靠性。鉴于支持向量机分类结果依赖于训练数据的情况，为了增强算法的泛化性能，根据自动提取的人工标注沙粒的几何和颜色特征构建沙粒特征的统计模型。利用沙粒特征的统计模型对支持向量机的分类结果再确认以提高算法的可靠性。利用沙粒特征的方法对高斯混合模型识别出的目标进一步区分，以确保风沙颗粒识别的准确性。

新算法的实施流程如下（图1）：

a） 初始化背景模型：分别选择6个视频段（E，L1，L2，M1，M2，M3）的前25帧图像作为高斯混合模型的资源库，利用Python 平台上的内置函数CV2c［简称CV2c，CV2.createBackground Subtractor MOG2（history=25， varThreshold=2.6，参数设置来源于敏感性试验）］提取初始化视频中每个像素点对应的灰度均值、方差和权重；在此基础上，选择权重最高的正态分布模型作为该像素点的背景模型，其对应的灰度均值作为视频段中每一个像素点的背景均值。

b） 更新背景模型：计算第26帧图像中每一个像素点的灰度值与对应背景模型的灰度均值之差的绝对值；如果该绝对值大于2.6，则该像素点被确定为目标对应的像素点，否则确定为背景对应的像素点，并把该像素点对应的灰度值输入背景资源库，并用内置函数CV2c更新该像素点的灰度均值、方差和权重，并以最大权重对应的正态分布模型为更新的背景模型。

c） 后续图像分类：对后续图像中每个像素点分类按照a）、b）流程进行操作直至完成对所有图像的所有像素点的目标/背景进行归类。

d） 沙粒分割：利用Python 平台上的内置函数cv2m （简称cv2m， cv2.morphologyEx， 经试验，对风沙图像分割，这里的参数history 设置为25， varThreshold设置为2.6） 对c） 步识别的所有沙粒进行开运算（腐蚀-膨胀）及闭运算操作（膨胀-腐蚀）以便把沙粒从背景中分割出来。此过程需要根据风沙颗粒的特征，通过敏感性试验确定合适的开、闭运算中内核（kernel）的几何形状和尺寸（见1.3节）。

e） 优化支持向量机：采用基于随机梯度下降的单类支持向量机（Stochastic Gradient Descent OneClassSVM，核函数为线性核，模型中惩罚参数nu为0.5），选取L3视频段，人工标注10 529个沙粒作为数据集，以其中的70%作为训练样本而其余的30%用于支持向量机的性能测试。

f） 支持向量机对沙粒分类：把d）流程中识别的所有跃移颗粒输入到优化的支持向量机，以便对初始识别跃移颗粒进行再分类。

g） 构建沙粒特征统计模型：利用作者在Python 平台上的编写代码提取人工标注的L3视频段10529个沙粒的几何特征和颜色分布，据此构建沙粒几何和颜色特征的统计模型（见1.4节）。

h） 确认沙粒：利用g）步获得的统计模型对f）步分类结果进行再确认，即将符合概率超过1%的几何和颜色特征参数的沙粒确认为真实沙粒。

i） 新算法验证：对h）步获得的分类结果进行人工核对以验证新算法的表现，这里采用确定召回率、准确率和F1分数来评价。

式中：Ra，Ar， Nic，Ni分别是沙粒识别算法的召回率、准确率、沙粒识别算法识别的正确颗粒数及识别的总颗粒数。

新算法运行的硬件环境如下：在Python平台进行编程，运行设备为Aspire vx5-591g，处理器为Intel（R） Core（TM） i5-7300HQ CPU @2.50 GHz。

为了获得运行高斯混合模型所需要的沙粒分割的最佳开闭运算的内核形状及大小，分别选择了十字形、方形和椭圆形等3类内核，对应的开运算核大小分别选择为5×5像素²（记作（5，5））、5×3像素²（记作（5，3））、3×3像素²（记作（3，3））；对应的闭运算核大小可分别选择为（ 5，5 ）、（ 3，3 ）、（ 5，3 ）。通过组合形成开闭运算核尺度的4个组合，其分别表示为：（5，5）（ 3，3 ）、（5，3）（ 5，3 ）、（5，5）（ 5，5 ）和（3，3）（ 3，3 ）。这里以极稀浓度下的1 000帧图像中共计4 263个沙粒进行敏感性试验，以确定合适的沙粒分割中开闭运算对应的内核形状及大小（图2）。

图2为沙粒分割中不同开闭运算内核形状和大小组合下沙粒检测的结果。由此可知，核大小为（5，5）（ 3，3 ）和（5，5）（ 5，5 ）时检测结果相似，故仅考虑核大小为（5，5）（ 3，3 ）的情况。当开闭运算核形状为十字形、核大小分别为（5，5）（ 3，3 ）时，高斯混合模型检测到的沙粒数量为4 106个，与样本数量最接近。而当核大小分别为（3，3）（ 3，3 ）时，高斯混合模型则出现严重的高估。当开闭运算核形状为十字形、核大小分别为（5，3）（ 5，3 ）时，高斯混合模型检测到沙粒数量略高于样本数量；当开闭运算核形状为方形且与上述不同核大小组合时，高斯混合模型对沙粒数量的检测存在显著的低估或高估。当开闭运算核形状为椭圆形时，高斯混合模型也会出现与方形核相似的情况。据此，当开闭运算核形状为十字形、核大小分别可以设置为（5，5）（ 3，3 ）或者（5，3）（ 5，3 ）；当开闭运算核为方形，核大小可设置为（3，3）（ 3，3 ）。

考虑到高斯混合模型的表现不仅受沙粒分割过程中开闭运算核形状及核大小的影响，也与沙粒浓度有关，这里进一步根据高斯模型检测沙粒的F1分数随沙粒浓度的变化来进一步确定开闭运算合适的核形状及核大小。高斯模型检测沙粒的F1分数随着沙粒浓度（50%~85%）剧烈波动，这意味着开闭运算核形状及核大小设置时，应根据沙粒的浓度而调整。该设置对识别的影响可能是由于极稀浓度下颗粒几乎没有表观重叠，故选用方形核可以最大程度保留沙粒信息；随着浓度的升高，表观重叠现象增多，选用较大的十字形核能有效避免颗粒重叠对沙粒分割的干扰。总之，设置极稀浓度E下开闭运算核形状为方形，大小为（3，3）（ 3，3 ）；在低浓度（~20颗粒/帧）及中浓度（~100颗粒/帧）下，核形状为十字形，核大小分别为（5，3）（ 5，3 ）。

本文选取人工标注的低浓度L3数据集10 529个沙粒用作者代码对沙粒的几何和颜色特征进行分别提取，以建立了相应的统计模型。其中，几何特征如图3所示。可以发现，沙粒宽度在2~18像素、沙粒高度在3~8像素、沙粒周长在10~50像素、沙粒面积在0~140像素。同时满足上述4个参数要求，即被几何统计模型确认为真实的沙粒。

需要说明的是，这里选择1%作为筛选真实沙粒的依据，具有一定的经验性。即1%以下的目标可能是噪音或其他的环境干扰。但后续结果表明，此1%的筛选能轻微地提高沙粒识别的准确率，因而统计模型的筛选是有意义的。

同样地，可以分别建立基于RGB（Red Green Blue）的沙粒三色概率分布模型和基于HSV模型（Hue Saturation Value）的概率分布模型。在这两个模型中，仍然选择概率超过1%的参数作为认识沙粒的依据。具体颜色参数见表1。

图4为新算法与已有算法的召回率和准确率。和Matlab平台的灰度阈值分割算法^［10］与动态灰度阈值分割算法^［14］相比较，新算法在低浓度下的召回率分别提高了30%和40%，而在中等浓度下分别提高了30%和25%，这意味着新算法在提高沙粒识别的召回率方面比灰度阈值方案有明显的优势。和YOLOv8相比较^［11］，新算法在极稀浓度和低浓度情景下与其相似而在中等浓度下比该模型提升10%以上，这进一步显示了高斯混合模型在提高沙粒识别的召回率的作用。另一方面，和改进的YOLOv5相比较^［12］， 新算法在极稀浓度和低浓度情景下召回率高10%而在中等浓度下低20%，这意味着改进深度神经网络在中等浓度下对提升沙粒召回率有优势。另外，由于新算法涉及的支持向量机分类和沙粒特征的统计模型再确认等环节，和高斯混合模型相比较，新算法的召回率分别下降了0.7%-5.5%；和高斯混合模型-支持向量机相比较，新算法召回率下降了0.2%~3.9%。

和动态灰度阈值分割算法^［14］相比较，新算法准确率提升了50%以上，这意味着新算法包含的支持向量机分类和沙粒统计模型的确认环节有效地提升了算法的可靠性。不仅如此，和改进的YOLOv5相比较^［12］， 新算法的准确率也提升了5%~20%，这表明新算法包含的支持向量机分类和沙粒统计模型的确认环节在提升沙粒识别准确率方面扮演着重要的角色。和高斯混合模型相比较，新算法准确率提高了1.2%~17.4%；与高斯混合模型-支持向量机相比较，除M2浓度下准确率降低1.2%外，其余浓度下准确率提高了0.5%~3.7%。

综上，在目前所讨论的跃移沙粒识别算法中，新算法具有中等召回率，其数值变化在60%~90%，其准确率保持在80%~95%，考虑到目前所有算法的召回率和准确率强烈地依赖于沙粒浓度的特点^［10-12］，而新算法的训练测试是在低浓度和中等浓度下进行的，因而我们认为其主要适于中、低浓度环境下的风沙颗粒识别。至于高浓度下环境下的问题仍需进一步研究。

表2比较了新算法与发表的风沙颗粒识别算法的推理时间。这里的推理时间是指处理每帧图像所需要的用时。YOLOv5系列^［12］是在网上租赁GPU进行图像处理的，因而其推理时间仅有新算法的一半左右。但使用GPU平台意味着高的经济成本。而新算法可直接在个人电脑上运行，显示了其和复杂神经网络的显著优势。和能在个人电脑上运行的YOLOv8模型相比较，其推理时间仅是YOLOv8算法的十分之一。和高斯混合模型、高斯混合模型-支持向量机相比较，其运行时间成本略高，但考虑到它在准确率等方面的优势，略高的时间成本是可以接受的。

针对目前风沙颗粒识别算法准确率低或算法运行成本高的特点，本文提出了一个集成高斯混合模型、支持向量机及统计模型的风沙颗粒识别算法。该集成算法包括高斯混合模型分割、支持向量机分类及沙粒统计模型确认等环节。风沙图像检测结果表明，新的集成算法具有中等召回率和最高的准确率和较低的时间成本，是一个轻量化的适合中低浓度的风沙图像识别算法。目前工作为风沙图像识别与追踪提供了新的手段。