“华为杯”竞赛论文模板

竞赛论文模板Realthought

题目降低汽油精制过程中的辛烷值损失模型

摘要：

本文主要针对汽油辛烷值进行量化分析，依据从催化裂化汽油精制装置采集的325个数据样本（每个数据样本都有354个操作变量），通过降维分析法找到主要变量，并且借助数据挖掘技术来建立汽油辛烷值（RON）损失的预测模型，找出每个样本的优化操作条件，在汽油产品脱硫效果得到有效保证，并且满足产品硫含量不大于5μg/g的前提要求下，尽量降低汽油辛烷值损失在30% 以上。

针对问题1，原始数据采集来自于中石化高桥石化实时数据库（霍尼韦尔PHD）及LIMS 实验数据库,其中操作变量数据来自于实时数据库，特点是变量较多且数据采集时间长，导致数据精确处理难大，特别地会有空值或者是部分为空值，所以首先要对原始数据做一些处理后方可使用。

针对问题2，由于催化裂化汽油精制过程是连续的，虽然操作变量每3 分钟就采样一次，但辛烷值（因变量）的测量比较麻烦，一周仅2次无法对应。但根据实际情况可以认为辛烷值的测量值是测量时刻前两小时内操作变量的综合效果，因此预处理中取操作变量两小时内的平均值与辛烷值的测量值对应，从367个操作变量中通过降维的方法筛选出建模主要变量，使之尽可能具有独立性、代表性，特别地将原料辛烷值纳入变量，从而找到建模地主变量。

针对问题3，采用上述样本和建模主要变量，将数据分为训练集和测试集，通过深度学习技术进行建模，进而构建深层神经网络并进行模型验证。

针对问题4，为寻找最优操作变量，我们通过构建粒子群算法，找到可接受的最优辛烷值的操作变量参数。

针对问题5，将问题4的过程进行可视化即可。

_

关键词：汽油辛烷值；操作变量；分析降维法；深层神经网络；预测模型；操作变量优化；模型可视化展示；粒子群算法。

目录

一，问题重述 (2)

1、背景叙述 (2)

2、数据采集 (3)

3、问题解析 (3)

3、1问题1，处理数据 (3)

3、2问题2，寻找建模主变量 (6)

3、3 问题3，建立辛烷值（RON）损失预测模型 (8)

3、4 问题4，主要操作变量方案的优化 (9)

3、5问题五，模型的可视化展示 (13)

二、问题结论 (13)

1、问题1、2：数据处理结果 (13)

2、问题3，建立辛烷值（RON）损失预测模型 (13)

3、问题4、5，操作方案的优化及模型可视化 (13)

三、参考文献 (14)

一，问题重述

1、背景叙述

在世界的各个角落，小型车辆的主要燃料是汽油，然而汽油燃烧产生的尾气排放大气中，对大气环境有不可忽略的重要影响。因此，在各个国家切合本国国情并且日益严格的汽油质量标准相继被拟定出来。汽油清洁化重点是尽量保持其辛烷值的同时，降低汽油中的硫、烯烃含量。

图1

中国作为全球主要石油消费大国之一，其中超过70%原油要依赖进口，且从中东地区进口的大部分原油都不同程度含有较高的硫。原油中的重油通常占比40-60%，这部分重油（以硫为代表的杂质含量也高）很难充分直接利用。为了加大重油资源的利用效率，我国

大力发展了以催化裂化为核心的重油轻质化工艺技术，以此将重油转化为相对环保的汽油、柴油和低碳烯烃等可使用油。在我国，由催化裂化生产得到超过70% 的汽油，因此催化裂化汽油提供了成品汽油中95% 以上的硫和烯烃。因而有必要对催化裂化汽油进行进一步的精制化处理，从而达到对汽油辛烷值达到国家标准或者更优。

辛烷值（以RON表示）是对于汽油燃烧性能最重要的指标，并且汽油的商品牌号（例如89#、92#、95#）以辛烷值作为标准。然而，现有以催化裂化汽油的方法进行脱硫和降烯烃过程中，汽油辛烷值普遍偏低。根据目前市场价来计算，每降低1个单位的辛烷值，相当于约150 元/吨的损失。以一个100 万吨/年催化裂化汽油精制装置为例，若能降低RON损失0.3个单位，其经济效益损失将达到四千五百万元。

一般通过数据关联或机理建模的方法来实现化工过程的建模，虽然已经取得了不错的成果，但是出于炼油工艺过程的复杂性以及设备的多样性的原因，使得它们的操作变量（控制变量）之间有高度非线性和相互强耦联的关系存在，然而传统的数据关联模型中存在变量相对较少存在，使得机理建模对原料的分析要求较高，同时过程优化的响应并不及时，因而并未达到理想的效果。

某石化企业运行了4年催化裂化汽油精制脱硫装置，积累了大量可用的原始数据，有平均为1.37个单位的汽油产品辛烷值损失，然而同类装置的最小损失值仅仅只有0.6个单位，因而优化空间比较大。本文将探索利用数据挖掘技术来解决化工过程建模问题，通过模型可视化方法，从而达到找出优化汽油辛烷值的目的。

2、数据采集

从LIMS实验数据库及中石化高桥石化实时数据库（霍尼韦尔PHD），我们得到大量的原始数据，其中实时数据库为我们提供了大量的操作变量数据，采集时间从2017年4月至2020年5月，时间跨度大约为3年，采集操作位点数共354 个。从2017 年4月到2019年9月，数据采集频次为3 分钟/次；自2019 年10 月至2020 年5 月，数据采集频次为6 分钟/次。原料、催化剂和产品数据来源于LIMS的实验数据库，数据时间范围为3年，即从2017 年4 月至2020年5 月。我们以每周2/次的采集频次，采集到重要的建模变量为原料及产品的辛烷值这两样数据。

3、问题解析

3、1问题1，处理数据

在原始数据中，大部分变量数据是正常的，但每套装置的数据难免都有部分位点会存在问题，部分变量只含有某时间段的数据，部分变量的数据部分为空值或全部数据为空值。因此，处理后的原始数据方可以使用，我们对数据进行剔除或者补全两种办法进行预处理。数据处理方法如下：

步骤一：假如325个样本中数据全部为空值的位点，则表明对该问题毫无帮助，古我们只能选择删除。

步骤二：如果某个点位的缺失个数占总数据量的比值达到给定阈值，则将此点位删除。在此假设下，情况（2）当归属于情况（1）。根据数据取值范围信息，绝大部分操作变量取值范围不包含0 值，但亦有个别操作变量可取0值。有鉴于此，我们假定，如果某个操作变量随时间的变化量取0值的个数占总量的比值也超过给定阈值，就把该操作变量当作不可置信的数据；另外，如果某个点位超出范围的数据个数占总数据量的比值达到给定阈值，则认为该点位置信度较差，将其删除。

此外，对于缺失值，我们也可以用缺失值比率的方法进行量化。该方法基本的思想就是基于包含太多缺失值的数据列，则该数列包含有用信息的可能性也就应当越低。所以，我们可以去掉数据缺失值大于某个阈值的列。阈值越大，降维方法效果更佳。一般地，数据缺失值大于90%的时候，我们可以考虑把这列去掉。在某些特定的行业如信用评分建模时，我们可以考虑把缺失值大于90%的变量设置成特殊规则，这里我们考虑使用pandas 包的dropna 方法。

步骤三：对于部分数据为空值的位点，这里我们可以利用前后面的数值，然后用拉格朗日插值法，给空值位点赋值。这里我们介绍一下拉格朗日值法在处理操作变量上的应用：

命题1：假设有函数 是互不相等的一系列点这里的为有限个数列我们将称为关于的一阶差商（也称均差），可以记为即j i j i j i x x x f x f x x f --=

)()(],[。下面的式子，我们称为一阶

差商的差商，即

为)(x f 关于k j i x x x ,,的二阶差商，记为],,[k j i x x x f 。一般地，称

是)(x f 关于k x x x ,...,,10的k 阶差商。定义了差商之后，我们发现了如下性质，有点像对称矩阵：

此时我们要开始对公式进行预测，发现一阶是线性插值，然后可以逐次往高阶迭代。线性插值表示为：

称为一次Newton公式。一般地，由各阶差商的定义，我们可的

把以上式子从下往上依次带入，我们就得到一下多项式，

。这样我们就可以得到插值多项式和误差函数，及插值多项式为

误差函数为

步骤四：根据质量要求与操制经验，原始数据变量的操作范围需要被提取出来，然后采用最大最小的限幅方法找出一部分不在此范围的样本,然后予以剔除，具体操作就是如果某个点位超出范围的数据个数占总数据量的比值达到给定阈值，则认为该点位置信度较差，

将其删除。

步骤五；去除异常值，我们用拉依达准则（3σ准则）。

具体操作如下：首先对被测量的数据进行处理，这里我们用精准测量的方法。得到n x x x ,...,,21，把它们的算术平均值算出来，将i x x -=i v 作为剩余误差，并且σ由贝塞尔公式算出。如果某个位点的剩余误差b v 与其测量值b x 达到σ3||||≥-=b b x x v 的条件，则认为b x 是含有较大误差，视为不好的值，应将其删除，或者用平均值代替。贝塞尔公式如下：

。

3、2问题2，寻找建模主变量

由于催化裂化汽油精制过程是连续的，虽然操作变量每3 分钟就采样一次，但辛烷值（因变量）的测量比较麻烦，一周仅能完成至多完成两次测量。但根据实际情况可以认为辛烷值的测量值是测量时刻前两小时内操作变量的综合效果，因此预处理中取操作变量两小时内的平均值与辛烷值的测量值对应这样产生了325个样本。

建立降低辛烷值损失模型涉及包括7个原料性质、2个待生吸附剂性质、2个再生吸附剂性质、2个产品性质等变量以及另外354个操作变量（共计367个变量），工程技术应用中经常使用先降维后建模的方法，这有利于忽略次要因素，发现并分析影响模型的主要变量与因素。因此，根据提供的325个样本数据，通过降维的方法从367个操作变量中筛选出建模主要变量，使之尽可能具有 独立性 、代表性，为了工程应用方便，建议降维后的主要变量会控制在在30 个以下，并且将原料的辛烷值作为建模变量之一。

PCA 主成分分析方法可以实现把数据降维到30个维度以内，然后用回归分析的方法，进行曲线拟合，虽然能实现降维处理，但有可能拟合效果较差。

因此我们用数据挖掘技术，采用最小冗余最大相关（mRMR ） 这种滤波特征选择方法使之能够达到选择30 个以内的特征的目的。这里我们将介绍mRMR 的基本原理。

有一种常用的特征选择方法就是分类变量与最大化特征之间的相关度，就是选择与分类变量拥有最高相关度的前k 个变量，然而，在特征选择中，单个好的特征的组合不能达到增加分类器的性能的目的，原因是有可能特征之间是高度相关的，这就导致了特征变量的冗余。这就是Peng et.al 说的“the m best features are not the best m features ”。因而找到一种实现最小化特征与特征之间的相关性或者最大化特征与分类变量之间的相关性的方法显得尤为重要，诞生mRMR 的背景就是如此。

互信息：假设y x ,为随机变量，然后给出他们的概率密度因变量为),(),(),(y x p y p x P ，则互信息为

mRMR 的目的就是找出含有m 个特征子集的S 。

因为本文仅提供离散数据，这里我们仅考虑离散变量，其最大相关性为：),(max c S D ，这里

这里i x 为第i 个特征，c 为类别变量，S 为特征子集。

下面我们介绍最小冗余度：),(min c D R ，满足

然后我们整合最大相关性和最小冗余度，这里我们有两种方式进行整合，第一种为剑法的方式进行整合，即为

第二种以除法的方式进行整合，即为

在具体操作过程中，为达到寻找近似最优的特征的目的，可以通过增量搜索方法来实现。假设1S -m 是我们已有的特征集,我们的目标就是在剩下的特征1S -X -m 中找出第m 个特征，通过特征选择使得Φ最大。下面是增量算法优化的条件：

其算法的复杂度为)

S

O。

(M

,

3、3 问题3，建立辛烷值（RON）损失预测模型

在问题二中，我们已经成功地寻找到了建模主变量，并且实现了建模主变量不超过30个的目标，此时我们将数据分成测试集和训练集，这里可用深度学习技术构建深层神经网络，然后可以成功建模。然而，构建深层学习网络，技术难度较大且耗时较长。为了方便起见，我们先用训练集做回归模型，然后用测试集去测试这个模型。虽然这种效果可能不如深度神经网络，但是短期内代码比较容易实现。

要确定不同变量间的彼此关系，回归分析是一种经常见到的统计方法。我们可以直接俄通过Excel中数据分析菜单中的回归功能快速完成回归建模的过程。这里我们简述一下在paython中借助机器学习库sklearm完成建立简单回归模型的过程。

步骤一;创建回归模型，以下几个库文件是我们需要用到的，分别为，numpy库，pandas 库，sklearn库和matplotlib库。对于数据导入，创建数据表和一些基础的计算工作，我们会用到numpy和pandas库。对于计算线性回归模型中的关键参考数据，我们会用到上述Learn库，而且还可以对模型进行检测。至于matplotlib库，则用于绘制散点图。

步骤二：完成准备工作之后，开始数据的读取，这里我们的数据已经预处理和降维处理，可以使用。将这组数据读取到python中并命名为cost and click，此时就通过head 函数至少可以查看数据表中前5行的内容。

步骤三：在cost_and_click数据表中，我们将已经经过降维处理的操作变量cost设置为自变量X，将汽油辛烷值click设置为因变量Y。希望通过回归模型发现主要操作变量对于汽油辛烷值造成的影响，以及两者间的关系。下面我们将X设置为cost，Y设置为click。并通过shape函数查看了两个变量的行数，可查看25行左右的数据。

步骤四：在输入好自变量X和因变量Y后，为了更清楚的看到两者间的关联，我们可以使用散点图来观察现象，看到不同操作变量对汽油辛烷值的影响。可能操作变量A，B 单独发生更改，都对汽油辛烷值产生比较积极的作用，但A，B同时变化时，可能会产生汽油辛烷值下降的现象，所以就需要对操作变量进行优化，从而达到优化汽油辛烷值的目的。至于怎么优化，优化到何种程度，我们需要从回归模型中找到答案。

步骤五：将这降维后的操作变量的数据分为两部分，回归模型的创建用训练集当中的数据，对回归模型的结果进行检测用测试集当中的数据。例如，将降维后的操作变量数据，使用随机方式抽取40%作为测试集，另外60%作为训练集数据。可以使用shape函数对训练集的数据行数进行查看，查看到的数据将被用于建立回归模型。

步骤六：将训练集数据代入到线性回归模型中，以计算关键的参数值，特别注意，因为自变量较大，这里我们建立的是多元回归模型。

步骤七：将测试集中的自变量代入到多元回归模型中，用以评估回归模型，将结果与要评估模型的准确性，我们将与测试集的因变量进行对比。

步骤八：为了衡量回归模型的计算结果与测试集中因变量的数据差异，我们使用误差平方和对两者的差异进行评估。误差平方和并不能说明回归模型的结果与测试集数据间的差异大小，它只能用来进行模型间准确率的对比。以下是误差平方和的计算公式，为了防止两组数据间差异的正负相互抵消，故取平方和，取绝对值亦可。

3、4 问题4，主要操作变量方案的优化

利用模型获得325个数据样本中，在保证产品硫含量不大于5μg/g的讲求前提下，辛烷值（RON）损失降幅大于30%的样本对应的主要变量优化后的操作条件（优化过程中原料、待生吸附剂、再生吸附剂的性质保持不变，以它们在样本中的数据为准）。这里我们通过粒子群法（PSO），来实现优化方案的抓取。这里我们简述一下PSO的基本原理以及如何实现算法。

在1995年，由于对鸟群捕食行为的研究，以实现获得群体最优解，Eberhart博士和Kennedy博士首次提出了粒子群算法。粒子群算法的基本思想：通过对群体中的个体对信息的共享的使用，使从无序到有序的演化过程在整个群体的运动的求解空间中产生，从而达到获得问题的最优解的目的。

假设这么一个场景：远处有一片森林，而此时一群鸟进行觅食，森林的位置所有的鸟都不知道，但是鸟儿们当前的位置距离森林有多远它们时知道的。那么最简单有效的策略就是搜寻目前距离森林最近的鸟群的周围区域。

在PSO中，每个搜索空间中的优化问题的解是每一只鸟，这里搜索空间每一只鸟称之为"粒子"，即忽略质量和体积，只关心位置和速度，而问题的最优解就对应于鸟群要寻找到的"森林"。所有的粒子都具有速度向量（决定下次飞行的方向和绝对速度）和一个位置向量（粒子在解空间的位置），这里距"森林"的距离用当前的所在位置的适应值（fitness value）来表示，至于适应值可以通过目标函数来计算。在迭代的进程中，鸟群除了机器学习自身的经验（历史位置）以外，还会学习种群中最优粒子的"经验"，这样下一次迭代就有了下一次迭代的根据，即需要如何调整和改变飞行的速度和方向。就会使得整个种群的粒子逐步趋于最优解就，加入如以上逐步迭代下去。为说明名如何寻找最优解，下面将举一个简单的例子。

在一个湖中有两个人他们之间可以通信，并且可以探测到自己所在位置鱼群的密度。

初始位置如上图所示，这里有一条曲线，离水面的深度度代表鱼群的密度，由于右边比较深，因此左边的人会往右边移动一下小船。

现在左边比较深，因此右边的人会往左边移动一下小船。

一直重复该过程，最后两个小船会相遇，就会得到一个局部最优解。

将单个个体视为为粒子，x(t)为每一个单个个体在某一时刻的位置,v(t)为每一个单个个体在某一时刻的速度。

P（t）为在t时刻x个体的自己的最优解，g(t)为在t时刻所有个体的最优解，v(t)为个体在t时刻的方向，x(t)为个体在t时刻的位置。

如上图所示，下一个位置由x,p,g共同决定。为达到找到问题的最优解的目的，种群中的粒子必须不断地向自身和种群的历史信息取得经验。

但是，上述原始的公式中存在一个问题，即公式中V的更新太具有随机性，在后来进行的研究也表明这一点，从而使得整个PSO算法的局部的搜索能力不强，但全局优化能力很强。然而在实际中，我们只需要在算法迭代后期需要PSO有着较强的局部优化能力，而在算法的前期期则整个种群在全局部搜索能力显得尤为重要。由于上述的弊端，shi和Eberhart提出了PSO的惯性权重模型，即是通过引入惯性权重将公式改进，。

迭代公式表示如下：

这里我们称w是PSO的惯性权重，它的取值不超出[0,1]范围，自适应的取值方法在一般应用中会被采用，即一开始令w=0.9，使得PSO较强具有较强的全局优化能力，随着迭代的继续进行，对参数w进行递减，从而使得PSO的局部优化能力得到加强，当迭代终止时，w=0.1。这里参数c1和c2可视为学习因子，一般可设置为1.4961.r1和r2则视为随机概率值且介于[0,1]之间。

我们可以得到粒子群优化算法，其算法框架如下：

注：该图片引用于网络。

Step1：种群初始化，可以进行随机初始化或者根据被优化的问题设计特定的初始化方法，然后计算个体的适应值，从而选择出个体的局部最优位置向量和种群的全局最优位置向每一个单个个体在某一时刻的位置Step2：迭代设置，变换个体的速度向量

Step4：位置变换，变换每个个体的位置向量

Step5：局部位置和全局位置向量变换，变换种群的全局最优解和每个个体的局部最优解。

Step6：终止条件判别，如果满足判别迭代次数时都达到最大迭代次数，输出全局最优解，否则继续进行迭代，跳转至step 3。

对位置向量和速度迭代算子的设计，这就是对于粒子群优化算法的运用。如何设计PSO 的迭代算子是PSO算法应用的研究重点和难点，这是因为迭代算子是否有效将决定整个PSO 算法性能的优劣。

3、5问题五，模型的可视化展示

工业装置为了平稳生产，优化后的主要操作变量（即：问题2中的主要变量）往往只能逐步调整到位，我们将对133号样本（原料性质、待生吸附剂和再生吸附剂的性质数据保持不变，以样本中的数据为准），以图形展示其主要操作变量优化调整过程中对应的汽油辛烷值和硫含量的变化轨迹。这里我们用粒子群算法，实现模型可视化。

二、问题结论

1、问题1、2：数据处理结果

在问题解析3、1中，我们首先对数据进行了预处理，问题解析3、2中，通过数据降维分析法（LASSO）来寻找建模主变量。以下是先进行预处理，然后进行降维的数据结果：

30个操作变量：

2、问题3，建立辛烷值（RON）损失预测模型

采用上述样本和建模主要变量，通过数据挖掘技术建立辛烷值（RON）损失预测模型，并进行模型验证。

很多图片。。。。。。

3、问题

4、5，操作方案的优化及模型可视化

以下是参数参数变量的调整幅度：

很多图片。。。。。。

含硫轨迹图

汽油辛烷值轨迹图

结论分析：从硫含量轨迹图可知，在迭代过程中，硫含量最开始总体平稳，然后开始急剧上升，最后趋于平稳。而辛烷值，随着迭代进行，辛烷值的上升速率稳定在以在一个固定的值附近，所以该模型是优化操作变量，使达到增加辛烷值的的目的。优质汽油要求有较高的汽油甲烷值，同时要有相对较低的硫含量。从可视化图像可知，对于问题4，要求在

保证产品硫含量不大于5μg/g，辛烷值（RON）损失降幅大于30%的样本对应的主要变量

优化后的操作条件可大致在10-16号操作变量数据库中找到。

三、参考文献

[1]Jorgensen P E T, Pedersen S. Spectral pairs in Cartesian coordinates. J Fourier Anal Appl, 1999, 5: 285--302.

[2]kexinxin，粒子群算法，https://https://www.360docs.net/doc/10745078.html,/kexinxin/ p/9858664.html，2020年9月20日。

[3]littlely_ll,特征选择之最小冗余最大相关性(mRMR),https:// https://www.360docs.net/doc/10745078.html,/littlely_ll/article/details/71749776,2020年9月20日。

[4]Pedersen S. Spectral theory of commuting self-adjoint partial differential operators. J Funct Anal, 1987, 73: 12 2--134.

[5]zhangyue_lala,,一些插值法及理解,https://https://www.360docs.net/doc/10745078.html,/zh angyue_lala/article/details/64437547,2020年9月20日。