2021年电工杯数学建模
B题 光伏建筑一体化板块指数发展趋势分析及预测
原题再现：
国家《第十四个五年规划和 2035 年远景目标纲要》中提出，将 2030 年实现“碳达峰”与 2060 年实现“碳中和”作为我国应对全球气候变暖的一个重要远景目标。光伏建筑一体化(BIPV)是充分利用工业建筑、公共建筑屋顶等资源实施分布式光伏发电工程，它对我国实现“碳达峰”和“碳中和”起到重要作用。目前已有北京、天津、上海、重庆、内蒙古、浙江等 31 地发布光伏建筑一体化未来三至五年相关政策，这将对光伏建筑一体化相关上市企业的发展带来一定影响。在股票交易软件中，每支股票都有日 K 线（由开盘价、收盘线、最高价、最低价构成），还有移动平均线（5 日、10 日、20 日等），通过对日 K 线和移动平均线的分析，可知每支股票的走势。股票市场中有板块指数（将同一板块的个股按不同的权重方式生成相关指数），它是对该板块走势的整体反映。现统计沪深股市中 37 家光伏建筑一体化相关企业的股票数据（见附件 1），将这 37 家企业股票看作一个整体称为光伏建筑一体化板块。光伏建筑一体化板块作为新兴板块，通过对板块指数发展趋势做出预测，可以了解板块相应行业未来的发展趋势。
请根据所给数据资料，解决以下问题：
1.利用附件 1 中数据，给出光伏建筑一体化板块指数的移动平均线（5 日、10 日、20 日等）模型，并绘制 2019 年 4 月 1 日至 2021 年 4 月 30 日该板块指数的移动平均线。
2.利用 2021 年 5 月 6 日至 5 月 28 日数据，对所建立模型进行误差分析并修正模型，根据修正后的模型对该板块未来发展趋势做出预测，给出 5 月 28 日后20 个交易日的日移动平均线、3 周的周移动平均线、2 个月的月移动平均线。
3.利用 2019 年 4 月 1 日至 2021 年 5 月 28 日数据，以 2 个月为一个时间段，对上证指数和光伏建筑一体化板块指数进行相关性分析。
4.对光伏建筑一体化板块个股投资风险进行评估，给出该板块 37 支股票2021 年 6 月份投资风险由低到高的排序结果。假定投资者持有资金 100 万元人民币，欲全部用于该板块的投资，请给出 2021 年 6 月份的最优投资方案（当日可用于投资的资金为上一个交易日结束后投资者所持有的资金，且每日持股数不超过 5 支）。
5.请给相关部门写一份不少于 1500 字关于我国光伏建筑一体化行业未来发展趋势的报告。
附件 1：沪深股市中 37 家光伏建筑一体化相关企业的股票数据
附件 2：http://www.csindex.com.cn/zh-CN/indices/index-rules中证指数有限公司
整体求解过程概述(摘要)
光伏建筑一体化指的是利用公共建筑、工业建筑的顶层布置光伏发电板，利用这些闲置资源来进行分布式光伏发电，以此实现“碳中和”目标。因此，对光伏建筑一体化的研究显得极为重要，本文针对光伏建筑一体化板块指数的发展趋势，建立了相关模型进行分析与预测。
针对问题一，首先建立了光伏建筑一体化板块指数的数学模型，通过平均加权求和的方式，将该板块 37 只股票的加权平均和作为板块指数。为了计算出板块指数的移动平均线，首先将板块指数的收盘价数据进行整理，如图 2 所示。据此绘制了板块指数的 5 日、10 日、20 日的移动平均线，如图 3 所示。然后以 5日移动平均线为例进行分析，建立了多幂次多项式数学模型，见式(3)。最后采用递推最小二乘法进行迭代求解，结果如图 4 所示，得到最终的模型见式(7)。
针对问题二，首先通过问题一模型对 2021 年 5.6~5.28 的移动平均线进行预测，然后与实际数据进行对比分析，结果如图 5 所示，此时平均误差为 1.2。然后通过 2021 年 5.6~5.28 时间段的数据对模型进行修正，同时修改模型的最大幂次，结果如图 6 所示，此时的平均误差为 0.2329，相比模型修正前的误差，修正后的平均误差同比下降了 80.6%。然后通过修正后的模型对未来数据进行预测，结果如图 7 所示。
针对问题三，由于上证指数和光伏建筑一体化指数的数值差异较大，因此首先对其进行归一化处理，然后利用加权平均的方式，建立了上证指数和光伏建筑一体化指数各自的综合指数模型，并对其进行求解，结果如图 8 所示。最后采用pearson 系数对二者指数进行相关性分析，结果如图 9 所示，其中 13 个时间段中，有 2 个时间段高度正相关；3 个时间段中度正相关；3 个时间段低度正相关；4个时间段不相关；1 个时间段中度负相关。
针对问题四，首先将将开盘价、最高价、最低价、收盘价数值与上一交易日的数值的差作为评价指标，然后以指标平均和的负数作为标准建立了综合评价模型，见式(12)。通过该模型对光伏一体化板块的 37 只股票进行了风险平均，结果如图 10 所示。其中，风险最低的 5 只股票依次是隆基股份、森特股份 、南玻A、中来股份 、 苏美达。最后通过风险程度，给出了投资者的最优投资方案，如表 2 所示。
针对问题五，根据本文研究结果和资料查询，给出了关于我国光伏建筑一体化行业未来发展趋势的报告。
模型假设：
1、假设搜集到的上证指数数据和附件所给数据真实有效。
2、假设数据预测时，不考虑其他突变因素的影响
3、假设 BIPV 板块的个股企业在预测时间内不会退市。
问题分析：
问题一分析
针对问题一，需建立出光伏建筑一体化板块指数移动平均线模型，首先，根据已知数据采用平均加权求和的方式得到板块指数中各变量的模型，进而可以得到板块指数中的收盘价，由题需绘制 5、10、20 日移动平均线，故分别以 5、10、20 为采样周期，得到对应的移动平均线，在此基础上，进一步建立多幂次多项式数学模型来表示板块指数移动平均线模型，通过求解得到该模型系数结果。
问题二分析
针对问题二，需利用 2021 年 5 月 6 日至 5 月 28 日的数据，来分析问题一中所建模型的误差并进行修正，由于该时间段内的数据较少，利用 10 日、20 日移动平均线模型进行预测时可验证的数据不足，因此，本问主要考虑利用问题一中所建立的 5 日移动平均线模型进行预测，通过与实际值进行对比，从而可以得到模型的预测误差，完成误差分析，根据题意得，需利用 5 月 6 日至 5 月 28 日的数据对模型进行修正，因此，考虑将此段时间内实际数据代入问题一所建立的多幂次多项式数学模型中，可以利用递推最小二乘法重新求得最优系数，完成对模型的修正，进而利用修正后的模型重新进行预测。
问题三分析
针对问题三，需分析上证指数和光伏建筑一体化板块指数的相关性，根据相关性定理，需先分别求得两者指数的一个代表性指标，再利用相关系数计算公式在求得相关性，考虑到两者指数中包括多个指标且每个指标差异较大，因此，先分别将两者的各个指标进行归一化处理并加权求和，来得到两者的综合指数，根据题目要求，需以两个月为一个时间段，故将总时间分成 13 段，分别进行相关性分析。
问题四分析
针对问题四，需对板块中 37 支股票进行风险评估，本问考虑当股票波动较大时，其风险也会较大，因此选取个股指数中的开盘价、最高价、最低价以及收盘价作为评价指标，但考虑到每个指标差异较大，故重新选取该指标当前时刻数值与上一时刻数值的差值作为新的评价指标，在此基础上进行加权平均，得到个股的风险评价模型，由于该模型主要包括个股指数两个时刻的差值，即变化量，因此，判断模型计算值越大，则风险越高，反之越低，最后选取风险评价值最低的 5 支按比例分配投资资金即可。
问题五分析
针对问题五，通过查阅相关资料，了解光伏建筑一体化的概念以及发展现状，根据本文对于该板块股市的分析，进一步分析该板块未来发展趋势，并提出相关建议。
模型的建立与求解整体论文缩略图
全部论文请见下方“ 只会建模 QQ名片” 点击QQ名片即可
程序代码：
部分程序如下:
%先求时延和嵌入维，带入本程序，Lyapunov_rosenstein_2.dll调用需要安装完整Simulink功能
% 计算时间序列 Lyapunov 指数
clc
clear
close all
data = textread('深证.txt'); % 时间序列，列向量
data = data(end-2000:end);
%-----------------------------------------------------------------
% Logistic 入口参数
tau = 2; % 时延
m = 25; % 嵌入维
[xn] = PhaSpaRecon(data,tau,m); % 每列为一个点
xn = xn'; % 每行为一个点
P =fix(period_mean_fft(data)); % 序列平均周期
if P>500P=500;
end
taumax = 500; % 最大离散步进时间
fs = 1; % 采样频率
Y = Lyapunov_rosenstein_2(xn,fs,tau,m,taumax,P);
Y(1)=[];
figure
plot(log(Y),'b-'); grid; xlabel('i'); ylabel('y(i)');
linear_zone=[1:length(Y)]';
F = polyfit(linear_zone,log(Y),1);
y=linear_zone.*F(1)+F(2);
hold on
plot(linear_zone,y,'r-')
legend('data','线性')
str = ['y=',num2str(F(1)),'x+(',num2str(F(2)),')'];
text(50,max(y),str)
clear
clc
%cao法求嵌入维
data = textread('深证.txt'); % 时间序列，列向量
data = data(end-100:end);
min_m=1;
max_m=30;
tau=2;
k=1;
N=length(data);
En=[];
for m=min_m:max_mNm=N-tau*(m-1);Y=[];D=[];ad=[];Y=reconstitution(data,N,m,tau);for i=1:N-m*taud=[];for j=1:N-m*taud(j)=norm(Y(:,i)-Y(:,j),inf);endtemp=[];temp=sort(d);D(i,1)=i; temp1=[];temp2=[];temp1=find(temp>0);temp2=find(d==temp(temp1(1)));D(i,2)=temp2(1); D(i,3)=temp(temp1(1));%计算a(i,m)Y1=[];Y2=[];Y1=[Y(:,i);data(m*tau+i)];Y2=[Y(:,D(i,2));data(D(i,2)+m*tau)];ad(i)=norm(Y1-Y2,inf)/D(i,3);end%求E(d)E(k,1)=m;E(k,2)=sum(ad)/(N-m*tau);%求E*(d)En(k,1)=m;dd=[];for kk=1:N-m*taudd(kk)=abs(data(D(kk,1)+m*tau)-data(D(kk,2)+m*tau));endEn(k,2)=sum(dd)/(N-m*tau);k=k+1;
end
%求E1(d)
E1=[];
for i=1:(max_m-min_m)E1(i,1)=E(i,1);E1(i,2)=E(i+1,2)/E(i,2);
end
%求E2(d)
E2=[];
for i=1:(max_m-min_m)E2(i,1)=En(i,1);E2(i,2)=En(i+1,2)/En(i,2);
end
figure(1)
plot(E1(:,1),E1(:,2),'-bs',E2(:,1),E2(:,2),'-r*');xlabel('嵌入维数');ylabel('E1(m)&E2(m)');
grid on
function [xn,dn] = PhaSpaRecon(s,tau,m,T)
% 混沌序列的相空间重构 (phase space reconstruction)
% [xn, dn, xn_cols] = PhaSpaRecon(s, tau, m)
% 输入参数： s 混沌序列(列向量)
% tau 重构时延
% m 重构维数
% T 直接预测步数
% 输出参数： xn 相空间中的点序列(每一列为一个点)
% dn 一步预测的目标(行向量)[rows,cols] = size(s);
if (rows>cols)len = rows;s = s';
elselen = cols;
endif (nargin < 4)T = 1;
endif (nargout==1)if (len-(m-1)*tau < 1)disp('err: delay time or the embedding dimension is too large!')xn = [];elsexn = zeros(m,len-(m-1)*tau);for i = 1:mxn(i,:) = s(1+(i-1)*tau : len-(m-i)*tau); % 相空间重构，每一行为一个点 endendelseif (nargout==2)if (len-T-(m-1)*tau < 1)disp('err: delay time or the embedding dimension is too large!')xn = [];dn = [];elsexn = zeros(m,len-T-(m-1)*tau);for i = 1:mxn(i,:) = s(1+(i-1)*tau : len-T-(m-i)*tau); % 相空间重构，每一行为一个点 enddn = s(1+T+(m-1)*tau : end); % 预测的目标endend
