How do I plot 14th and 86th percentile?
Here I am currently running Bayesian time varying model in matlab and this is my whole code. I am not sure how can I create 14th and 86th percentile of my posterior distribution/credible interval that I have found using gibbs sampling and carter kohn algorithm. Furthermore, how do I plot the 14th and 86th percentile along with my mean posterior distribution in one graph?
clear;
% Load data from Excel file
Cyclicalitydata = readtable("Could be final data set.xlsx");
Countrydata = Cyclicalitydata(strcmp(Cyclicalitydata.Country, ‘Netherlands’), :);
% Prepare variables
Fiscalvariable = Countrydata.CAPB1; % Dependent variable, original series
LaggedFiscalvariable = [NaN; Fiscalvariable(1:end-1)]; % Creating one-period lag
Economicactivity = Countrydata.Outputgap; % Independent main variable
Govsize = Countrydata.Govsize; % Control variable 1
Debtp = Countrydata.DebtGDP; % Control variable 2
CPIlevel = Countrydata.CPIlevel; % Additional control variable 3
Years = Countrydata.Year; % Year for plotting and analysis
CPIgrowth = [NaN; CPIlevel(1:end-1)];
% Define the state-space model components
T = numel(LaggedFiscalvariable); % Number of observations, using the lagged variable now
C = [Economicactivity, Govsize, Debtp, CPIgrowth]; % Each row of C corresponds to an observation time
N = size(C,2);
A = eye(N); % Assuming independent evolution of each coefficient
Q = 0.01 * eye(N); % Small variance to assume minor changes over time
R = 1; % Observation noise variance
% Calculate OLS estimates for initialization
validIdx = ~isnan(LaggedFiscalvariable) & ~any(isnan(C), 2);
C_valid = C(validIdx, :);
Y_valid = LaggedFiscalvariable(validIdx);
beta_OLS = (C_valid’ * C_valid) (C_valid’ * Y_valid);
residuals = Y_valid – C_valid * beta_OLS;
sigma2_OLS = residuals’ * residuals / (length(Y_valid) – size(C_valid, 2));
% Initialize beta_0 and V_0 with OLS estimates
beta_0 = beta_OLS; % Use OLS estimates as initial states
V_0 = sigma2_OLS * eye(N); % Scale identity matrix by OLS residual variance
% Number of Gibbs sampling iterations
numIterations = 1000;
beta_samples = zeros(N, T, numIterations); % Store samples of all beta coefficients
% Gibbs Sampling Loop
for iter = 1:numIterations
% Run the Carter-Kohn algorithm for the current iteration
[beta_t, ~] = carterKohn(LaggedFiscalvariable, A, C, Q, R, beta_0, V_0, T);
% Store the sample
beta_samples(:, :, iter) = beta_t;
end
% Post-processing after Gibbs sampling
% Focus on the OutputGap coefficient’s estimates and its posterior distribution
burn_in = 100;
posterior_beta_outputgap = squeeze(mean(beta_samples(1, :, burn_in:end), 3));
% Plot only the OutputGap coefficient estimates over years
figure;
plot(Years, posterior_beta_outputgap, ‘LineWidth’, 2);
title(‘Posterior Estimates of OutputGap Coefficient over Years for the Netherlands’);
xlabel(‘Year’);
ylabel(‘Coefficient for OutputGap’);
function [beta_t, V_t] = carterKohn(Y, A, C, Q, R, beta_0, V_0, T)
n = size(A, 1); % Number of coefficients
beta_t = zeros(n, T);
V_t = zeros(n, n, T);
beta_pred = beta_0;
V_pred = V_0;
for t = 1:T
if isnan(Y(t))
continue; % Skip iterations where Y is NaN
end
% Observation update
C_t = C(t, :); % Current row of C
y_pred = C_t * beta_pred;
S = C_t * V_pred * C_t’ + R;
K = V_pred * C_t’ / S;
beta_t(:, t) = beta_pred + K * (Y(t) – y_pred);
V_t(:, :, t) = (eye(n) – K * C_t) * V_pred;
% Time update
if t < T
beta_pred = A * beta_t(:, t);
V_pred = A * V_t(:, :, t) * A’ + Q;
end
end
endHere I am currently running Bayesian time varying model in matlab and this is my whole code. I am not sure how can I create 14th and 86th percentile of my posterior distribution/credible interval that I have found using gibbs sampling and carter kohn algorithm. Furthermore, how do I plot the 14th and 86th percentile along with my mean posterior distribution in one graph?
clear;
% Load data from Excel file
Cyclicalitydata = readtable("Could be final data set.xlsx");
Countrydata = Cyclicalitydata(strcmp(Cyclicalitydata.Country, ‘Netherlands’), :);
% Prepare variables
Fiscalvariable = Countrydata.CAPB1; % Dependent variable, original series
LaggedFiscalvariable = [NaN; Fiscalvariable(1:end-1)]; % Creating one-period lag
Economicactivity = Countrydata.Outputgap; % Independent main variable
Govsize = Countrydata.Govsize; % Control variable 1
Debtp = Countrydata.DebtGDP; % Control variable 2
CPIlevel = Countrydata.CPIlevel; % Additional control variable 3
Years = Countrydata.Year; % Year for plotting and analysis
CPIgrowth = [NaN; CPIlevel(1:end-1)];
% Define the state-space model components
T = numel(LaggedFiscalvariable); % Number of observations, using the lagged variable now
C = [Economicactivity, Govsize, Debtp, CPIgrowth]; % Each row of C corresponds to an observation time
N = size(C,2);
A = eye(N); % Assuming independent evolution of each coefficient
Q = 0.01 * eye(N); % Small variance to assume minor changes over time
R = 1; % Observation noise variance
% Calculate OLS estimates for initialization
validIdx = ~isnan(LaggedFiscalvariable) & ~any(isnan(C), 2);
C_valid = C(validIdx, :);
Y_valid = LaggedFiscalvariable(validIdx);
beta_OLS = (C_valid’ * C_valid) (C_valid’ * Y_valid);
residuals = Y_valid – C_valid * beta_OLS;
sigma2_OLS = residuals’ * residuals / (length(Y_valid) – size(C_valid, 2));
% Initialize beta_0 and V_0 with OLS estimates
beta_0 = beta_OLS; % Use OLS estimates as initial states
V_0 = sigma2_OLS * eye(N); % Scale identity matrix by OLS residual variance
% Number of Gibbs sampling iterations
numIterations = 1000;
beta_samples = zeros(N, T, numIterations); % Store samples of all beta coefficients
% Gibbs Sampling Loop
for iter = 1:numIterations
% Run the Carter-Kohn algorithm for the current iteration
[beta_t, ~] = carterKohn(LaggedFiscalvariable, A, C, Q, R, beta_0, V_0, T);
% Store the sample
beta_samples(:, :, iter) = beta_t;
end
% Post-processing after Gibbs sampling
% Focus on the OutputGap coefficient’s estimates and its posterior distribution
burn_in = 100;
posterior_beta_outputgap = squeeze(mean(beta_samples(1, :, burn_in:end), 3));
% Plot only the OutputGap coefficient estimates over years
figure;
plot(Years, posterior_beta_outputgap, ‘LineWidth’, 2);
title(‘Posterior Estimates of OutputGap Coefficient over Years for the Netherlands’);
xlabel(‘Year’);
ylabel(‘Coefficient for OutputGap’);
function [beta_t, V_t] = carterKohn(Y, A, C, Q, R, beta_0, V_0, T)
n = size(A, 1); % Number of coefficients
beta_t = zeros(n, T);
V_t = zeros(n, n, T);
beta_pred = beta_0;
V_pred = V_0;
for t = 1:T
if isnan(Y(t))
continue; % Skip iterations where Y is NaN
end
% Observation update
C_t = C(t, :); % Current row of C
y_pred = C_t * beta_pred;
S = C_t * V_pred * C_t’ + R;
K = V_pred * C_t’ / S;
beta_t(:, t) = beta_pred + K * (Y(t) – y_pred);
V_t(:, :, t) = (eye(n) – K * C_t) * V_pred;
% Time update
if t < T
beta_pred = A * beta_t(:, t);
V_pred = A * V_t(:, :, t) * A’ + Q;
end
end
end Here I am currently running Bayesian time varying model in matlab and this is my whole code. I am not sure how can I create 14th and 86th percentile of my posterior distribution/credible interval that I have found using gibbs sampling and carter kohn algorithm. Furthermore, how do I plot the 14th and 86th percentile along with my mean posterior distribution in one graph?
clear;
% Load data from Excel file
Cyclicalitydata = readtable("Could be final data set.xlsx");
Countrydata = Cyclicalitydata(strcmp(Cyclicalitydata.Country, ‘Netherlands’), :);
% Prepare variables
Fiscalvariable = Countrydata.CAPB1; % Dependent variable, original series
LaggedFiscalvariable = [NaN; Fiscalvariable(1:end-1)]; % Creating one-period lag
Economicactivity = Countrydata.Outputgap; % Independent main variable
Govsize = Countrydata.Govsize; % Control variable 1
Debtp = Countrydata.DebtGDP; % Control variable 2
CPIlevel = Countrydata.CPIlevel; % Additional control variable 3
Years = Countrydata.Year; % Year for plotting and analysis
CPIgrowth = [NaN; CPIlevel(1:end-1)];
% Define the state-space model components
T = numel(LaggedFiscalvariable); % Number of observations, using the lagged variable now
C = [Economicactivity, Govsize, Debtp, CPIgrowth]; % Each row of C corresponds to an observation time
N = size(C,2);
A = eye(N); % Assuming independent evolution of each coefficient
Q = 0.01 * eye(N); % Small variance to assume minor changes over time
R = 1; % Observation noise variance
% Calculate OLS estimates for initialization
validIdx = ~isnan(LaggedFiscalvariable) & ~any(isnan(C), 2);
C_valid = C(validIdx, :);
Y_valid = LaggedFiscalvariable(validIdx);
beta_OLS = (C_valid’ * C_valid) (C_valid’ * Y_valid);
residuals = Y_valid – C_valid * beta_OLS;
sigma2_OLS = residuals’ * residuals / (length(Y_valid) – size(C_valid, 2));
% Initialize beta_0 and V_0 with OLS estimates
beta_0 = beta_OLS; % Use OLS estimates as initial states
V_0 = sigma2_OLS * eye(N); % Scale identity matrix by OLS residual variance
% Number of Gibbs sampling iterations
numIterations = 1000;
beta_samples = zeros(N, T, numIterations); % Store samples of all beta coefficients
% Gibbs Sampling Loop
for iter = 1:numIterations
% Run the Carter-Kohn algorithm for the current iteration
[beta_t, ~] = carterKohn(LaggedFiscalvariable, A, C, Q, R, beta_0, V_0, T);
% Store the sample
beta_samples(:, :, iter) = beta_t;
end
% Post-processing after Gibbs sampling
% Focus on the OutputGap coefficient’s estimates and its posterior distribution
burn_in = 100;
posterior_beta_outputgap = squeeze(mean(beta_samples(1, :, burn_in:end), 3));
% Plot only the OutputGap coefficient estimates over years
figure;
plot(Years, posterior_beta_outputgap, ‘LineWidth’, 2);
title(‘Posterior Estimates of OutputGap Coefficient over Years for the Netherlands’);
xlabel(‘Year’);
ylabel(‘Coefficient for OutputGap’);
function [beta_t, V_t] = carterKohn(Y, A, C, Q, R, beta_0, V_0, T)
n = size(A, 1); % Number of coefficients
beta_t = zeros(n, T);
V_t = zeros(n, n, T);
beta_pred = beta_0;
V_pred = V_0;
for t = 1:T
if isnan(Y(t))
continue; % Skip iterations where Y is NaN
end
% Observation update
C_t = C(t, :); % Current row of C
y_pred = C_t * beta_pred;
S = C_t * V_pred * C_t’ + R;
K = V_pred * C_t’ / S;
beta_t(:, t) = beta_pred + K * (Y(t) – y_pred);
V_t(:, :, t) = (eye(n) – K * C_t) * V_pred;
% Time update
if t < T
beta_pred = A * beta_t(:, t);
V_pred = A * V_t(:, :, t) * A’ + Q;
end
end
end percentile, quantile, bayesian time varying coefficient, plot, credible interval MATLAB Answers — New Questions
