Рассмотрим модель множественной линейной регрессии, которая предсказывает реальный валовой национальный продукт США (GNPR) с использованием линейной комбинации индекса промышленного производства (IPI), общая занятость (E) и реальная заработная плата (WR).

Для всех $$t$$$${}_{\mathrm{\delta \; t}}$$ - это ряд независимых гауссовых возмущений со средним значением 0 и дисперсией $${}^{\mathrm{start2}}$$.

Предположим, что эти предыдущие распределения:

Эти допущения и правдоподобие данных подразумевают нормально-обратно-гамма-сопряженную модель.

Создайте нормально-обратно-гамма-сопряженную предыдущую модель для параметров линейной регрессии. Укажите количество предикторов p и имена переменных.

p = 3;
VarNames = ["IPI" "E" "WR"];
PriorMdl = bayeslm(p,'ModelType','conjugate','VarNames',VarNames);

PriorMdl является conjugateblm Объект байесовской модели линейной регрессии, представляющий предварительное распределение коэффициентов регрессии и дисперсии возмущений.

Постройте график предыдущих распределений.

plot(PriorMdl);

plot строит графики предельных предварительных распределений перехвата, коэффициентов регрессии и дисперсии возмущений.

Предположим, что среднее значение коэффициентов регрессии равно $$[-\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}2040.0012{]}^{}$$′ и их масштабированная ковариационная матрица равна

Кроме того, предшествующая шкала дисперсии возмущений составляет 0,01. Укажите предшествующую информацию с помощью точечной нотации.

PriorMdl.Mu = [-20; 4; 0.001; 2];
PriorMdl.V = diag([1 0.001 1e-8 0.01]);
PriorMdl.B = 0.01;

Запросите новый рисунок и постройте график предыдущего распределения.

plot(PriorMdl);

plot заменяет текущий показатель распределения графиком предыдущего распределения дисперсии возмущений.

Загрузите набор данных Нельсона-Плоссера и создайте переменные для данных предиктора и ответа.

load Data_NelsonPlosser
X = DataTable{:,PriorMdl.VarNames(2:end)};
y = DataTable.GNPR;

Оцените апостериорные распределения.

PosteriorMdl = estimate(PriorMdl,X,y,'Display',false);

PosteriorMdl является conjugateblm объект модели, который содержит апостериорные распределения $$\beta$$ и $${}^{\mathrm{start2}}$$.

Постройте график задних распределений.

plot(PosteriorMdl);

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

plot(PosteriorMdl,PriorMdl);

Рассмотрим регрессионную модель в сюжетных предварительных и задних распределениях.

Загрузите набор данных Нельсона-Плоссера, создайте предыдущую сопряженную модель по умолчанию, а затем оцените заднюю модель, используя первые 75% данных. Выключите отображение оценки.

p = 3;
VarNames = ["IPI" "E" "WR"];
PriorMdl = bayeslm(p,'ModelType','conjugate','VarNames',VarNames);

load Data_NelsonPlosser
X = DataTable{:,PriorMdl.VarNames(2:end)};
y = DataTable.GNPR;

d = 0.75;
PosteriorMdlFirst = estimate(PriorMdl,X(1:floor(d*end),:),y(1:floor(d*end)),...
    'Display',false);

Постройте график предыдущего распределения и заднего распределения дисперсии возмущения. Верните дескриптор фигуры.

[~,~,~,h] = plot(PosteriorMdlFirst,PriorMdl,'VarNames','Sigma2');

h - маркер перемещения фигуры для графика распределения. Если изменить имя марки фигуры, изменив Tag свойство, затем следующее plot вызов помещает все новые участки распределения на другой показатель.

Измените имя маркера перемещения фигуры на FirstHalfData с использованием точечной нотации.

h.Tag = 'FirstHalfData';

Оцените апостериорное распределение, используя остальные данные. Укажите апостериорное распределение на основе последних 25% данных в качестве предыдущего распределения.

PosteriorMdl = estimate(PosteriorMdlFirst,X(ceil(d*end):end,:),...
    y(ceil(d*end):end),'Display',false);

Постройте график задней части дисперсии возмущений на основе половины данных и всех данных на новый рисунок.

plot(PosteriorMdl,PosteriorMdlFirst,'VarNames','Sigma2');

Этот тип графика показывает эволюцию заднего распределения при включении новых данных.

Рассмотрим регрессионную модель в сюжетных предварительных и задних распределениях.

Загрузите набор данных Нельсона-Плоссера и создайте предыдущую сопряженную модель по умолчанию.

p = 3;
VarNames = ["IPI" "E" "WR"];
PriorMdl = bayeslm(p,'ModelType','conjugate','VarNames',VarNames);

load Data_NelsonPlosser
X = DataTable{:,PriorMdl.VarNames(2:end)};
y = DataTable.GNPR;

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

[pointsUsedPrior,priorDensities1] = plot(PriorMdl);

pointsUsedPrior является массивом ячеек 5 на 1 из 1 на 1000 числовых векторов, представляющих значения параметров, которые plot используется для построения графика соответствующих плотностей. Первый элемент соответствует перехвату, следующие три элемента соответствуют коэффициентам регрессии, а последний элемент соответствует дисперсии возмущений. priorDensities1 имеет те же размеры, что и pointsUsed и содержит соответствующие значения плотности.

Оцените апостериорное распределение. Выключите отображение оценки.

PosteriorMdl = estimate(PriorMdl,X,y,'Display',false);

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

[pointsUsedPost,posteriorDensities1] = plot(PosteriorMdl);

pointsUsedPost и posteriorDensities1 имеют те же размеры, что и pointsUsedPrior. pointsUsedPost выходные данные могут отличаться от pointsUsedPrior. posteriorDensities1 содержит значения задней плотности.

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

[pointsUsedPP,posteriorDensities2,priorDensities2] = plot(PosteriorMdl,PriorMdl);

Все выходные значения имеют те же размеры, что и pointsUsedPrior. posteriorDensities2 выход содержит задние значения плотности. priorDensities2 выходные данные содержат предыдущие значения плотности.

Подтвердить, что pointsUsedPP равно pointsUsedPost.

compare = @(a,b)sum(a == b) == numel(a);
cellfun(compare,pointsUsedPost,pointsUsedPP)

ans = 5x1 logical array

   1
   1
   1
   1
   1

Используемые точки эквивалентны.

Подтвердите, что задние плотности одинаковы, но что предыдущие плотности не одинаковы.

cellfun(compare,posteriorDensities1,posteriorDensities2)

ans = 5x1 logical array

   1
   1
   1
   1
   1

cellfun(compare,priorDensities1,priorDensities2)

ans = 5x1 logical array

   0
   0
   0
   0
   0

При печати только предыдущего распределения, plot оценивает предыдущие плотности в точках, которые создают четкий график предыдущего распределения. При построении графика как предшествующего, так и заднего распределения, plot предпочитает рисовать заднюю чётко. Поэтому график может определять другой набор точек для использования.

Рассмотрим регрессионную модель в сюжетных предварительных и задних распределениях.

Загрузите набор данных Нельсона-Плоссера и создайте по умолчанию сопряженную предыдущую модель для коэффициентов регрессии и дисперсии возмущений. Затем оцените апостериорное распределение и получите сводную таблицу оценки из summarize.

p = 3;
VarNames = ["IPI" "E" "WR"];
PriorMdl = bayeslm(p,'ModelType','conjugate','VarNames',VarNames);

load Data_NelsonPlosser
X = DataTable{:,PriorMdl.VarNames(2:end)};
y = DataTable.GNPR;

PosteriorMdl = estimate(PriorMdl,X,y);

Method: Analytic posterior distributions
Number of observations: 62
Number of predictors:   4
Log marginal likelihood: -259.348
 
           |   Mean      Std          CI95        Positive       Distribution      
-----------------------------------------------------------------------------------
 Intercept | -24.2494  8.7821  [-41.514, -6.985]    0.003   t (-24.25, 8.65^2, 68) 
 IPI       |   4.3913  0.1414   [ 4.113,  4.669]    1.000   t (4.39, 0.14^2, 68)   
 E         |   0.0011  0.0003   [ 0.000,  0.002]    1.000   t (0.00, 0.00^2, 68)   
 WR        |   2.4683  0.3490   [ 1.782,  3.154]    1.000   t (2.47, 0.34^2, 68)   
 Sigma2    |  44.1347  7.8020   [31.427, 61.855]    1.000   IG(34.00, 0.00069)     
 

summaryTbl = summarize(PosteriorMdl);
summaryTbl = summaryTbl.MarginalDistributions;

summaryTbl - таблица, содержащая статистику, которая estimate отображается в командной строке.

Для каждого параметра определите набор из 50 равномерно разнесенных значений в пределах трех стандартных отклонений среднего значения. Поместите значения в ячейки вектора ячейки 5 на 1, следуя порядку параметров, которые составляют строки сводной таблицы оценки.

Points = cell(numel(summaryTbl.Mean),1); % Preallocation

for j = 1:numel(summaryTbl.Mean)
    Points{j} = linspace(summaryTbl.Mean(j) - 3*summaryTbl.Std(j),...
        summaryTbl.Mean(j) + 2*summaryTbl.Std(j),50);
end

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

plot(PosteriorMdl,'Points',Points)