介绍

研究了好一阵子马欢飞老师在PNAS上发的文章，下面附上个人的研究心得与代码与大家讨论。

在基于非线性系统的理论基础上，延迟嵌入理论以及广义嵌入理论等相空间重构的理论基础上，观察者便有可能从一个观察到的系统中的时间序列数据中重建一个低维吸引子。通过对N维时间序列进行相空间重构，三维的吸引子可以被重建出来，并且认为这种重建出来后的系统，保留了原始系统的动态特性。随机分布嵌入理论表明，在通过延迟嵌入和非延迟嵌入重构出的相空间之间，存在一定的映射关系，因此可以通过非延迟嵌入重构出的相空间转化为延迟嵌入的相空间，从而达成预测。

随机分布嵌入（RDE）框架是2018年Huanfei Ma在PNAS刊中提出的无模型框架，用于根据观察到的短期高维数据来准确预测未来动态。RDE框架可以被看做是观测到的高维变量之间的空间信息与时间动力学的时间相关概率分布之间的交换方案。在面临原始数据的维度灾难和数据量匮乏的窘境时，随机分布嵌入框架可以使用高维中各维度之间的联系，获得更多的信息量，用以弥补数据量不足所带来的困难，因此，其显著提高了目标变量之间的可预测性。通过对代表性模型和真实系统产生的短时高维数据使用RDE框架，表明即使在噪声扰动之下，高维特征也不再是短期数据准确预测的障碍，而是信息来源。

在算法实现层面，该理论认为，经过很多次迭代后，就会得到很过个预测值，虽然每一个预测值都有可能有误差，但对所有的预测值经过处理，去除掉一些错误非常明显的值之后，对剩余的可能正确值取数学期望，便可以得到最后的预测值。在预测的算法层面，其主要运用了高斯过程回归GPR，该函数在matlab里面可以用fitrgp来调用，详见matlab的技术文档。

由于CSDN代码段不支持matlab的代码块，所以在线阅读起来比较麻烦。

主函数

clear all;
% Data generation
Y=mylorenz();% coupled lorenz system
noisestrength=0; %noise strength
X=Y+noisestrength*rand(size(Y));% noise could be added
xx=X(2000+1:end,:)';% after transient dynamicstrainlength=30; % length of training data (observed data)
traindata=xx(:,1:trainlength);D=size(traindata,1); % number of varialbes in the system. s=2000;% number of non-delay embedding
L=6;% embedding dimension, which could be determined using FNN or set empiricallypredictions=zeros(1,s);
j=1; % the index of target variable
indexr=zeros(s,D);
score=zeros(1,s);
% making predictions with RDE using KDE schema;
for i=1:sindexr(i,:)=randperm(D);predictions(i)=myprediction_gp(traindata(indexr(i,1:L),1:trainlength-1),traindata(j,2:trainlength),traindata(indexr(i,1:L),trainlength));% other kinds of fitting method could be used here
end
pp=outlieromit(predictions);% exclude the outliers
[F,XI]=ksdensity(pp,linspace(min(pp),max(pp),10000));% use kernal density estimation to approximate the probability distribution
prediction=sum(XI.*F/sum(F)); % use expectation as the final one-step predicted value
ystd=std(pp);% % plot the result
% figure
% plot(xx(j,1:trainlength+1),'-*');% real data
% hold on;
% plot(trainlength+1,prediction,'ro','MarkerSize',8); %predicted data
%
% figure
% plot(XI,F); % probablity distribution generated by RDE framework
% hold on
% plot(xx(j,trainlength+1),max(F),'bo','MarkerFaceColor','r'); % true value of the target variable% making predictions with RDE using aggregation schema
for i=1:sindexr(i,:)=randperm(D);score(i)=LOOtest(traindata(indexr(i,1:L),1:trainlength-1),traindata(j,2:trainlength));% leave one out test to get training error, or simply one can use the training residul directlyi
end
[tempp,in]=sort(score,'ascend');
indexrecord=indexr(in(1:L),:);% r best tuples
s=size(indexrecord);
s=s(1);
for i=1:Lpredictions2(i)=myprediction_gp(traindata(indexrecord(i,1:L),1:trainlength-1),traindata(j,2:trainlength),traindata(indexrecord(i,1:L),trainlength));omega(i)=exp(-score(in(i))/score(1));
end
omega=omega./sum(omega);  %归一化
prediction2=omega*predictions2';% use the weighted mean as the final prediction  figure
plot(xx(j,1:trainlength+1),'-*');% real data
hold on;
plot(trainlength+1,prediction2,'ro','MarkerSize',8); %predicted data

去除错误值代码

function score=LOOtest(inputs,outputs)error=zeros(1,length(inputs));
for i=1:length(inputs)inputs1=inputs;inputs1(:,i)=[];outputs1=outputs;outputs1(:,i)=[];prediction=myprediction_gp(inputs1,outputs1,inputs(:,i));error(i)=outputs(i)-prediction;
end
score=mean(abs(error));

预测代码段

% using gaussian process to train the input-output relation and make
% prediction for inputs
% traininputs: n*m, n is dimension, m is length
% trainoutputs: 1*m, m is length
% inputs: n*l, n is dimension, l is the length of inputsfunction outputs = myprediction_gp(traininputs,trainoutputs,inputs)
warning off;
s=size(traininputs);
n=s(1);     %n 训练数据的维数 即嵌入维数
powers=mypower(n,2);
s=size(powers);
l=s(1);
gprMdl = fitrgp(traininputs',trainoutputs','Basis',@mypowerseries,'Beta',zeros(l,1),'FitMethod','exact','PredictMethod','exact');
outputs=predict(gprMdl,inputs');

最后的预测结果如图：

随机分布嵌入（RDE）框架详解附代码相关推荐

动态规划---01背包问题--Dp（详解附代码）
一.动态规划代表一类问题(最优子结构或子问题最优性)的一般解法,是设计方法或者策略,不是具体算法本质:递推,核心是找到状态转移的方式,写出dp方程. 解决问题:交叉,重叠子问题(最优子问题) 形式 ...
【排序】堆排序详解附代码
按照国际惯例,开篇前先简单介绍(吹一波)堆排序(Heapsort).Heapsort是一种优秀的排序算法(个人感觉基本排序算法中仅次于快速排序),时间复杂度为O(nlgn),同时,Heapsort具有 ...
目标检测模型的评估指标mAP详解(附代码）
https://zhuanlan.zhihu.com/p/37910324 对于使用机器学习解决的大多数常见问题,通常有多种可用的模型.每个模型都有自己的独特之处,并随因素变化而表现不同. 每个模型在 ...
c++实现贪吃蛇详解(附代码)
文章目录前言一.运行界面二.类的大致抽象三.关于一些问题的思考四.最后一些想说的五.代码前言经过一个多月的学习,又加深了对c++的理解,所以接下来,就和大家分享一下,一个月学习c++的 ...
各种进制转换（二，八，十，十六进制间转换）详解附代码
进制转换原理进制转换是人们利用符号来计数的方法.进制转换由一组数码符号和两个基本因素"基数"与"位权"构成. 基数是指,进位计数制中所采用的数码(数制中用来 ...
前序遍历、中序遍历、后序遍历层序遍历详解附代码（数据结构C语言）
目录 (1)前序遍历 (DLR) 递归算法 (2)中序遍历 (LDR) 递归算法 (3)后序遍历 (LRD) 递归算法 (4)层序遍历队列实现方法层序遍历的定义: 实现方法: 代码实现结果截图 ...
【大道至简】机器学习算法之EM算法(Expectation Maximization Algorithm)详解(附代码)---通俗理解EM算法。
☕️ 本文来自专栏:大道至简之机器学习系列专栏
Numpy学习笔记(二):argmax参数中axis=0,axis=1,axis=-1详解附代码
文章目录 1.argmax和max函数区别 2.axis=0/axis=1/axis=-1的区别 3.具体代码分析 ---3.1一维数组 ---3.2二维数组 ---3.3三维数组 1.argmax和 ...
MATLAB生成伯努利图分布,C++ - 随机生成器伯努利分布(bernoulli distribution) 的详解及代码...
随机生成器伯努利分布(bernoulli distribution) 的详解及代码本文地址: http://blog.csdn.net/caroline_wendy/article/deta ...

随机分布嵌入（RDE）框架详解附代码

介绍

主函数

去除错误值代码

预测代码段

随机分布嵌入（RDE）框架详解附代码相关推荐

最新文章

热门文章