改进多目标萤火虫优化的软子空间聚类算法及在短期负荷预测中的应用

张 曦 康 平 付雪峰,2 叶 军,2 赵 嘉,2,3*

1(南昌工程学院信息工程学院 江西 南昌 330000) 2(南昌工程学院江西省水信息协同感知与智能处理重点实验室 江西 南昌 330000) 3(南昌工程学院鄱阳湖流域水工程安全与资源高效利用国家地方联合工程实验室 江西 南昌 330000)

0 引 言

在大数据时代，高维数据聚类是一项重要且具挑战性的工作。高维数据中大量的无关属性使目标簇所在的子空间只与某些维度有关；稀疏的数据分布使数据之间的距离几乎相等，导致传统聚类方法中的全维空间距离度量失效。为解决以上问题，软子空间聚类(Soft Subspace Clustering,SSC)为每个维度分配一个0到1之间的权值，以权重来表示维度对目标簇的贡献程度[1]。

许多传统的软子空间聚类将若干标准合并成一个目标函数进行优化，由于这些标准之间可能会相互冲突，其在目标函数中引入加权参数来平衡这些标准以提高性能，但是如何设置适当的加权参数仍然是一个未解决的难题，而且聚类效果对加权参数的设置极为依赖，同一参数无法适用于不同的数据集。因此，相关学者将软子空间聚类问题转化为多目标优化问题(Multi-objective Optimization Problem,MOP)以减少这种加权参数，并运用多目标进化算法对多个目标函数同时进行优化，以提高软子空间聚类的性能[2]。Zhu等[3]提出了一种基于多目标进化算法的软子空间聚类(Soft subspace clustering with a multi-objective evolutionary approach,MOSSC)。Xia等[4]考虑到权重的分布，将负权值熵和簇间分离度整合成一项，并将其与加权簇内紧凑度视为两个目标同时进行优化，提出了一种基于多目标进化方法的软子空间聚类(A novel soft subspace clustering with a multi-objective evolutionary approach,MOEASSC)。在MOEASSC的基础上，Liu等[5]采用基于分解的多目标进化算法(A multi-objective evolutionary algorithm based on decomposition,MOEA/D)[6]对簇内紧凑度、簇间分离度、负权值熵三个目标同时进行优化，提出了一种对高维数据聚类的改进多目标进化方法(An improved multi-objective evolutionary approach for clustering high dimensional data,MOEC)。虽然MOEC比一些软子空间聚类算法的性能要好，但是其采用的MOEA/D算法寻优精度不高，获得的最优解集可能分布不均匀[7]，影响了聚类效果。

2013年剑桥学者Yang[8]提出了一种多目标萤火虫算法(Multi-Objective Firefly Algorithm,MOFA)，相对于几种常用的多目标进化算法而言能得到较均匀的最优解集。由于MOFA易提前收敛、寻优精度较低，许多MOFA的改进算法相继被提出。谢承旺等[7]采用正交实验方法产生初始群体，利用外部档案中的精英解指导萤火虫移动，并运用3点最短路径方法保持档案群体多样性，提出了一种混合型多目标萤火虫算法。Lv等[9]在萤火虫学习公式中加入补偿因子进行迭代，并随机选取一个外部档案粒子作为精英粒子参与种群的进化，提出了一种基于补偿因子与精英学习的多目标萤火虫算法(Multi-objective firefly algorithm based on compensation factor and elite learning,CFMOFA)。现有的MOFA能较好地解决一些MOP问题，但其性能仍有很大的提升空间[9]。

基于上述思想，本文采用多目标萤火虫算法优化簇内紧凑度、簇间分离度及负权值熵三个目标函数，同时，为弥补MOFA易提前收敛、寻优精度较低的缺陷，将萤火虫位置更新公式中的步长α、吸引力β0进行动态定义，并设计一种萤火虫单行随机学习的机制来提高最优解集分布的均匀性，提出一种改进多目标萤火虫优化的软子空间聚类算法(IMOFASSC)。IMOFASSC将簇心和权重拟化成萤火虫种群，通过萤火虫位置的更新来搜索最优簇心和权重，并发掘子空间中隐藏的簇。在UCI标准数据集和癌症基因表达数据集上的实验结果表明，IMOFASSC算法在处理低维和高维数据时较其他算法具有更好的聚类效果。

短期负荷预测对于电力系统的可靠性至关重要，因为负荷预测结果是为确定系统是否易受攻击而进行离线网络分析的基础。如果预测不充分，则会导致分配不足的备用容量和使用昂贵的调峰装置，或者导致不必要的大备用容量，这两者都与增加系统运行成本有关。此外，在公开市场环境中，短期负荷的准确预测是电能交易和现货价格确定的基础，旨在获得最低的购电成本。因此，开展短期电力负荷预测研究工作，提高其预测精度，具有重大的现实意义[10]。基于该背景，本文将IMOFASSC算法应用到短期电力负荷预测中，仿真实验结果表明IMOFASSC算法在短期电力负荷预测中具有一定的推广应用价值。

1 相关工作

1.1 软子空间聚类

软子空间聚类算法的最小化目标函数可以概括为：

J(V,U,W)=Jwc(V,U,W)+γ·Jadd(V,U,W)

(1)

Jwc通常表示为：

(2)

式中：Jwc表示加权簇内紧凑度，其值越小簇内样本点越紧凑；Jadd是一个从其他角度评估聚类结果的附加项；γ是用户设定的参数，用于平衡Jwc和Jadd；V=[vik]C×D为簇心矩阵；U=[uij]C×N是划分矩阵，表示N个样本分别属于C个簇的程度；W=[wik]C×D为权重矩阵，表示C个簇分别在D个维度上的权重；X=[xjk]N×D是D维空间中N个样本的数据矩阵；m(m≥1)和β(β≥1)是两个参数。

最早提出的软子空间聚类算法不包含Jadd附加项[11-12]，为进一步提高算法性能，许多学者在目标函数中添加一些附加项。为避免权重无限大，Jing等[13]提出了一种模糊加权K均值(Fuzzy weighting k-Means,FWKM)算法，其中：Jadd是全部样本所在簇的模糊权重的总和；γ是整个样本的平均离散度。Gan等[14]提出了模糊子空间聚类(Fuzzy Subspace Clustering,FSC)算法，其中：Jadd是模糊权重总和；γ是一个恒定参数。为使算法选择更多的维度来搜索子空间，熵加权K均值(Entropy weighting k-Means,EWKM)算法[15]和局部自适应聚类(Local Adaptive Clustering,LAC)算法[16]都采用负权值熵作为附加项。Deng等[17]在最小化EWKM目标函数的同时最大化加权簇间分离度，提出了增强的软子空间聚类(Enhanced Soft Subspace Clustering,ESSC)算法，其目标函数如下：

J(V,U,W)=Jwc(V,U,W)+γ·Jadd(V,U,W)-

η·Jbs(V,U,W)

(3)

式中：Jbs表示加权簇间分离度；η是用户预设的控制Jbs影响的参数。

在式(1)和式(3)对应的软子空间聚类算法中，γ和η值难以设置，因为式中项与项之间可能会冲突，没有确切的方法定义适当的γ和η值来平衡这些项[4,17]。例如，对于Jadd项，EWKM、LAC和ESSC中的负权值熵应被最小化以在更多维度上搜索子空间，可能导致某些簇的子空间被一些不太相关的特征错误地构造，虽然最后也能得到较小的Jwc，但是结果不准确；对于Jbs项，在ESSC中它也可能与Jwc冲突，例如，具有较小Jwc的某些簇可能在样本点密集的区域中挤在一起，使它们之间的分离度较小。

1.2 多目标萤火虫算法

自然界中，萤火虫利用自身的发光特性来相互吸引、传递信息，实现求偶交配和繁殖的目的。MOFA包含亮度和吸引力两个寻优要素，亮度决定萤火虫的移动方向，吸引力体现移动距离，迭代中萤火虫不断被较亮萤火虫吸引并向其移动。在M维目标优化问题中，萤火虫i适应值函数表示为：

f(xi)=(f1(xi),f2(xi),…,fM(xi))

(4)

其亮度由f(xi)决定，xi表示萤火虫i的位置。萤火虫的吸引力是相对的，它随萤火虫i与萤火虫j的距离的变化而变化，同时吸引力也与吸收因子有关，由于传播媒介对光的吸收，光的亮度随着与光源的距离增大而变小。吸引力β可定义为：

(5)

(6)

式中：β0为初始吸引力，即在rij=0处的吸引力，一般取为1；γ为光吸收系数，标志吸引力的变化，一般取γ=[0.01,100]；rij为萤火虫i到萤火虫j的欧氏距离；d是数据维度。

MOFA以Pareto支配[18]的概念确定萤火虫之间的吸引关系，即：如果萤火虫j支配萤火虫i(记为j≻i)，萤火虫i就被更亮的萤火虫j吸引，位置更新公式如下：

xi(t+1)=xi(t)+βji(xj(t)-xi(t))+α·εi

(7)

式中：t为当前迭代次数；xi(t)、xj(t)为萤火虫i和j的第t代空间位置；βji表示萤火虫j对萤火虫i的吸引力；α为步长因子，一般取[0,1]上的常数；εi是服从均匀分布或高斯分布或者其他分布的随机数向量，一般用rand-0.5来表示，rand为[0,1]上服从均匀分布的随机数。

若萤火虫不被其他任何萤火虫支配则采用式(8)更新它的位置。

xi(t+1)=g*+α·εi

(8)

标准MOFA通常通过式(9)将多个目标函数以随机加权和的方式转换为单目标函数，g*是ψ(x)值最优时萤火虫的位置。

(9)

式中：M为目标函数个数；ωm∈[0,1]，并且每一代ωm都不相同。

每一代萤火虫完成进化后，将所有最优解保存在外部档案(External Archive,EA)中，并对外部档案进行维护，以确保外部档案中所有的解互不支配。为节省计算资源，一般会设定外部档案的最大容量。当最优解的个数超过外部档案数时采用自适应网格删除法[19]确保算法所得的Pareto最优解分布更均匀。

2 改进多目标萤火虫算法优化的软子空间聚类(IMOFASSC)

2.1 目标函数

本文采用文献[5]提出的一种对高维数据聚类的改进多目标进化方法MOEC的目标函数作为IMOFASSSC的目标函数，见式(10)。

(10)

(11)

式中：xjk为第j个数据点第k维的值；Awi表示第i个簇重要权重的平均值，即该簇中值大于1/D的权重的平均值，Awi的最小化可以抑制一些非常大的权重；Sepi表示第i个簇心与其他簇心之间的欧几里得距离的总和；σ一般设为0.000 1，可避免分母为零。

2.2 多目标萤火虫算法改进

根据我们之前的研究[20]，萤火虫算法的性能很大程度上取决于其参数，如式(7)中固定步长α和初始吸引力β0。当萤火虫算法收敛时，α应满足当t→∞时limα(t)=0。根据式(5)，β由两只萤火虫之间的距离决定，随着迭代的进行，算法不断收敛，萤火虫之间的距离逐渐减小到零，在后期的迭代中吸引力β恒等于β0，对局部搜索没有帮助。基于以上思路，本文对MOFA中的α和β0做出同样的修改。

α(t)=(1-t/Maxiter)α0

(12)

(13)

式中：rand1为第rand1代萤火虫的步长；rand1为初始步长；rand1为最大迭代次数；rand1是第rand1代萤火虫的初始吸引力；rand1是第rand1代萤火虫rand1和rand1之间的吸引力；rand1和rand2表示范围在j服从均匀分布的两个随机因子。

MOFA外部档案中的最优解又称精英解，精英解包含种群的优质信息，而在种群进化时精英解没有得到利用，浪费了计算代价。此外，标准MOFA没有考虑到萤火虫彼此不存在支配关系时该如何移动的问题，以上原因导致了萤火虫种群进化缓慢，多样性差。针对这些问题，本文设计一种萤火虫单行随机学习的机制。

(1) 当萤火虫j支配萤火虫i时，即j≻i时，则i利用式(14)向j移动。

xi(t+1)=xi(t)+β(t)·(xj(t)-xi(t))+α·εi

(14)

(2) 当萤火虫i不被萤火虫j支配时，则以萤火虫i的每一行为研究对象，在外部档案EA中随机选取一个精英解g*，利用式(15)对i的当前行进行位置更新并用一个新粒子i′来记录更新结果。所有行全部更新之后得到xi′，判断此时的i′和i之间的支配关系，若i′≻i，则xi(t+1)=xi′，否则不作处理。

(15)

2.3 编码方式及初始化

萤火虫个体采用混合实数编码，如图1所示，其中i=(1,2,…,C)，C和D的定义同式(2)。例如：在萤火虫第i行的位置中，前半部分代表第i个簇的簇心，后半部分代表第i个簇的权重。在计算目标函数之前，应将各个权重除以它们所在簇的权重总和，以确保每个簇的权重和等于1。

图1 萤火虫个体编码

C

/D

V

W

U

(16)

式中：uij表示第j个样本属于第i个簇的程度，其确定方法是计算第j个样本到各个簇心的加权距离，若j到某个簇心的距离最小，则u=1，否则u=0；1≤p≤C。

2.4 IMOFASSC算法步骤

IMOFASSC的流程如算法1所示。

算法1IMOFASSC算法

1.初始化萤火虫种群并设置参数，如：萤火虫数量Npop，外部档案最大容量Epop，最大迭代次数Maxiter；

2.利用式(10)计算每个萤火虫的三个目标函数值；

3.while(t

4.fori=1:Npop

5.forj=1:Npopd

6.ifxj≻xi

7.萤火虫j利用式(14)向j移动；

8.else

9.forh=1:C

11.end for

12.ifxi′≻xi

13.xi=xi′；

14.end if

15.end if

16.利用式(10)更新萤火虫i的三个目标函数值；

17.end for

18.end for

19.更新并保存外部档案中的Pareto最优解；

20.t=t+1；

21.end while

22.输出外部档案中最优的V和W

3 实验与结果分析

为检验IMOFASSC算法的优良性能，本文使用6个UCI标准数据集和4个癌症基因表达数据集[21-22]进行测试，并将IMOFASSC与FSC[14]、EWKM[15]、ESSC[16]、MOEASSC[4]和MOEC[5]等5种先进的软子空间聚类算法进行比较。数据集基本信息在表1中给出。

表1 数据集信息

本文采用Rand Index(RI)[23]和Normalized Mutual Information(NMI)[24]两个指标评价聚类效果。

(17)

(18)

式中：K为标签数量；C为聚类簇数；N为样本总数；f00表示具有不同标签的点被分配到不同簇的样本对数；f11表示具有相同标签的点被分配到同簇的样本对数；ni表示具有第i个标签的样本数；nj表示簇j中的样本数；nij表示具有第i个标签的样本点被分配到簇j的样本对数。RI和NMI的值在[0,1]区间内，两指标越接近1，聚类效果越好。

实验所有数据集归一化处理，参数设置见表2。为公平比较，设置所有算法的最大迭代次数为5 000，每种算法运行30次，记录最佳结果见表3、表4。由表3、表4中的RI、NMI评价指标可知，IMOFASSC算法在6个UCI数据集上的RI、NMI最好值均明显优于其他5种算法，表明IMOFASSC在处理低维数据集时具有较高的精度。其次，在CNS_tumors、DLBCL_B、StJude_Leukemia三个癌症基因表达数据集上，IMOFASSC取得的RI、NMI最好值均最优，这说明在高维数据维度灾难的影响下，IMOFASSC不仅性能更好而且具有一定的稳定性。最后，从所有数据集上来看，IMOFASSC一共在9个数据集上取得了最优的RI、NMI值。此外，还可发现MOEASSC和MOEC的整体性能也比单目标的FSC、EWKM和ESSC要好，这印证了第1节中的观点，即式(1)和式(3)中的项彼此冲突，γ和η参数的设置对聚类效果影响极大，多目标进化算法实现的软子空间聚类可以同时优化这些项有效地解决这个问题。

表2 算法参数设置

表3 每种算法运行30次的最佳结果(RI)

表4 每种算法运行30次的最佳结果(NMI)

综上所述，对于低维UCI数据集，IMOFASSC算法能够有效地提高聚类精度，同时在高维癌症基因表达数据集上也取得了较好的聚类效果。

4 短期负荷预测模拟实验

4.1 IMOFASSC-SVR模型

短期负荷预测是一个受时间、经济、季节、天气和随机效应多因素影响的复杂的非线性问题。在众多预测方法中，支持向量回归(Support vector machine for regression,SVR)在解决有限样本、凸二次规划及非线性高维的问题时具有优异的性能，已成为短期负荷预测领域的研究热点问题之一[25]。但是，SVR预测性能对训练样本数据的相关性极为敏感。传统SVR预测方法未对训练样本进行筛选，训练样本数据量大，相关性较弱，导致预测性能不佳。针对上述问题，本文采用IMOFASSC实现SVR训练样本的自动选取，建立IMOFASSC-SVR模型。

本模型利用IMOFASSC对历史日和预测日的影响因子数据进行聚类，并找到预测日所在的类，则该类中的历史日即为预测日的相似日，这些历史日的负荷数据与相应的影响因子数据共同构成SVR预测模型的训练样本。为检验模型预测性能，采用平均绝对百分误差(Mean Absolute Percentage Error,MAPE)作为评价指标，计算公式如下：

(19)

具体预测步骤如下：

(1) 寻找历史负荷数据，确定影响因子，并进行归一化处理。

(2) 利用改进的软子空间聚类选取与预测日相似度较高的历史负荷数据，与影响因子一起作为训练样本。

(3) 将训练样本的影响因子作为SVR的输入向量，对应的负荷数据作为输出向量，训练SVR预测模型。

(4) 在训练好的SVR预测模型中输入预测日的影响因子，输出值反归一化后即为相应的预测结果。

(5) 计算预测误差MAPE。

4.2 仿真实验

本节实验数据引自第九届电工数学建模竞赛试题A，即以某地区2014年10月25日至12月24日每天96个时刻点的负荷数据和气象数据作为的训练样本，2014年12月25日至12月31日每天96个时刻点的实际负荷值作为测试样本进行实验，并将IMOFASSC-SVR模型与SVR模型进行比较。

为公平起见，IMOFASSC算法中的聚类簇数为3，种群规模为20，最大迭代次数为100，其他参数与第3节相同。IMOFASSC-SVR与SVR模型的惩罚系数C、不敏感损失系数ε和核参数σ的取值分别设为默认值，即[0.01,1 000]、[0,1]和[0.1,100]，SVR核函数设为径向基核函数。实验中IMOFASSC-SVR和SVR模型独立运行20次，记录两种模型得到的最佳预测结果，绘制出两种模型的预测结果，见图2。

图2 两种模型的预测结果

图2表示这两种模型7天672个采样点的负荷预测值与实际负荷值的比较。可以看出，两种模型的预测负荷曲线与实际负荷曲线的走势基本相同，但在前3天和第4至第5天期间，两种模型的预测负荷曲线与实际负荷曲线基本重合，说明两种模型的预测效果相当；在第4天预测阶段，两种模型的预测负荷曲线均与实际负荷曲线存在较大的间距，表明两种模型的预测效果不佳；在第7天的预测中，SVR模型的预测曲线与实际负荷曲线存在明显的间距，IMOFASSC-SVR模型的预测曲线与实际负荷曲线间距较小，说明IMOFASSC-SVR模型的预测效果较好。

为进一步比较模型性能，将IMOFASSC-SVR模型的预测误差MAPE与标准SVR预测模型进行对比，结果见表5。

表5 四种模型的预测误差(MAPE)(%)

续表5

表5显示了两种模型分别对12月25日至12月31日的电力负荷进行预测时产生的误差MAPE值，以及两种模型在7天预测中的MAPE误差的平均值。可以看出，在12月25日、26日、29日、31日的预测中，由于IMOFASSC-SVR采用IMOFASSC自动选取预测日的相似日训练样本使训练出的SVR模型更准确，有利于后续预测，因此IMOFASSC-SVR模型的MAPE值均小于SVR模型，尤其31日的误差较SVR模型明显降低了4.37百分点，表明IMOFASSC-SVR模型在这几天的预测效果优于SVR模型；对于27日、28日、30日，SVR模型取得了最优的MAPE值，这是因为当预测日的气象情况差别很大时，在气象数据的聚类过程中仅设置一个固定的聚类簇数变得不准确；从两种模型7天的平均MAPE值来看，IMOFASSC-SVR模型7天的平均MAPE值为2.49%，SVR模型为3.34%，表明IMOFASSC-SVR模型在这7天的整体预测效果最优；从全局来看，IMOFASSC-SVR模型在7天的预测中有4天取得了最优的预测效果，并且7天的整体预测效果也是最优的，而SVR模型只在3天中的预测效果较好，以上分析表明IMOFASSC-SVR模型的预测性能优于SVR模型。

5 结 语

本文提出一种改进多目标萤火虫优化的软子空间聚类算法。IMOFASSC首先对MOFA的步长α、吸引力β0进行了动态定义；设计一种萤火虫单行随机学习机制来提高最优解集分布的均匀性；将改进的MOFA运用到软子空间聚类问题中，同时优化簇内紧凑度、簇间分离度及负权值熵三个目标函数。IMOFASSC将簇心矩阵和权重矩阵拟化成萤火虫种群，通过萤火虫位置的更新来搜索最优簇心和权重，并发掘子空间中隐藏的簇。在UCI标准数据集、癌症基因表达数据集和短期电力负荷预测的实验结果表明，IMOFASSC不仅在低维数据聚类中有较高的精度，而且在处理高维数据时也表现出良好的性能，可以有效地解决软子空间聚类中的MOP问题。在短期电力负荷预测模拟实验中，IMOFASSC有效地去掉训练样本中一些无关数据，加强数据间的相关性，使IMOFASSC-SVR模型预测精度较SVR更高，表明了IMOFASSC在短期负荷预测中具有一定推广应用价值。下一步的研究方向是继续研究改进IMOFASSC，使其能同时识别子空间和簇的数量。