时间序列预测和预测预后

打开脚本

这个例子显示了如何创建一个时间序列模型，并使用模型进行预测，预测和状态估计。测得的数据是从一个感应炉，其时隙大小溶蚀随时间。所述时隙大小不能直接测量，但在炉的电流和消耗的功率进行测量。已知的是，作为所述时隙大小的增加，槽电阻降低。平方以测得的功率测量的电流的比率是正比于所述槽的大小。可以使用测量的电流 - 功率比（电流和功率测量是有噪声的）来创建一个时间序列模型，并使用模型来估计当前时隙大小和预测未来的时隙大小。通过物理检查感应炉时隙大小在时间上在某些点是已知的。

加载和绘制测量数据

所测量的电流 - 功率比数据被存储在iddata_TimeSeriesPredictionMATLAB文件。数据在每小时间隔和显示测量为表示炉槽的侵蚀随时间的比率增加。您开发使用该数据的时间序列模型。由数据分离成识别和验证段开始。

加载iddata_TimeSeriesPredictionN = numel（Y）;NS =地板（N / 2）;y_id = Y（1：NS，:);y_v = Y（（NS + 1：结束），:);data_id = IDDATA（y_id，[]中，Ts，'TIMEUNIT'，'小时'）;data_v = IDDATA（y_v，[]中，Ts，'TIMEUNIT'，'小时'，“T开始”，NS + 1）;图（data_id，data_v）传说（“标识数据”，“验证数据”，'位置'，'东南'）;

型号识别

槽侵蚀可以被建模为具有噪声输入和测量的电流 - 功率比作为输出的状态空间系统。所测量的电流 - 功率比正比于系统的状态，或

$X_ {N + 1} = Ax_n + Ke_n$

$ y_n = Cx_n + e_n $

哪里 $ x_n $ 状态矢量，包含时隙大小; $ y_n $ 是所测量的电流 - 功率比; $ e_n $ 噪声和 $ A，C，K $ 将被识别。

使用ssest（）命令到从测量数据识别离散的状态空间模型。

SYS = ssest（data_id，1，'TS'，TS，'形成'，'典范'）

SYS =离散时间识别出的状态空间模型：X（t + TS）= A X（t）+ K E（t）的Y（吨）= C X（t）+ E（t）的A = X1 X1 1.001Ç= X1 Y1 1 K = Y1 X1 0.09465采样时间：1小时参数化：规范形式具有索引：1.干扰成分：自由系数的估计数：2使用 “idssdata”， “getpvec”， “getcov” 为参数和它们的不确定度。状态：使用SSEST时域数据“data_id”估算。适合于估计数据：67.38％（预测焦点）FPE：0.09575，MSE：0.09348

所识别的模型最小化了1步提前预测。使用10领先一步预测器验证该模型，即，给定 $ y_0，...，y_n $ 使用模型预测 $Y_ {N + 10}$ 。需要注意的是所测量的值和预测值之间的误差， $y_0- \ {帽子} y_0，...，y_n- \ {帽子} y_n$ ，是用来做 $Y_ {N + 10}$ 预测。

使用10提前一步预测的识别数据和独立的验证数据。

n步= 10;比较（SYS，data_id，n步）10步预测的％比较估计数据网格（'上'）;

数字;比较（SYS，data_v，n步）％比较验证数据网格（'上'）;

上述两种运动数据组显示，预测器所测量的数据相匹配。

预测用于进一步验证模型。预测使用测量的数据记录 $y_0，Y_1，...，y_n- \ {帽子} y_n$ 以计算在时间步骤n中的模型状态。该值被用作初始条件为预测为将来的时间跨度中的模型响应。我们预计在验证数据的时间跨度模型响应，然后比较这两个。我们也可以计算在预报的不确定性和绘制+/- 3 SD它们的值。

MeasuredData = IDDATA（Y，[]中，Ts，'TIMEUNIT'，'小时'）;（％）= [data_id; data_v]T0 = MeasuredData.SamplingInstants;地平线=尺寸（data_v，1）;％预测地平线[YF，〜，〜，yFSD] =预测（SYS，data_id地平线）;％注：YF是IDDATA对象而yFSD是双载体T = yF.SamplingInstants;％提取时间样本yFData = yF.OutputData;％提取响应如双载体图（MeasuredData）保持上图（T，yFData，'r.-'，T，yFData + 3 * yFSD，'R--'，T，yFData-3 * yFSD，'R--'）保持离标题（的时间跨度“在验证数据预测的响应”）网格上

该图显示，与置信区间模型响应（由红色着色虚曲线表示）重叠的验证数据的测量值。将合并的预测和预测结果表明，该模型代表测得的电流 - 功率比。

预测结果还表明，在大视野模型的方差较大，为了实用的目的未来的预测应限于短期的视野。用于感应炉模型200小时地平线是适当的。

最后，我们用模型来预测响应为502-701小时的时间跨度200步到未来。

地平线= 200;％预测地平线[yFuture，〜，〜，yFutureSD] =预测（SYS，MeasuredData地平线）;T = yFuture.SamplingInstants;％提取时间样本yFutureData = yFuture.OutputData;％提取响应如双载体图（T0，Y，...T，yFutureData，'r.-'，...吨，yFutureData + 3 * yFutureSD，'R--'，...吨，yFutureData-3 * yFutureSD，'R--'）标题（“预测响应（200步）”）网格上

蓝色曲线表示所测量的数据跨越1-501小时。红色曲线为超出测量数据的时间范围200小时预测响应。红色虚曲线示出了基于所识别的模型的随机抽样在预测响应中的3 SD的不确定性。

状态估计

所识别的模型中的测量的电流 - 功率比相匹配，但我们感兴趣的是在炉槽尺寸，该尺寸在模型中的状态。所识别的模型具有可转化的，使得状态已经意，在这种情况下，槽的大小的任意的状态。

创建任意状态的预测。辨识模型协方差需要被翻译使用的预测模型translatecov（）命令。该createPredictor（）函数简单地提取的第三输出参数预测（）与其一起使用的功能translatecov（）。

类型createPredictorEST = translatecov（@（S）createPredictor（S，data_id），SYS）

预解码功能= createPredictor（MDL，数据）％CREATEPREDICTOR返回1-领先一步预测器。％％SYS = createPredictor（MDL，数据）％％创建指定的模型MDL％和测量的数据的1领先一步预测器模型SYS。该功能使用％| TimeSeriedPredictionExample |和| translatecov（）|命令％转化辨识模型协方差预测。％版权所有2015 MathWorks公司[〜，〜，预解码值] =预测（MDL，数据，1）;EST =离散时间识别出的状态空间模型：X（t + TS）= A X（t）+ B u（t）的Y（吨）= C X（t）+ d U（t）的A = X1 X1 0.9064B = Y1 X1 0.09465 C = X1 Y1 1 d = Y1 Y1 0采样时间：1小时参数化：规范形式具有索引：1.穿心：无干扰成分：自由系数的无数：2使用 “idssdata”， “getpvec”“getcov”为参数及其不确定性。状态：创建者直接建设或改造。未估计。

该模型美东时间是1领先一步预测器在相同的状态坐标原始模型表达SYS。我们如何改变状态坐标，使模型的状态对应于（与时间有关）插槽的大小？该解决方案是依靠间歇截取的槽尺寸的实际的，直接的测量。这是在实践中没有地方采取直接测量成本高，只有定期进行（如组件被替换时）少见。

具体而言，变换预测器状态， $ x_n $ ，至 $ z_n $ ，以便 $ y_n = Cz_n $ 哪里 $ y_n $ 所测量的电流 - 功率比，和 $ z_n $ 是在炉槽的大小。在这个例子中，炉槽尺寸的四个直接测量，sizeMeasured和炉当前功率比，ySizeMeasured，用于估计 $ C $ 。在转型中预测预测协方差也需要进行改造。因此，我们使用translatecov（）命令执行国家坐标变换。

Cnew = sizeMeasured \ ySizeMeasured;EST = translatecov（@（S）ss2ss（S，S.C / Cnew），EST）

EST =离散时间识别出的状态空间模型：X（t + TS）= A X（t）+ B u（t）的Y（吨）= C X（t）+ d U（t）的A = X1 X1 0.9064B = Y1 X1 0.9452 C = X1 Y1 0.1001 d = Y1 Y1 0采样时间：1小时参数化：规范形式具有索引：1.穿心：无干扰成分：自由系数的无数：2使用 “idssdata”， “getpvec”“getcov”为参数及其不确定性。状态：创建者直接建设或改造。未估计。

该预测现在是在理想状态坐标来表示。它有一个输入端，是测得的系统输出（炉电流 - 功率比）和一个输出，该输出的预测系统输出（炉时隙大小）。预测器是模拟来估计系统输出和系统状态。

OPTS = simOptions;opts.InitialCondition = sizeMeasured（1）;U = IDDATA（[]，[data_id.Y; data_v.Y]，TS，'TIMEUNIT'，'小时'）;[你们，ye_sd，XE] = SIM（EST，U，OPTS）;

比较测量和已知值估计的输出和位尺寸。

yesdp =你们;yesdp.Y = ye.Y + 3 * ye_sd;yesdn =你们;yesdn.Y = ye.Y-3 * ye_sd;N = numel（XE）;图中，图（[data_id; data_v]，咋，yesdp，'G'，yesdn，'G'）图例（“测量输出”，“预计产出”，'99 0.7％结合的”，'位置'，'东南'）网格（'上'）图，图（tSizeMeasured，sizeMeasured，'R *'，1：N，XE，1：N，yesdp.Y / est.C，'G'，1：N，yesdn.Y / est.C，'G'）;传说（“测量的状态”，“估计状态”，'99 0.7％结合的”，'位置'，'东南'）xlabel（“时间（小时）”）ylabel（'振幅'）;网格（'上'）