Загрузка и проверка данных

Загрузите набор данных, содержащий ставки дефолта кредита.

load Data_CreditDefaults

Для получения дополнительной информации о наборе данных отобразите Description.

Таблица DataTable содержит ежегодные измерения, но набор данных является агностическим для временной основы. Чтобы применить временную основу к данным, преобразуйте DataTable в расписание.

Tbl = table2timetable(DataTable,"RowTimes",datetime(DataTable.Dates,"Format","yyyy"));

Econometrics Toolbox™ arima программа удаляет отсутствующие наблюдения (обозначенные NaNs) из данных при помощи спискового удаления. Это действие влияет на то, как вы разделяете данные.

Определите, содержит ли какая-либо строка данных хотя бы одно отсутствующее наблюдение.

idxmiss = any(ismissing(Tbl),2);
anymissing = sum(idxmiss)

anymissing = 0

Отсутствующих наблюдений в данных нет.

Определите общий размер выборки.

T = size(Tbl,1)

T = 21

Временная шкала разбиения

Предположим, что все последовательности в данных являются стационарными, и рассмотрим эту модель ARX (1) для скорости IGD:

где ${\epsilon }_{\mathit{t}}$ - iid-ряд Гауссовых случайных переменных со средним значением 0 и конечным отклонением.

Чтобы оценить модель ARX (1), условное среднее значения первого отклика требует предыдущего отклика. Поэтому период предварительного образца для оценки является первой строкой данных. Также рассмотрите прогнозирование IGD на 8-летний прогнозный горизонт и сравнение прогнозов с наблюдаемыми ответами. Поэтому прогнозный период является последними восемью строками данных. Подбор модели к остальным данным.

Создайте векторы индекса, которые идентифицируют предварительную выборку, оценку и прогнозные периоды.

fh = 8;                          % Forecast horizon
idxpre = 1;                      
idxest = 2:(T - fh);   
idxf = (T - fh + 1):T; 

fprintf('\nPresample period: %d\nEstimation Period: %d - %d\nForecast period: %d - %d\n',...
    year(Tbl.Time([idxpre idxest(1) idxest(end) idxf(1) idxf(end)])))

Presample period: 1984
Estimation Period: 1985 - 1996
Forecast period: 1997 - 2004

Подгонка модели к данным

Создайте arima объект модели, который представляет модель ARX (1).

Mdl = arima(1,0,0)

Mdl = 
  arima with properties:

     Description: "ARIMA(1,0,0) Model (Gaussian Distribution)"
    Distribution: Name = "Gaussian"
               P: 1
               D: 0
               Q: 0
        Constant: NaN
              AR: {NaN} at lag [1]
             SAR: {}
              MA: {}
             SMA: {}
     Seasonality: 0
            Beta: [1×0]
        Variance: NaN

Mdl является шаблоном для оценки. Свойства Mdl которые имеют значение NaN соответствуют неизвестным параметрам модели ARX (1), которые должны быть оценены estimate. Хотя свойство коэффициента регрессии Beta пуст, estimate определяет его размер по количеству столбцов в указанной матрице экзогенных данных.

Подгонка модели к реакции выборки оценки и экзогенным данным. Задайте предварительный образец отклика.

X = Tbl{:,["AGE" "CPF" "SPR"]};
Mdl = estimate(Mdl,Tbl.IGD(idxest),'Y0',Tbl.IGD(idxpre),...
    'X',X(idxest,:));

 
    ARIMAX(1,0,0) Model (Gaussian Distribution):
 
                  Value       StandardError    TStatistic     PValue 
                __________    _____________    __________    ________

    Constant      0.022517        0.20278        0.11104      0.91158
    AR{1}         -0.18486        0.43861       -0.42146      0.67342
    Beta(1)      0.0091356       0.016249        0.56224      0.57395
    Beta(2)     -0.0098351      0.0053577        -1.8357     0.066402
    Beta(3)      0.0039044       0.057031       0.068461      0.94542
    Variance     0.0027082      0.0014648         1.8489     0.064468

Mdl является полностью заданным arima объект модели, представляющий подобранную модель ARX (1). Сводные данные оценок предполагают, что все переменные, кроме CPF (Beta(2)) незначительны на 10% уровне значимости.

Сгенерируйте прогнозы MMSE

Прогнозируйте отклики MMSE на 8-летний прогнозный горизонт и вычисляйте соответствующие прогнозные отклонения ошибки. Потому что прогноз на 1 период вперед $${\underset{}{\overset{ˆ}{y}}}_{1997}$$ требуется предварительный пример отклика ${\mathit{y}}_{1996}$ и $\mathit{t}$-прогнозирование с опережением требует экзогенных данных ${\mathit{X}}_{\mathit{t}}$, задайте конечный ответ в выборке оценки как предварительный пример отклика и укажите экзогенные данные прогнозного периода.

[yf,ymse] = forecast(Mdl,fh,Tbl.IGD(idxest(end)),...
    'XF',X(idxf,:));

Вычислите интервалы предсказания для прогнозов, принимая, что они обычно распределены.

ci = yf + 1.96*[-sqrt(ymse) sqrt(ymse)];

Построение прогнозов MMSE

Постройте график серии ответов с помощью их прогнозов MMSE и интервалов предсказания.

figure
h1 = plot(Tbl.Time,Tbl.IGD,'ko-','LineWidth',2);
hold on
h2 = plot(Tbl.Time(idxf),yf,'k:','LineWidth',2);
h3 = plot(Tbl.Time(idxf),ci,'r:','LineWidth',2);
h = gca;
px = Tbl.Time(idxf([1 end end 1]));
py = h.YLim([1 1 2 2]);
hf = fill(px,py,[0.9 0.9 0.9]);
uistack(hf,"bottom");
legend([h1 h2 h3(1) hf],["Observation" "Forecast"...
    "95% Interval" "Forecast period"],'Location','NorthWest')
title('IGD Rates and MMSE Forecasts')
xlabel('Year')
ylabel('IGD (%)')
axis tight
hold off

Прогнозы кажутся разумными, несмотря на крайние наблюдения, имевшие место в 2001 году.

Сгенерируйте прогнозы Монте-Карло

Прогнозные отклики путем симуляции 100 путей из подобранной модели ARX (1). Затем агрегируйте случайные пути путем вычисления среднего и 2,5% и 97,5% процентильных интервалов в каждой временной точке в прогнозном периоде. Аналогично оценке прогноза MMSE, simulate требуется предварительный пример отклика ${\mathit{y}}_{1996}$ произвести ${\mathit{y}}_{1997}$ из подобранной модели ARX (1), и она требует экзогенных данных в прогнозном периоде.

nsim = 100;
rng(1);
YMC = simulate(Mdl,fh,'NumPaths',nsim,...
    'Y0',Tbl.IGD(idxest(end)),'X',X(idxf,:));
YMCbar = mean(YMC,2);
YMCci = quantile(YMC',[0.025 0.975]);

Сюжет Прогнозы Монте-Карло

Постройте график серии ответов с прогнозами Монте-Карло и 95% интервалами процентиля.

figure
h1 = plot(Tbl.Time,Tbl.IGD,'ko-','LineWidth',2);
hold on
h2 = plot(Tbl.Time(idxf),YMC,'Color',[0.7 0.7 0.7]);
h3 = plot(Tbl.Time(idxf),YMCbar,'c:','LineWidth',3);
h4 = plot(Tbl.Time(idxf),YMCci,'r:','LineWidth',3);
h = gca;
px = Tbl.Time(idxf([1 end end 1]));
py = h.YLim([1 1 2 2]);
hf = fill(px,py,[0.9 0.9 0.9]);
uistack(hf,"bottom");
legend([h1 h2(1) h3 h4(1) hf],["Observation" "Monte Carlo path"...
    "Monte Carlo mean" "95% Interval" "Forecast period"],'Location','NorthWest')
title('IGD Rates and Monte Carlo Forecasts')
xlabel('Year')
ylabel('IGD (%)')
axis tight
hold off

Прогнозы Монте-Карло и интервалы предсказания напоминают прогнозы MMSE.