diagnostics
Class: HamiltonianSampler
Markov Chain Monte Carlo diagnostics
Syntax
tbl = diagnostics(smp,chains)
tbl = diagnostics(smp,chains,'MaxLag',maxlag)
Description
returns
Markov Chain Monte Carlo diagnostics for the chains in tbl
= diagnostics(smp
,chains
)chains
.
specifies
the maximum number of autocorrelation lags to use for computing effective
sample sizes.tbl
= diagnostics(smp
,chains
,'MaxLag',maxlag
)
Input Arguments
smp
— Hamiltonian Monte Carlo sampler
HamiltonianSampler
object
Hamiltonian Monte Carlo sampler, specified as a HamiltonianSampler
object.
Use the hmcSampler
function
to create a sampler.
chains
— MCMC chains
matrix | cell array
MCMC chains, specified as one of the following:
A matrix, where each row is a sample and each column a parameter.
A cell array of matrices, where the chain
chains{i}
is a matrix where each row is a sample and each column a parameter.
The number of parameters (that is, matrix columns) must equal
the number of elements of the StartPoint
property
of the smp
sampler.
maxlag
— Maximum number of autocorrelation lags
100
(default) | positive integer
Maximum number of autocorrelation lags for computing effective sample sizes, specified as a positive integer.
The effective sample size calculation uses lags of 1,2,...,maxlag
for
each chain in chains
that has more than maxlag
samples.
For chains with maxlag
or fewer samples,
the calculation uses Ni -
1 lags, where Ni is the
number of samples of chain i.
Example: 'MaxLag',50
Output Arguments
tbl
— MCMC diagnostics
table
MCMC diagnostics, computed using all the chains in chains
and
returned as a table with these columns.
Column | Description |
---|---|
Name | Variable name |
Mean | Posterior mean estimate |
MCSE | Estimate of the Monte Carlo standard error (the standard deviation of the posterior mean estimate) |
SD | Estimate of the posterior standard deviation |
Q5 | Estimate of the 5th quantile of the marginal posterior distribution |
Q95 | Estimate of the 95th quantile of the marginal posterior distribution |
ESS | Effective sample size for the posterior mean estimate. Larger effective sample sizes lead to more accurate results. If the samples are independent, then the effective sample size is equal to the number of samples. |
RHat | Gelman-Rubin convergence statistic. As a rule of thumb,
values of |
Examples
Compute Markov Chain Monte Carlo Diagnostics
Create MCMC chains using a Hamiltonian Monte Carlo (HMC) sampler and compute MCMC diagnostics.
First, save a function on the MATLAB® path that returns the multivariate normal log probability density and its gradient. In this example, that function is called normalDistGrad
and is defined at the end of the example. Then, call this function with arguments to define the logpdf
input argument to the hmcSampler
function.
means = [1;-2;2]; standevs = [1;2;0.5]; logpdf = @(theta)normalDistGrad(theta,means,standevs);
Choose a starting point. Create the HMC sampler and tune its parameters.
startpoint = randn(3,1); smp = hmcSampler(logpdf,startpoint); smp = tuneSampler(smp);
Draw samples from the posterior density, using a few independent chains. Choose different, randomly distributed starting points for each chain. Specify the number of burn-in samples to discard from the beginning of the Markov chain and the number of samples to generate after the burn-in.
NumChains = 4; chains = cell(NumChains,1); Burnin = 500; NumSamples = 1000; for c = 1:NumChains chains{c} = drawSamples(smp,'Burnin',Burnin,'NumSamples',NumSamples,... 'Start',randn(size(startpoint))); end
Compute MCMC diagnostics and display the results. Compare the true means in means
with the column titled Mean
in the MCMCdiagnostics
table. The true posterior means are within a few Monte Carlo standard errors (MCSEs) of the estimated posterior means. The HMC sampler has accurately recovered the true means. Similarly, the estimated standard deviations in the column SD
are very near the true standard deviations in standev
.
MCMCdiagnostics = diagnostics(smp,chains)
MCMCdiagnostics=3×8 table
Name Mean MCSE SD Q5 Q95 ESS RHat
______ _______ ________ _______ ________ ______ ______ ____
{'x1'} 1.0038 0.016474 0.96164 -0.58601 2.563 3407.4 1
{'x2'} -2.0435 0.034933 1.999 -5.3476 1.1851 3274.5 1
{'x3'} 1.9957 0.008209 0.49693 1.2036 2.8249 3664.5 1
means
means = 3×1
1
-2
2
standevs
standevs = 3×1
1.0000
2.0000
0.5000
The normalDistGrad
function returns the logarithm of the multivariate normal probability density with means in Mu
and standard deviations in Sigma
, specified as scalars or columns vectors the same length as the startpoint
. The second output argument is the corresponding gradient.
function [lpdf,glpdf] = normalDistGrad(X,Mu,Sigma) Z = (X - Mu)./Sigma; lpdf = sum(-log(Sigma) - .5*log(2*pi) - .5*(Z.^2)); glpdf = -Z./Sigma; end
Tips
After creating an HMC sampler using the
hmcSampler
function, you can compute MAP (maximum-a-posteriori) point estimates, tune the sampler, draw samples, and check convergence diagnostics using the methods of theHamiltonianSampler
class. For an example of this workflow, see Bayesian Linear Regression Using Hamiltonian Monte Carlo.
Version History
Introduced in R2017a
See Also
Functions
Classes
Comando de MATLAB
Ha hecho clic en un enlace que corresponde a este comando de MATLAB:
Ejecute el comando introduciéndolo en la ventana de comandos de MATLAB. Los navegadores web no admiten comandos de MATLAB.
Select a Web Site
Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .
You can also select a web site from the following list:
How to Get Best Site Performance
Select the China site (in Chinese or English) for best site performance. Other MathWorks country sites are not optimized for visits from your location.
Americas
- América Latina (Español)
- Canada (English)
- United States (English)
Europe
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)