%% Estimate Conditional Posterior Distribution
% Consider the linear regression model in Create Custom Multivariate t Prior Model for Coefficients.
% Create an anonymous function that operates like priorMVTIG, but accepts the parameter values only and holds the hyperparameter values fixed.

dof = 50;
C = eye(4);
ct = [-25; 4; 0; 3];
st = ones(4,1);
a = 3;
b = 1;
logPDF = @(params)priorMVTIG(params,ct,st,dof,C,a,b);

%
% Create a custom joint prior model for the linear regression parameters. Specify the number of predictors p. Also, specify the function handle for priorMVTIG and the variable names.

p = 3;
PriorMdl = bayeslm(p,'ModelType','custom','LogPDF',logPDF,...
    'VarNames',["IPI" "E" "WR"])

%
% Load the Nelson-Plosser data set. Create variables for the response and predictor series.

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

%
% Estimate the conditional posterior distribution of \beta given the data and \sigma^2 = 2, and return the estimation summary table to access the estimates. Specify a width for the slice sampler that is close to the posterior standard deviation of the parameters assuming a diffuse prior model. Reduce serial correlation by specifying a thinning factor of 10, and reduce the effective default number of draws by a factor of 10.

width = [20,0.5,0.01,1];
thin = 10;
numDraws = 1e5/thin;
rng(1) % For reproducibility
[Mdl,Summary] = estimate(PriorMdl,X,y,'Sigma2',2,...
    'Width',width,'Thin',thin,'NumDraws',numDraws);

% estimate displays a summary of the conditional posterior distribution of \beta. Because \sigma^2 is fixed at 2 during estimation, inferences on it are trivial.
% Extract the mean vector and covariance matrix of the conditional posterior of \beta from the estimation summary table.

condPostMeanBeta = Summary.Mean(1:(end - 1))
CondPostCovBeta = Summary.Covariances(1:(end - 1),1:(end - 1))

%
% Display Mdl.

Mdl

% Because estimate computes the conditional posterior distribution, it returns the original prior model, not the posterior, in the first position of the output argument list. Also, estimate does not return the MCMC sample. Therefore, to monitor convergence of the MCMC sample, use simulate instead and specify the same random number seed.
% Copyright 2012 The MathWorks, Inc.