1. 引言及文献综述
2011年10月29日，国家发展改革委批准北京、天津、上海等7省市开展碳排放权交易试点工作，开启了我国碳排放权交易市场(以下简称碳市场)建设的重要一步。2021年7月，全国性的碳市场在上海正式启动交易，持续多年的碳市场试点工作迎来统一，全国碳市场的开市标志着政府主管部门开始依靠市场机制推动碳减排。
然而，与股票市场不同，碳市场交易参与者相对较少，市场交易稀薄，碳交易价格过分偏离碳减排成本，碳价格波动时有异常，市场价格发现功能弱化。为了完善碳排放交易机制，更好利用市场的价格发现功能，完成碳排放控制目标，同时，降低碳金融市场的交易风险，有必要对影响碳价格波动的因素以及碳价格预测问题展开更深入探索。
(一) 碳市场交易价格影响因素研究
基础性的研究以分析碳市场价格波动的影响因素为主，通过不同类型的分解模型或分解方法，抽取影响碳市场价格波动的各种因素，这些因素包括政府双碳目标政策、最终需求总量变化、产业结构调整、等宏观层面因素 [1] ，及人均财富、技术水平、碳信息披露质量、产业集聚、能源价格、环境政策和对外贸易度等微观层面因素 [2] ，通过分析具体因素与研究对象之间关系，给出经济发展的政策建议。
进一步的研究则在因素分析的基础上，获取因素变量的分布特征、高阶矩以及碳期货价格的隐含波动率等信息 [3] ，研究因素变量的波动规律以及时变特征 [4] 。或者，从国际碳价以及国内外经济指标等多维度构建指标，筛选并识别碳价格的影响因素 [5] 。
(二) 碳市场价格预测研究
更为深入的研究则在时域分析基础上，充分考虑到碳价格的非线性和非平稳性以及尖峰厚尾和偏斜等特点，使用人工智能以及频域分析等理论和工具，综合构建具有自适应特征的频域分析模型或者基于数据自身结构特点的人工智能类分析方法，研究碳排放权交易价格的预测问题，可概括为如下三类基本模式。
第一类为“时域模型 + 人工智能”模式，时域模型多以GARCH或ARIMA为基础，人工智能模型多以机器学习和神经网络为主，如朱帮助等构建的基于数据分组处理方法(GMDH)–粒子群算法(PSO)–最小二乘支持向量机(LSSVM)的国际碳市场价格预测模型 [6] ，George S等构建的混合神经模糊控制模型 [7] ，Lei Ji等构建的自回归移动平均–卷积神经网络–长短期记忆网络(ARIMA-CNN-LSTM)组合模型 [8] 。
第二类为“频域模型 + 时域模型或人工智能”模式，频域模型多以经验模态分解(以下EMD)为基础，如朱帮助等构建的基于EMD的进化最小二乘支持向量回归多尺度集合预测模型 [9] 。Jiaming Zhu等构建的针对非平稳、非线性碳价格序列的变分模态分解(VMD)、模态重构(MR)和最优组合预测模型(CFM)相结合的混合预测模型 [10] 以及Wei Sun等构建的EMD、VMD和LSTM的混合碳价格预测模型 [11] 。
第三类为“频域模型 + 时域模型 + 人工智能”模式，是上述两类模式的综合，如周坤等构建的EMD、极端梯度提升(XGBoost)、极限学习机(ELM)、分数阶季节性灰色系统(FSGM)组合模型 [12] 以及Boqiang Lin等考虑结构断点而构建的结合Bai & Perron检验、迭代累积平方和算法(ICSS)、小波变换和长短时记忆神经网络(LSTM)的混合预测模型 [13] 。
本文的创新之处主要表现在：(1) 构建了碳交易价格序列的频域解析表达式，函数解析式是变量之间关系的最为理想地表达，但在时域分析框架下，碳价格影响因素不可能全部明确列出，导致无法构建碳价格影响因素与碳价格之间的函数解析式，本文则在频域分析框架下，通过现代信号分析与数据处理技术，完备且清晰地提取不同频域的影响因子，进而构建碳价格波动的频域解析式，预测碳价格，预测思想和理念上有新意。(2) 将原始碳价格序列经过完全集合经验模态分解后得到的各IMF分量重构为高频、低频和趋势项，探究不同时间尺度下碳市场交易价格的影响因素差异。(3) 当前对碳市场价格预测工作的研究主要集中于区域性碳市场价格预测和欧盟等西方经济体碳价格预测，以全国性碳市场为研究对象的碳价格预测鲜有着墨。考虑到区域性市场和欧盟市场与全国性市场在市场规模、市场交易者、运行机制等方面的差异，因此，有必要基于全国碳市场交易数据预测碳价格的波动规律。
2. 研究方法介绍
(一) 完全集合经验模态分解
经验模态分解(EMD)通过将原始信号分解成若干个本征模态分量(IMF)和一个代表长期趋势的残差项，得到的不同分量可以更好的呈现出数据在不同时间尺度下的特征。近年来，EMD在处理非平稳非线性信号上得到了较为普遍的应用。
完全集合经验模态分解(以下CEEMDAN)在EMD和EEMD两种方法上加以改进，通过自适应加入多组独立同分布的高斯白噪声，并通过多次迭代和平均以抵消白噪声的影响。相较于EEMD方法，CEEMDAN方法不仅大幅提高了分解效率，同时克服了其在加入白噪声后迭代次数不足导致较差的完备性和出现重构误差的问题。
(二) 傅里叶变换
根据傅里叶理论，任意一个周期函数都可以通过谐波函数的线性组合来表示，非周期函数可以视作周期为无穷大的周期函数，只要找到这些基本谐波函数，将其加总，在满足一定误差精度条件下，可以得到原始价格数据序列的频域解析表达式，时域分析框架下诸多无法列出的影响因素就可以通过频率谱的形式表现，因此，从频域角度构建碳交易价格数据序列的频域解析式，能够成为解决原始序列解析表达问题的可行之道。
3. 实证分析
(一)数据选择与描述性统计分析
本文选取上海环境能源交易所于2021年7月份上线正式交易的全国碳市场碳排放配额(CEA)收盘价数据作为样本，样本期间为2021年7月16日至2022年10月27日共310个交易日的昨日收盘价数据。数据来源于上海环境能源交易所，数据处理主要通过Matlab2019软件进行处理。本文将样本数据分为估计样本和预测样本，估计样本区间为2021/7/16~2022/10/13，预测样本区间为2022/10/14~2022/10/27。
Figure 1. Original series of CEA prices in the national carbon market
图1. 全国碳市场碳排放配额(CEA)价格原始序列图
根据图1可以发现，全国碳交易价格在2022年2月份之前出现暴涨暴跌的波动特征，2月份之后序列波动开始相对平稳，异质性特征显著。同时，由图1发现，在2022年2月份之后多处连续的局部时间跨度内，全国碳市场交易价格没有变化，呈现直线模式，如图1中箭头指向标示。相对于股票市场中股价的波动模式，碳市场中碳交易价格的直线模式意味着市场机制的价格发现功能失灵，“失灵”的原因以及“失灵”如何影响碳价格的预报工作有待研究。进一步给出碳交易价格序列的描述性统计结果如下表1。
Table 1. Descriptive statistical results of carbon trading price series
表1. 碳交易价格序列的描述性统计结果
根据表1可知，最大值和最小值相差较大，同时出现较大的标准差，说明碳交易价格的波动存在较大程度的不稳定性。偏度系数值为−0.915，小于0，数据呈现出明显的左偏分布。峰度系数值为2.094，小于3，相较于正态分布而言，序列呈现扁平厚尾特征，样本序列具备典型金融数据的一般特点。
(二) 碳价格波动特征和影响因素构成分析
1、CEEMDAN分解
此处运用CEEMDAN方法对全国碳市场交易价格序列进行分解，最大程度揭示原始数据的真实波动特征。经过反复比较，最终报告了分解层级为8的分解结果，如图2所示，并且计算了IMF分量的方差贡献率以及IMF分量与原序列的Pearson相关系数，分析价格波动特征及其微观影响因素构成。
Figure 2. CEEMDAN decomposition results of national carbon market trading prices
图2. 全国碳市场交易价格CEEMDAN分解结果
根据图2，IMF分量基本呈现出围绕零均值线、局部极大值和极小值对称的振荡模式。除RES以外的分量的均值基本保持为零，不随时间变化。各IMF分量的振幅和频率的变化较大，波形比原始序列更加规整、简洁，单位根检验表明所有IMF分量序列基本平稳。
与IMF7、IMF8以及RES序列不同，IMF1、IMF2、IMF3、IMF4以及IMF5和IMF6并非在整个样本区间始终呈现出规律性周期波动特征，如图2所示，IMF1和IMF2分量的波动行为在部分期间存在显著的近似直线模式。把IMF1、IMF2、IMF3、IMF4、IMF5和IMF6序列中含有直线的波动模式定义为间歇式波动模式，把IMF7、IMF8以及RES序列表现出的有规律的周期性波动称为周期式波动模式。
产生间歇式波动模式的原因可归结为碳价格数据在连续的多个交易日保持不变，反映出投资者在碳市场的交易行为极度不活跃，导致分解后的序列呈直线模式。从序列微观因子之间相互作用的角度分析，直线模式可归结为序列微观因子之间相互对冲的结果，内部微观因子之间的对冲效应导致综合效应为零，致使特定样本区间内波动状态呈现出直线模式。
定义相邻两个极大值点之间的时间间隔为序列波动的特征时间尺度，根据图2可知，IMF1、IMF2、IMF3和IMF4序列波动的时间尺度特征值均未超过10个交易日，具有随机扰动的特点，反映了短期的投机性交易行为。IMF5以及IMF6序列波动的时间尺度特征值均超过30个交易日，IMF7和IMF8序列的时间尺度特征值均超过60个交易日，IMF8序列的时间尺度特征值超过180个交易日，反映了重大事件冲击以及长期经济政策对碳交易价格的影响效果。
Table 2. IMF component variance contribution rate and Pearson correlation coefficient
表2. IMF分量方差贡献率和Pearson相关系数
根据表2可知，IMF6、IMF7以及IMF8中低频分量的方差贡献之和为69.87%，是序列波动的主要影响因素，并且IMF6、IMF7以及IMF8与原始序列的Pearson相关系数分别为0.373、0.549、0.867，此外，RES的方差贡献率为26.58%，与原始序列的Pearson相关系数为0.720，均呈现出显著的相关性。
2、重构IMF
为进一步分析序列波动的主要影响因子，对经过CEEMDAN分解得到的8个IMF序列进行高、低频区分，分析方差贡献率及Pearson相关系数，结果见表3。
Table 3. Variance contribution rate and Pearson correlation coefficient of reconstructed high, low frequency and trend terms
表3. 重构后的高、低频和趋势项的方差贡献率、Pearson相关系数
由表2和表3分析可知，由IMF1、IMF2、IMF3、IMF4、IMF5以及IMF6构成的高频分量的时间尺度特征值不超过30个交易日，其中，IMF1、IMF2、IMF3和IMF4序列波动的时间尺度特征值未超过10个交易日，同时，高频分量对原始序列的方差贡献率只有11.96%，Pearson相关系数未超过0.4，这说明由IMF1、IMF2、IMF3、IMF4、IMF5及IMF6构成的高频分量代表的影响碳市场价格的因素发生的周期较短，对全国碳市场原始价格序列的波动和影响程度较小。
究其原因，一方面市场交易主体较少，考虑到当前全国性碳市场尚处于早期建设完善的发展现状，首批纳入全国碳市场的交易单位只有2162家发电行业的重点碳排放企业，交易主体稀薄，并且可供交易的产品只有碳排放配额，导致供需双方对交易的产品和数量存在矛盾。另一方面市场机制有待完善，与证券市场相比，碳市场在信息披露、制度机制建设和金融化程度等方面都不具备投机条件，大大削弱了碳市场的投机属性。与国外发展程度相对成熟的碳市场以及国内九个区域性碳市场相比，全国碳市场的交投活跃度有待提高，因此，高频分量所代表的投机、短期供需不平衡等短周期因素对全国碳市场价格的影响较小。
低频分量和趋势项对原始序列的方差贡献率分别达到了65.21%和22.84%，Peason相关系数达到了0.9259和0.7200，结合低频分量中的IMF7和IMF8的时间尺度均超过了60个交易日，IMF8的时间尺度达到了180个交易日的周期特征，说明低频分量和趋势项与碳市场原始价格序列存在着高度相关性，其代表的影响碳市场价格的因素发生周期更长、对全国碳市场价格的波动和影响程度较大。结合当下全国碳市场的发展情况分析，这可能是由于低频分量代表的全国碳市场现行的内在交易机制以及趋势项代表的长期价值的影响导致的，这些因素的影响程度更深、持续时间更久。
(三) 构建基于CEEMDAN和频域解析分析组合方法的预测模型
快速傅里叶变换(FFT)是基于离散傅里叶变换(DFT)的一种快速算法，在信号分析领域得到了广泛的应用。傅立叶变换实现了信号从时域到频域的转换，在频域视角下，可以发现信号在时域视角下隐藏的某些特征。依据这种思想，把碳价格影响因素与碳价格波动规律之间的内在关系，通过傅里叶变换转换到频域视角下清晰的表现出来。
根据傅里叶理论，任意一个周期函数都可以通过谐波函数的线性组合来表示，非周期函数可以视作周期为无穷大的周期函数，在满足一定误差精度条件下，可以将非周期的观测数据分解为正弦或余弦谐波，通过确定主要谐波分量的振幅、频率和相位，得到三角函数的线性组合作为基本谐波函数，并将其加总，可以得到原始价格数据序列的频域解析表达式并通过频域解析式分析预测碳价格。
这里对全国碳市场交易价格原始序列经过CEEMDAN分解重构后得到的高频分量、低频分量以及趋势项分别做快速傅里叶变换，估计各分量的振幅、角频率和相位参数获取分量的频域解析式，并将其加总，构建全国碳市场交易价格的综合频域解析式并进行预测。
Figure 3. Actual amplitude and phase frequency curves of the high-frequency component after fast Fourier transform
图3. 高频分量快速傅立叶变换后的实际振幅和相位频率曲线图
对重构后得到的高频分量图3(a)，由于FFT变换的对称性，图3(b)仅绘出了前半部分的图像，根据图3(b)可以发现高频分量经过快速傅立叶变换以后在第10个频率上出现峰值。根据每个点的频率计算公式(1)：
(1)
其中，Fn表示采样频率，此处为300，N表示采样点，此处为300，计算得到频率分量10对应点11。MATLAB工作区数据分析表明，点11的值比附近的值都要大，这表明11频率点上对应的信号振幅较大，计算出11点的模值为124.845，将模值乘以2/N得到该频率下的真实振幅为0.8323。第一个点表示直流分量，值为0.0445。然后，确定相位信息，根据图3(c) MATLAB工作区数据可知，10 Hz频率的相位为−2.10150。
依据以上信息得到高频分量离散信号的解析表达式(2)，为计算简便，相位参数均做了取整，直流分量和振幅参数保留2位小数。
(2)
同理，可以分别得到低频分量和趋势项的离散信号解析表达式X( )、X( )，见如下表达式(3)、(4)。
(3)
(4)
对得到的三个分量的函数解析式进行加总，获得全国碳市场交易价格原始序列的函数解析式X(t)，见如下表达式(5)。
(5)
(四) 碳价格频域解析预测结果分析
式(5)滤除了时域分析框架下难以量化的随机扰动因素，清晰反映了碳价格波动规律与碳价格影响因素之间的关系，用于预测碳价格，具体结果如下表4。
Table 4. Prediction results of the frequency domain analytical model of the national carbon market trading price
表4. 全国碳市场交易价格频域解析模型预测结果
根据表4可以看出，就短期预测而言，频域解析模型的预测值与真实值比较接近，相对误差较小，预测结果较为理想。
为进一步论证频域解析模型的预测能力，选取常用的GARCH/ARIMA时域模型、EEMD-GARCH时频域模型以及BP/EEMD-BP神经网络与频域集成模型预测碳市场价格，选用均方根误差(RMSE)和平均绝对百分比误差(MAPE)两项指标衡量模型的预测能力，选取相同的估计样本区间和预测样本区间，比较结果如表5所示。
Table 5. Comparison of prediction results of recent research models
表5. 新近研究模型预测结果比较
由表5的各模型预测对比结果可以看出，相较于其他预测模型，在RMSE和MAPE两个评价指标下，频域解析模型预测精度显著提升。就BP神经网络模型而言，BP神经网络等价于表达了从n个自变量(输入节点数)到m个因变量(输出节点数)的函数映射关系，预测能力高度依赖于数据样本，如何确定合适的自变量，如何构建网络并训练出合适的网络参数较为复杂；GARCH和ARIMA模型更多考虑了数据的时变特质，缺乏对变化原因的回溯；EEMD通过分解的方法获得序列的影响因子，在一定程度上弥补了单纯GARCH和ARIMA模型以及BP神经网络模型的不足，因此，EEMD-GARCH模型以及EEMD-BP模型的预测精度较独立使用GARCH或ARIMA或BP预测结果要好，表5数据可以看出，集成模型的RMSE和MAPE两项指标均优于独立的GARCH或ARIMA或BP模型。频域解析模型在EEMD分解的基础上，通过FFT变换，综合了上述GARCH、ARIMA和BP模型的特点，因此，预测结果较为理想，预测能力有提升，是现有预测理论和方法的必要补充。
4. 研究结论
尽管已有大量的关于碳交易价格的研究，但是当前对碳市场价格预测工作的研究主要集中于区域性碳市场价格预测和欧盟等西方经济体碳价格预测，以全国性碳市场为研究对象的碳价格预测鲜有着墨。本文利用CEEMDAN分解方法将全国碳市场原始价格序列分解为8个IMF序列和趋势项，分解后的各序列在不同时间尺度上反映了影响全国碳市场价格波动的内在因素。分析发现，由IMF1、IMF2、IMF3、IMF4、IMF5和IMF6构成的高频序列所代表的短期投机、短期供需不平衡以及天气状况等周期较短的因素对碳市场价格的波动影响不显著。由IMF7和IMF8构成的低频序列以及趋势项所代表的碳市场内在交易机制和长期价值等因素对碳市场价格的影响占主导地位。
对高、低频序列和趋势项分别进行快速傅里叶变换，估计每个分量的振幅、角频率和相位参数，获取分量的频域解析式，并将其加总，构建全国碳市场交易价格的综合频域解析式，通过解析式预测碳排放权交易价格。实证结果表明，频域解析模型预测法能够有效滤除时间序列数据中的扰动因素，提高了预测精度。
经济环境中的多种因素都会影响碳价格，较大的冲击会导致碳价格发生突变，形成价格波动变点，如何利用人工智能方法，准确预报碳价格序列的变点，有效管理碳金融产品的市场风险将是未来关注的主要问题。
基金项目
教育部人文社会科学研究规划基金项目“频域分析视角下金融系统极值风险防范与化解机制研究”(19YJA630023)。
NOTES
*通讯作者。