В этом примере показано, как предсказать наблюдения за известным, независимым от времени, моделью в пространстве состояний.
Предположим, что скрытый процесс является AR (1). Уравнение состояния
где является Гауссовым со средним значением 0 и стандартным отклонением 1.
Сгенерируйте случайную последовательность 100 наблюдений от , предположение, что ряд запускается в 1,5.
T = 100; ARMdl = arima('AR',0.5,'Constant',0,'Variance',1); x0 = 1.5; rng(1); % For reproducibility x = simulate(ARMdl,T,'Y0',x0);
Предположим далее, что скрытый процесс подвергается аддитивной погрешности измерения. Уравнение наблюдения
где является Гауссовым со средним значением 0 и стандартным отклонением 0.75. Вместе, скрытые уравнения процесса и наблюдения составляют модель в пространстве состояний.
Используйте случайный скрытый процесс состояния (x
) и уравнение наблюдения, чтобы сгенерировать наблюдения.
y = x + 0.75*randn(T,1);
Задайте четыре содействующих матрицы.
A = 0.5; B = 1; C = 1; D = 0.75;
Задайте модель в пространстве состояний с помощью содействующих матриц.
Mdl = ssm(A,B,C,D)
Mdl = State-space model type: ssm State vector length: 1 Observation vector length: 1 State disturbance vector length: 1 Observation innovation vector length: 1 Sample size supported by model: Unlimited State variables: x1, x2,... State disturbances: u1, u2,... Observation series: y1, y2,... Observation innovations: e1, e2,... State equation: x1(t) = (0.50)x1(t-1) + u1(t) Observation equation: y1(t) = x1(t) + (0.75)e1(t) Initial state distribution: Initial state means x1 0 Initial state covariance matrix x1 x1 1.33 State types x1 Stationary
Mdl
ssm
модель. Проверьте, что модель правильно задана с помощью отображения в Командном окне. Программное обеспечение выводит, что процесс состояния является стационарным. Впоследствии, программное обеспечение устанавливает среднее значение начального состояния и ковариацию к среднему значению и отклонению стационарного распределения модели AR (1).
Предскажите наблюдения 10 периодов в будущее и оцените их отклонения.
numPeriods = 10; [ForecastedY,YMSE] = forecast(Mdl,numPeriods,y);
Постройте прогнозы с ответами в выборке, и 95% интервалов прогноза вальдового типа.
ForecastIntervals(:,1) = ForecastedY - 1.96*sqrt(YMSE); ForecastIntervals(:,2) = ForecastedY + 1.96*sqrt(YMSE); figure plot(T-20:T,y(T-20:T),'-k',T+1:T+numPeriods,ForecastedY,'-.r',... T+1:T+numPeriods,ForecastIntervals,'-.b',... T:T+1,[y(end)*ones(3,1),[ForecastedY(1);ForecastIntervals(1,:)']],':k',... 'LineWidth',2) hold on title({'Observed Responses and Their Forecasts'}) xlabel('Period') ylabel('Responses') legend({'Observations','Forecasted observations','95% forecast intervals'},... 'Location','Best') hold off