0.   引言

镜质组反射率与挥发分产率、粘结系数等指标有较强相关性，是公认表征煤化程度（煤级）的重要指标。其中，油浸物镜下测得的镜质组平均最大反射率(Mean Maximum Vitrinite Reflectance, MMVR)应用最为广泛，是进行炼焦配煤、混煤鉴别的主要依据^[1-2]。目前MMVR测定主要有2种方式：① 基于显微光度计的人工测定。该方式准确度高，但耗时长，测定过程受人员的专业知识影响较大，测试结果主观性较强。② 基于计算机视觉技术的自动测定。该方式凭借其快速、准确和可重复性等优势，有望成为MMVR测定的主要方向^[3]。

MMVR自动测定主要包括煤岩显微组分识别（镜质组的识别）和MMVR估计2个部分^[4]。近年来基于计算机视觉技术的煤岩显微组分识别引起了研究者的极大关注^[5]。M. Mlynarczuk等^[6]采用形态学梯度、灰度特征和最近邻算法实现了三大显微组分和矿物的识别，平均准确率为97.23%。Wang Hongdong等^[7]提出了一种基于图像分析的显微组分识别方法，开发了基于图像分割和分类的显微组分识别策略，可识别7种显微组分或亚组分，准确率为90.44%。 Lei Meng等^[8]设计了具有增强注意力门的改进U−Net模型并评估了各种编码器的性能，在89张显微图像上的实验结果显示该方法的分割准确率达到91.56%。王培珍等^[9]根据各显微组分的纹理特点和亮度差异，提出了迁移学习的镜质组分类方法，实现了对结构镜质体、无结构镜质体和镜屑体3种显微组分的识别。

虽然前人对显微组分的智能识别做了大量尝试，但对MMVR估计的研究相对有限。B. M. England等^[10]利用图像分析仪测量MMVR的分布，获得了较好的效果。王洪栋^[11]提出一种基于机器学习的MMVR估计方法，并发布首款免费公开的镜质组反射率分析软件，其研究结果证明机器学习方法在估计MMVR上具有较大的潜力。

上述研究成果极大推动了基于计算机视觉技术的MMVR测定研究，但仍然存在一些问题和改进的空间。目前的研究仍多为半自动化方式，即图像中镜质组区域的识别仍然依赖于人工判定。虽然Wang Hongdong等^[12]实现了全自动MMVR测定，但由于未考虑镜质组和非镜质组样本的不均衡问题，镜质组识别率有待进一步提升。此外，MMVR回归多采用单一灰度特征，回归分析性能仍有一定提升空间。鉴此，本文提出了一种基于机器学习的MMVR估计方法。首先，使用K−Means算法将煤岩显微图像中的不同显微组分聚类，利用聚类后的标签将煤岩显微图像分割成不同显微组分。其次，采用随机森林(Random Forest，RF)对分割后的显微组分分类，得到镜质组区域的显微图像。然后，采用人工少数类过采样法(Synthetic Minority Over-Sampling Technique, SMOTE)对少数类样本进行过采样，以解决分类时镜质组区域与非镜质组区域样本不均衡的问题。最后，使用基于树突网络(Dendrite Net，DDNet)的回归算法估计MMVR。本文的主要创新点如下：① 通过引入SMOTE算法，解决了镜质组和非镜质组区域样本不均衡问题，提高了镜质组识别的准确率，为MMVR回归分析打下良好基础。② 采用DDNet进行MMVR回归分析，提升了MMVR估计的精度和鲁棒性。③ 将显微图像分割、镜质组识别和MMVR估计算法进行集成，开发了MMVR估计软件。

1.   研究方法

基于随机森林和树突网络的煤镜质组反射率估计方法主要由基于K−Means的煤岩显微图像分割、基于SMOTE和RF的镜质组区域识别和基于DDNet的MMVR估计3个部分组成。MMVR估计流程如图1所示。

图  1  MMVR估计流程

Figure  1.  Estimation process of mean maximum vitrinite reflectance

下载: 全尺寸图片 幻灯片

采用K−Means算法将显微图像中不同显微组分聚类，将显微图像分割成K个不同的区域，每个区域均对应一种显微组分。由于每张显微图像中显微组分的类别数不同，本文采用手肘法确定K值^[13-14]，分别采用不同的K值进行聚类，并计算其误差平方和（每个样本点与其聚类中心距离的平方和），误差平方和与K的关系图是一个手肘形，取肘部位置对应的K值作为最优类别数。从分割后的显微图像中提取包括纹理、灰度和几何特征在内的综合特征进行镜质组区域识别。考虑到较大的镜质组区域能更真实地反映煤岩的镜质组反射率，本文仅对面积大于100×100像素的区域进行识别。在MMVR估计阶段，从识别出的镜质组区域随机截取41×41像素大小的方形镜质组窗口，从中提取灰度特征，并建立基于DDNet的MMVR回归模型。

1.1   SMOTE算法

考虑到显微图像中非镜质组区域较多，为避免类别不平衡对模型性能的影响，本文使用SMOTE算法，通过合成新样本的方式对少数类样本进行过采样^[15]。与随机过采样不同，SMOTE算法是在特征空间进行采样，并非数据空间。该算法通过生成人工样本拓宽决策区域，添加到数据集中的新样本位于原始样本的附近，而不是样本本身，从而降低了过拟合的概率^[16]。

1.2   DDNet算法

Liu Gang等^[17-18]提出了一种只包含矩阵乘法和Hadamard乘积的机器学习算法DDNet。与传统的全连接神经网络(Full-connect Neural Network，FNN)相比，DDNet使用Hadamard乘积代替非线性激活函数。单层DDNet与FNN的输入输出关系分别为

式中：Y₁为单层DDNet的输出；W ^l，l−1为第l−1个模块到第l个模块的权值矩阵；X为输入数据；“$\circ $”为Hadamard乘积符号，表示2个矩阵对应元素相乘；Y₂为FNN的输出；f(·)为非线性激活函数。

DDNet的体系结构可表示为

式中L为DDNet模块的数量。

DDNet使用Hadamard乘积代替非线性激活函数，将展开式等价为特征之间的逻辑表达，实现了高次幂代替非线性映射的功能，各项前的权重矩阵则转换为泰勒展开式的系数。DDNet收敛后在类似于最优组合点处的泰勒展开，使得其非线性学习所受的限制更小。与传统的FNN相比，DDNet具有较快的有效收敛速度，且不易出现过拟合，泛化能力更好。本文使用DDNet构建神经网络，结构如图2所示。

图  2  DDNet网络结构

Figure  2.  The structure of dendrite net

下载: 全尺寸图片 幻灯片

1.3   特征提取

煤岩显微组分光学特征复杂，依靠单一类型的特征难以区分镜质组与其他显微组分^[19-20]。本文提取包括灰度、几何、纹理特征在内的101维特征用于镜质组识别。其中，灰度特征包括灰度均值、灰度最大值、灰度中值、灰度标准差、对比度、能量、灰度概率（不同灰度的像素在总像素中的占比）^[21]，共70维。为提高计算灰度概率的效率，实验对相邻的4个灰度级合并计算。同时选取7个Hu不变矩描述显微组分的几何特征^[22]。采用灰度共生矩阵描述显微组分的纹理特征，包含4个方向上的逆差矩、二阶矩、熵、对比度、差异性、相关性共24个特征^[23-24]。

MMVR与灰度分布之间存在密切关系，为了准确描述灰度分布并预测MMVR，从截取的镜质组方形窗口中提取14维灰度特征进行MMVR回归分析，包括灰度均值、灰度最大值、灰度中值、灰度标准差及对应像素点数量最多的10个灰度值。

1.4   评价指标

1.4.1   镜质组识别的性能评价指标

通过准确率A、查准率P、召回率R和F₁分数4个性能指标评价镜质组识别算法的性能。

式中：T_P为被预测为正类的正样本数；T_N为被预测为负类的负样本数；F_P为被预测为正类的负样本数；F_N为被预测为负类的正样本数。

1.4.2   MMVR估计的性能评价指标

采用均方误差M_se、平均绝对误差M_ae和决定系数S三个指标评价MMVR估计算法的性能。均方误差衡量MMVR预测值与真实值之间的差异程度。平均绝对误差为预测值与真实值的绝对差值的平均值。决定系数表明回归预测与真实值的接近程度，其取值范围为0~1。决定系数值越大，表示回归的性能越好。各指标计算公式为

式中：n为样本个数；$\mathop {{\hat y}}\nolimits_{i}$为预测值；$ {y}_{i} $为真值；$ \stackrel{-}{y} $为总体真实值的平均值。

2.   实验结果和分析

本文使用来源于美国科罗拉多和西弗吉尼亚的13个煤样进行实验。样本包含69张煤岩显微图像，其MMVR为0.7%~1.79%。10~15个相互独立的实验室严格遵守反射率测定标准ASTM D2798-21(Standard Test Method for Microscopical Determination of the Vitrinite Reflectance of Coal)，使用显微光度计测得13个煤样的MMVR，最终的MMVR数值为各实验室测得数据的平均值。13个煤样的煤岩显微图像均采用Leica DFC480数字照相机在白光照射和油浸物镜下拍摄，且均处于相同曝光条件下。可通过以下网址获取相关实验数据：https://energy.usgs.gov/PhotoAtlas/?aid=14。

2.1   图像分割结果

基于手肘法和K−Means算法的煤岩显微图像分割结果如图3所示。图3(a)−图3 (d)中显微图像所含有的显微组分类别数依次是1−4。分割结果显示，采用手肘法自动确定K−Means算法的参数K，具有良好的自适应能力，能够自动区分不同类别数的显微组分。

图  3  K−Means聚类结果

Figure  3.  K−Means clustering results

下载: 全尺寸图片 幻灯片

2.2   镜质组识别结果

使用K−Means算法分割13个煤样的69张显微图像后，去除小于100×100像素的区域，最终获得显微组分区域共891个，其中镜质组区域168个，非镜质组区域723个，类别分布不平衡。对每个区域提取几何、纹理和灰度共101维复合特征，形成891×101像素的矩阵，作为基础数据集。将数据集按照8∶2的比例随机划分为训练集和测试集进行实验。为了解决镜质组与非镜质组样本不均衡的问题，尝试了3种不同的过采样、下采样算法：① SMOTE算法，通过添加合成的少数类样本改善数据分布的不平衡。② 随机下采样(Random Under Sample，RUS)，通过随机选择对多数类样本进行下采样。③ SMOTE算法与RUS相结合的不均衡数据处理方法。在镜质组识别实验中，对比了4种经典分类算法（分类回归树（CART）、K近邻（KNN）、支持向量机（SVM）、RF）与不同数据处理方法结合得到的分类性能，实验结果见表1。同时，为研究采样后的少数类与多数类样本比例对各模型分类性能的影响，以使用SMOTE算法处理不平衡数据为例，分别将镜质组与非镜质组样本比例设置为0.3∶1，0.5∶1，0.7∶1，0.9∶1和1∶1。用SMOTE算法处理后，样本不平衡问题得到了缓解。随着镜质组与非镜质组样本比例的增大，分类模型整体性能逐渐变好，当该比例达到0.7时，模型性能趋于稳定。为简便起见，在后续实验中，将上述比例设置为1∶1。

表  1  过采样、下采样处理前后结果对比

Table  1.  Comparison of experimental results before and after oversampling and down-sampling

数据处理	分类算法	准确率	查准率	召回率	F1分数
处理前	CART	0.95±0.02	0.88±0.06	0.87±0.07	0.87±0.05
	KNN	0.95±0.01	0.90±0.05	0.86±0.05	0.88±0.04
	SVM	0.95±0.01	0.88±0.06	0.88±0.06	0.88±0.04
	RF	0.96±0.01	0.92±0.05	0.85±0.07	0.88±0.04
RUS	CART	0.93±0.02	0.77±0.08	0.92±0.05	0.83±0.05
	KNN	0.95±0.02	0.82±0.06	0.93±0.04	0.87±0.04
	SVM	0.95±0.02	0.81±0.08	0.95±0.04	0.87±0.05
	RF	0.95±0.02	0.82±0.06	0.96±0.04	0.89±0.04
SMOTE结合RUS	CART	0.94±0.02	0.83±0.06	0.87±0.06	0.85±0.04
	KNN	0.96±0.01	0.86±0.06	0.92±0.05	0.89±0.03
	SVM	0.96±0.01	0.88±0.05	0.90±0.04	0.89±0.03
	RF	0.97±0.01	0.91±0.05	0.93±0.04	0.92±0.03
SMOTE	CART	0.95±0.01	0.85±0.05	0.88±0.05	0.87±0.04
	KNN	0.96±0.01	0.87±0.05	0.91±0.05	0.89±0.03
	SVM	0.96±0.01	0.88±0.05	0.90±0.05	0.89±0.03
	RF	0.97±0.01	0.92±0.04	0.93±0.04	0.92±0.03

下载: 导出CSV 

| 显示表格

为了避免偶然因素影响，表1中的各项数据均为4种分类算法重复实验50次所得的平均值。由表1可知，RF算法在实验中表现最佳。此外，使用SMOTE算法对训练集少数类样本过采样后，4种分类算法在准确率和F₁分数略有提升的同时，召回率大幅增加，且各项指标的标准差有所降低。使用RUS算法处理后，镜质组与非镜质组样本的比例被提升为1∶1，数据不平衡问题得到了解决。相较于不做下采样处理，模型的平均召回率从0.85提升到0.96，平均F₁分数从0.88提升为0.89，但查准率明显下降，预测结果中可能包含较多的假阳性样本。实验结果表明，仅使用SMOTE算法即可有效避免模型因过度学习样本先验信息而导致对多数类识别好、少数类识别差的问题。

2.3   基于煤样的MMVR回归结果

受传统MMVR测定方法的启发，采用基于煤样的MMVR回归分析，单个煤样的MMVR估计值为其包含的所有煤岩显微图像的MMVR估计值的均值。在估计单幅煤岩显微图像的MMVR时，选取煤岩显微图像中大于100×100像素的镜质组区域进行MMVR估计。模仿传统取点测量的方法，从镜质组区域中截取多个方形窗口（41×41像素），所有方形窗口MMVR估计结果的平均值为该煤岩显微图像的MMVR。

煤岩显微图像的MMVR估计结果如图4所示，其中图4(b)为多个方形窗口的反射率预测结果分布。该煤岩显微图像的MMVR实际值为1.31%，多个方形窗口预测结果的平均值为1.285%，证明了对单幅显微图像截取窗口进行回归分析的有效性。

图  4  煤岩显微图像的MMVR分析结果

Figure  4.  MMVR analysis results of the coal photomicrograph

下载: 全尺寸图片 幻灯片

在煤样MMVR回归实验中，对比了7种回归分析算法，包括支持向量机回归(Support Vector Regression，SVR)、自适应增强(AdaBoost)、K 近邻(K-Nearest Neighbor，KNN)、梯度提升(Gradient Boosting)、RF、FNN和DDNet。实验采用留一法进行交叉验证，为了防止偶然因素影响，共进行了20次回归实验，最终得到的均方误差、平均绝对误差、决定系数结果对比见表2。

表  2  回归算法测试结果对比

Table  2.  Comparison of test results of regression algorithms

算法	均方误差	平均绝 对误差	决定系数
SVR	0.014	0.080	0.885
AdaBoost	0.010	0.078	0.919
KNN	0.010	0.071	0.925
Gradient Boosting	0.009	0.071	0.926
RF	0.009	0.070	0.926
FNN	0.033	0.064	0.735
DDNet	0.001	0.027	0.990

下载: 导出CSV 

| 显示表格

由表2可知，DDNet回归分析算法取得了最佳性能。DDNet作为FNN的改进算法，在将非线性激活函数替换为Hadamard乘积后，其回归分析的精度和泛化性能均有显著提升。

DDNet回归模型的预测结果如图5所示，直观地展示了13个煤样MMVR的实际值和预测值之间的相关性。图5中MMVR的预测值与实际值高度契合，证明了将DDNet回归算法用于MMVR预测估计具有较强的可行性。

图  5  DDNet模型的预测值与真实值对比

Figure  5.  Comparison of predicted value of DDNet and the ground truth

下载: 全尺寸图片 幻灯片

2.4   软件界面

为协助煤岩分析工作者进行MMVR测定，开发了一款煤岩MMVR估计软件。该软件集成了基于K−Means的图像分割算法、基于RF的镜质组识别算法和基于DDNet的MMVR估计算法。需要说明的是，该软件的适用对象需符合反射率测定标准ASTM D2798-21(Standard Test Method for Microscopical Determination of the Vitrinite Reflectance of Coal)。针对图像分割中超参数K的设定，软件提供手肘法自动确定和人工设定2种方式。

煤岩MMVR估计软件界面如图6所示。通过“导入煤样图片”按钮批量导入指定煤样的多幅煤岩显微图像，A区域显示导入的显微图像；通过“镜质组识别”按钮对导入的显微图像进行图像分割和镜质组识别，并将识别出的镜质组区域显示于B区域中；通过“切换图像”按钮切换该煤样的不同显微图像及其镜质组识别结果；通过“反射率估计”按钮估计镜质组区域的MMVR值，并将该煤样所截取窗口的反射率分布柱状图和MMVR估计结果显示在C区域和D区域中。

图  6  煤岩MMVR估计软件界面

Figure  6.  MMVR estimation software interface for coal and rock

下载: 全尺寸图片 幻灯片

3.   结论

1） 基于机器视觉技术开发了煤岩MMVR估计系统，主要包括煤岩显微图像分割、镜质组识别和MMVR回归3个部分。采用K−Means算法，分割显微图像中的不同显微组分，提取镜质组灰度、几何和纹理等特征；为改善样本不均衡问题，采用SMOTE算法对少数类样本过采样，构建了4种镜质组区域识别模型，其中RF方法性能最优，分类准确率为97.0%；建立了7种回归估计模型，其中DDNet回归算法取得了最优的结果，决定系数达到了0.990。

2） 实验结果表明，本文所提方法与传统基于显微光度计测定方法的测定结果高度契合，且克服了传统显微光度计测定方法对时间、精力、专业知识要求高的缺点，验证了机器学习在煤岩显微图像分析中应用的可行性。

3） 下一步拟对更多样本、更宽镜质组反射率范围的煤样进行分析，并尝试使用更复杂的神经网络完成图像分割、分类及回归任务，以减少对特征工程的依赖。