【MATLAB第117期】#源碼分享 | 基于MATLAB的SSM狀態空間模型多元時間序列預測方法（多輸入單輸出）
引言
本文使用狀態空間模型實現失業率遞歸預測，狀態空間模型（State Space Model, SSM）是一種用于描述動態系統行為的數學模型，通過狀態變量、輸入和輸出的關系來刻畫系統的時變特性。
本示例演示如何使用Econometrics Toolbox中的狀態空間模型，實現美國年度失業率的滾動窗口預測
% 模型結構：
% x(t) = A*x(t-1) + B*u(t) + w(t) (狀態方程)
% y(t) = C*x(t) + D*u(t) + v(t) (觀測方程)
% 其中：
% y(t) - 觀測失業率
% x(t) - 潛在狀態變量
% u(t) - 名義GNP增長率（外生變量）
% w(t) ~ N(0,Q) - 過程噪聲
% v(t) ~ N(0,R) - 觀測噪聲
一、關鍵思路
?1、滾動窗口機制?
使用31年的窗口進行參數估計，逐年滾動預測下一年失業率
窗口滑動方式：1940-1970年間共生成31個預測窗口
?2、狀態空間模型結構?
x(t)=[α1;01]?x(t?1)+[β;0]?ΔlnGNP(t)+w(t)y(t)=[10]?x(t)+γ?ΔlnGNP(t)+v(t)
包含AR(2)動態和GNP的外生影響
參數估計目標：α(自回歸系數), β(GNP對狀態的系數), γ(GNP對觀測的直接系數)
3、?數據預處理細節?
對原始GNP數據取對數差分處理：ΔlnGNP = diff(log(GNP))
失業率使用一階差分：ΔUR = diff(UR)
二、核心代碼
%% 1. 導入年度經濟數據
load Data_NelsonPlosser; % 加載Nelson-Plosser宏觀經濟數據集
% 數據集包含美國1909-1970年間的經濟指標，包括：
% GNPN - 名義國民生產總值
% UR - 失業率
% 其他指標：CPI、工資率等%% 2. 數據預處理
isNaN = any(ismissing(DataTable),2); % 標記包含缺失值的行
Z = DataTable.GNPN(~isNaN); % 提取有效期的名義GNP數據（61×1）
y = DataTable.UR(~isNaN); % 提取有效期的失業率數據（61×1）%% 3. 創建時間序列數組
WindowSize = 31; % 滾動窗口大小（31年）
ForecastPeriod = numel(y) - WindowSize + 1; % 預測期數（61-31+1=31）% 初始化存儲矩陣
ZZ = zeros(ForecastPeriod, WindowSize); % GNP窗口數據（31×31）
yy = zeros(ForecastPeriod, WindowSize); % 失業率窗口數據（31×31）% 創建滾動窗口數據集
m = 1;
for nYear = 1:ForecastPeriodZZ(nYear,:) = Z(m:m+WindowSize-1); % 當前窗口的GNP數據yy(nYear,:) = y(m:m+WindowSize-1); % 當前窗口的失業率數據m = m + 1;
end% 提取時間戳（1940-1970年）
Time = str2double(DataTable.Properties.RowNames(~isNaN));
Time = Time(end-ForecastPeriod+1:end); % 構建帶時間戳的數組
ObsUnemployRate = [Time, yy]; % 失業率時間序列（31×32）
nGNP = [Time, ZZ]; % GNP時間序列（31×32）%% 4. MATLAB狀態空間模型遞歸估計
eUR = zeros(numel(Time),1); % 存儲預測失業率
param0 = [0.5; 0.1; -20]; % 初始參數 [A; B; D]% 滾動窗口參數估計循環
for t = 0:numel(Time)-1% 數據準備dlZ = diff(log(ZZ(t+1,:)))'; % GNP對數差分（增長率）dy = diff(yy(t+1,:))'; % 失業率差分% 狀態空間模型定義（調用rwAR2ParamMap函數）Mdl = ssm(@(c)rwAR2ParamMap(c,dy,dlZ));% 參數估計（最大似然估計）[Mdl, param0] = estimate(Mdl, dy, param0, 'Display', 'off');% 1步超前預測dyhat = forecast(Mdl, 1, dy,...'Predictors0', dlZ,...'PredictorsF', dlZ(end),...'Beta', param0(end));% 整合預測結果eUR(t+1) = dyhat + yy(t+1,end); % 預測值 = 差分預測 + 當前水平值
end%% 5. 可視化MATLAB預測結果
figure;
axH = axes;
plot(axH, Time, y(end-numel(Time)+1:end),... % 實際失業率'Color', [0.9290 0.6940 0.1250], 'LineWidth', 1.2);
hold on; grid on;
plot(axH, Time, eUR,... % 預測失業率'Color', [0 0.4470 0.7410], 'LineWidth', 1.2);
axH.XLim(1) = Time(1);
axH.Color = [0.5020 0.5020 0.5020];
title('失業率預測對比 (%)');
legend(["實際值", "預測值"]);
三、代碼獲取
1.閱讀首頁置頂文章
2.關注CSDN
3.根據自動回復消息，回復“117期”以及相應指令，即可獲取對應下載方式。
