【太阳黑子预测】太阳黑子变化规律预测matlab仿真
1.软件版本
matlab2013b
2.本算法理论知识
太阳黑子是人们最早发现也是人们最熟悉的一种太阳表面活动。因为太阳内部磁场发生变化,太阳黑子的数量并不是固定的,它会随着时间的变化而上下波动,每隔一定时间会达到一个最高点,这段时间就被称之为一个太阳黑子周期。太阳黑子的活动呈现周期性变化是由施瓦贝首次发现的。沃尔夫 (R.Wolfer)继而推算出11年的周期规律。实际上,太阳黑子的活动不仅呈11年的周期变化,还有海耳在研究太阳黑子磁场分布时发现的22年周期;格莱斯堡等人发现的80年周期以及蒙德极小期等。由于太阳黑子的活动规律极其复杂,时至今日科学家们仍在努力研究其内在的规律和特性。事实上,对太阳黑子活动规律的研究不仅具有理论意义,而且具有直接的应用需求。太阳黑子的活动呈现周期性变化的,沃尔夫(R.Wolfer)根据在过去的288 年(1700年~1987 年)间每年太阳黑子出现的数量和大小的观测数据推算出11 年的周期规律。我们利用Matlab强大的数据处理与仿真功能,对Wolfer数进行功率谱密度分析从而可以得到对太阳黑子活动周期的结论。
使用的线性模型为
这个部分的MATLAB代码为:
现象模型的预测结构如下:
将上面的图放大之后,得到的局部估计效果如下图所示:
然后通过多变量的最小二乘法进行模型预测。引入多个变量作为参数估计,这里,将最小二乘法的估计函数为:
由于从历年的数据可知,虽然总体上数据呈现出周期波形,但是对于每个周期中的波形,其中的局部图中会出现部分高频周期量,所以在估计的时候,这里加上更高频率的分量。使用和上面所讲的最小二乘法进行上述参数的估计。
这个步骤,得到的仿真结果如下图所示:
3.部分源码
clc;
clear;
close all;
warning off;load Sunspot.txt
%YEAR MON SSN DEV%画出导入的活动数据
YEAR = Sunspot(:,1);
MON = Sunspot(:,2);
SSN = Sunspot(:,3); %sun spot number
DEV = Sunspot(:,4); %标准偏差表figure;
subplot(211);plot(SSN,'b.');
title('SSN');subplot(212);plot(DEV,'b.');
title('DEV');addpath 'funcs\'%step1:Devise and/or employ methods to ?nd the best frequency for the sunspot cycle.
%step1:Devise and/or employ methods to ?nd the best frequency for the sunspot cycle.
%step1:Devise and/or employ methods to ?nd the best frequency for the sunspot cycle.%找到一种较好的方法计算周期
K = 20;%前22年数据去除
Start = 12*K+1;
mainPeriod = func_cycle_cal(SSN(Start:end));
mainPeriod%Step2:Predict the time of the next solar maximum and minimum.
%Step2:Predict the time of the next solar maximum and minimum.
%Step2:Predict the time of the next solar maximum and minimum.
%构建预测模型,预测后面的最大值和最小值对应的时间
%构建预测模型,预测后面的最大值和最小值对应的时间
c0 = [1 1 1 1 1 1 1 1]';%直接给出被辨识参数的初始值,即一个充分小的实向量
p0 = 10^6*eye(8,8); %直接给出初始状态P0,即一个充分大的实数单位矩阵
t = 1/(12):1/(12):length(SSN)/(12);
yc = cos(2*pi*t/mainPeriod);
ys = sin(2*pi*t/mainPeriod);
yc2 = cos(2*pi*t/(2*mainPeriod));
ys2 = sin(2*pi*t/(2*mainPeriod));
yc3 = cos(2*pi*t/(3*mainPeriod));
ys3 = sin(2*pi*t/(3*mainPeriod));for k=1:length(SSN); %开始求K h1 =[1,k,yc(k),ys(k),yc2(k),ys2(k),yc3(k),ys3(k)]'; x = h1'*p0*h1 + 99; x1 = inv(x); %开始求K(k) k1 = p0*h1*x1; %求出K的值 d1 = SSN(k)-h1'*c0; c1 = c0+k1*d1; %求被辨识参数c e1 = c1-c0; %求参数当前值与上一次的值的差值 e2 = e1./c0; %求参数的相对变化 e(:,k) = e2; %把当前相对变化的列向量加入误差矩阵的最后一列 c0 = c1; %新获得的参数作为下一次递推的旧参数 c(:,k) = c1; %把辨识参数c 列向量加入辨识参数矩阵的最后一列 p1 = p0-k1*k1'*[h1'*p0*h1+1]; %求出 p(k)的值 p0 = p1; %给下次用
end
%估计得到的太阳黑子活动曲线
for k=1:length(SSN); %开始求K Y_predict2(k) = c(1,k) + c(2,k)*k + c(3,k)*yc(k) + c(4,k)*ys(k) + c(5,k)*yc2(k) + c(6,k)*ys2(k) + c(7,k)*yc3(k) + c(8,k)*ys3(k);
end
figure;
plot(Y_predict2,'b','LineWidth',2);hold on;
plot(SSN,'r');hold off;
legend('预测SSN','实际SSN');Predict_Len = 100;%定义预测时间长度t = 1/(12):1/(12):(length(SSN)+Predict_Len)/(12);
yc = cos(2*pi*t/mainPeriod);
ys = sin(2*pi*t/mainPeriod);
yc2 = cos(2*pi*t/(2*mainPeriod));
ys2 = sin(2*pi*t/(2*mainPeriod));
yc3 = cos(2*pi*t/(3*mainPeriod));
ys3 = sin(2*pi*t/(3*mainPeriod));
ind = 0;for k=1:length(SSN)+Predict_Len; %开始求K if k <= length(SSN)Y_predict3(k) = c(1,k) + c(2,k)*k + c(3,k)*yc(k) + c(4,k)*ys(k) + c(5,k)*yc2(k) + c(6,k)*ys2(k) + c(7,k)*yc3(k) + c(8,k)*ys3(k); else %Y_predict3(k) = c(1,end) + c(2,end)*k + c(3,end)*yc(k) + c(4,end)*ys(k) + c(5,end)*yc2(k) + c(6,end)*ys2(k) + c(7,end)*yc3(k) + c(8,end)*ys3(k); c0 = mean(c(1,end:end));c1 = mean(c(2,end:end));c2 = mean(c(3,end:end));c3 = mean(c(4,end:end));c4 = mean(c(5,end:end));c5 = mean(c(6,end:end));c6 = mean(c(7,end:end));c7 = mean(c(8,end:end));Y_predict3(k) = c0 + c1*k + c2*yc(k) + c3*ys(k) + c4*yc2(k) + c5*ys2(k) + c6*yc3(k) + c7*ys3(k); ind = ind + 1;Ys(ind) = Y_predict3(k);end
end
figure;
plot(Y_predict3,'b','LineWidth',2);hold on;
plot(SSN,'r');hold off;
legend('预测SSN','实际SSN');
grid on;%根据预测结果得到下次太阳黑子活动高峰和低峰的时间
%前一次高峰日期为XX = 59;
[Vmax1,Imax1] = max(Ys);
[Vmax2,Imax2] = max(SSN(length(SSN)-XX:length(SSN)));%3100~3160if Vmax1 > Vmax2II = Imax1;MM = Vmax1; time = (length(SSN) + II-3019);%原数据的最后一个月份+预测后的最大值 - 前一个高峰日期
elseII = Imax2;MM = Vmax2; time = (length(SSN) + (XX-II)-3019);%原数据的最后一个月份+预测后的最大值 - 前一个高峰日期
end
Years=time/12;fprintf('下次高峰期日期为:%d',round(2000 + Years));
fprintf('年\n\n');
fprintf('最大值为:%2.2f\n\n\n\n',MM);%计算下一次低谷值
[Vmin,Imin] = min(Ys);
fprintf('下次低峰期日期为:%d',round(2012 + Imin/12));
fprintf('年\n\n');
fprintf('最小值为:%2.2f\n\n',Vmin);
A16-13
【太阳黑子预测】太阳黑子变化规律预测matlab仿真相关推荐
- m基于ESN+BP神经网络的数据预测算法matlab仿真,测试数据为太阳黑子变化数据
目录 1.算法描述 2.仿真效果预览 3.MATLAB核心程序 4.完整MATLAB 1.算法描述 在人工神经网络的发展历史上,感知机(Multilayer Perceptron,MLP)网络曾对人工 ...
- 基于BP神经网络的多因素房屋价格预测matlab仿真
目录 一.理论基础 二.案例背景 1.问题描述 2.思路流程 三.部分MATLAB仿真 四.仿真结论分析 五.参考文献 一.理论基础 神经网络主要由处理单元.网络拓扑结构.训练规则组成.处理单元是神经 ...
- 风能matlab仿真_风能产量预测—深度学习项目
风能matlab仿真 DL DATATHON- AI4Impact DL DATATHON- AI4影响 Published by Team AI Traders - Suyash Lohia, Ng ...
- m分别使用BP神经网络和GRNN网络进行时间序列预测matlab仿真
目录 1.算法描述 2.仿真效果预览 3.MATLAB核心程序 4.完整MATLAB 1.算法描述 广义回归神经网络是径向基神经网络的一种,GRNN具有很强的非线性映射能力和学习速度,比RBF具有更强 ...
- 时间序列深度学习:状态 LSTM 模型预测太阳黑子
目录 时间序列深度学习:状态 LSTM 模型预测太阳黑子 教程概览 商业应用 长短期记忆(LSTM)模型 太阳黑子数据集 构建 LSTM 模型预测太阳黑子 1 若干相关包 2 数据 3 探索性数据分析 ...
- 极限学习机(Extreme Learning Machine, ELM)的训练与预测matlab仿真
目录 1.算法概述 2.仿真效果 3.MATLAB仿真源码 1.算法概述 极限学习机(ELM)是当前一类非常热门的机器学习算法,被用来训练单隐层前馈神经网络(SLFN
- 基于GA优化RBF神经网络(GA-RBF)数据预测的matlab仿真
目录 1.算法概述 2.仿真效果 3.MATLAB仿真源码 1.算法概述 遗传算法的基本运算过程如下: (1)初始化:设置进化代数计数器t=0,设置最大进化代数T,随机生成M个个体作为初始群体P(0 ...
- 基于IMM和UKF的三维路径预测跟踪matlab仿真,模型包括匀速模型CV,匀加速模型CA以及常速率协同转弯模型CSCT
目录 1.算法概述 2.仿真效果 3.MATLAB仿真源码 1.算法概述 IMM(交互式多模型)方法是Blom H.A.P.于1984年提出的.多模型方法主要用于特性随时问变化系统的状态估计,所以它特 ...
- 用matlab进行markov链预测,用MATLAB仿真markov链程序
<用MATLAB仿真markov链程序>由会员分享,可在线阅读,更多相关<用MATLAB仿真markov链程序(2页珍藏版)>请在人人文库网上搜索. 1.用MATLAB仿真ma ...
- 【负荷预测】基于神经网络的负荷预测和价格预测(Matlab代码实现)
目录 1 概述 2 基于神经网络的负荷预测(Matlab实现) 2.1 代码 2.2 结果 2.3 回归树模型的进一步改进 3 基于神经网络的价格预测(Matlab代码实现) 4 阅读全文(Ma ...
最新文章
- RateLimiter 的底层实现是啥?
- 关于html5和css3的新特性
- 【Uva 10934】Dropping water balloons
- 如何用Java创建不可变的Map
- lodash和debounce
- vue项目中更新element-ui版本
- [Swust OJ 188]--异面空间(读懂题意很重要)
- 帆软填报报表的使用教程
- C# ABB机器人上位机控制 .net PC SDK开发全流程(通信、控制、日志、二次开发)--Chapter 1
- 车辆管理系统实施方案
- HTML中font标签中size属性值对应的像素大小
- Robotics Toolbox :(1)建立机器人模型
- MoFlow:生成分子图的可逆流模型
- Windows server 2012 R2添加桌面图标(计算机、控制面板、网络等)
- 多语言适配分享会演讲稿
- 无需设置路由器,无需公网ip 实现永久免费内网穿透
- redis incr mysql_redis命令_INCR
- 电商app源码该如何布局流量模块(下)
- 毛球科技点评区块链在组织和商业生态系统中的未来
- Go语言实战 - 我需要站内搜索
热门文章
- 一道关于扔球/扔鸡蛋/摔手机的DP问题(蓝桥杯题目/面试题)
- javascript高级程序设计读书笔记----引用类型
- 样本分布不平衡处理策略(20210429)
- 【Java基础】NoClassDefFoundError 和 ClassNotFoundException的定义及其区别
- 蓝牙车载系统的组成结构和应用规范分析
- 3月23—3月27日三年级课程
- lemke算法 matlab,lemke是什么意思
- python之正则表达式及RE模块
- 西部矿业(601168):整合湖北铅锌资源
- https://gns3.com/community/discussion/gns3-doesn-t-work-on-vmware-play