# MCMC diagnostics

Explore the difference between MC and MCMC with a sample example.

A discrete variable $$\delta\sim{1,2,3}$$ and a continuous variable $$\theta\in\mathrm{I!R}$$.

```r
# discrete variable
delta = c(.45, .10, .45)

# continuous variable
mu = c(-3, 0, 3)
sigma2 = c(1/3, 1/3, 1/3)

# exact marginal density of theta
ext_margin_den <- function(x)
{
  dnorm(x, mu[1], sqrt(sigma2[1])) * delta[1] +
    dnorm(x, mu[2], sqrt(sigma2[2])) * delta[2] +
    dnorm(x, mu[3], sqrt(sigma2[3])) * delta[3]
}

theta = seq(-6, 6, length.out = 1000)
ptheta = ext_margin_den(theta)
plot(theta, ptheta, type = "l")
```

The marginal density of $$\theta$$ would be

![](/files/-LNYSZXgH7NGCDQh3PBD)

Firstly, generate 1000 Monte Carlo $$\theta$$-samples.

![](/files/-LNYSZXi4E7kFPCGk5s3)

```r
sample_delta <- function(x)
{
  denx = c(delta[1] * dnorm(x, mu[1], sqrt(sigma2[1])),
           delta[2] * dnorm(x, mu[2], sqrt(sigma2[2])),
           delta[3] * dnorm(x, mu[3], sqrt(sigma2[3])))
  condenx = denx / sum(denx)
  sample(1:3, 1, prob = condenx)
}
# gibbs sampler
gibbs <- function(m, delta = 3)
{
  DELTA = delta
  THETA = NULL
  for (i in 1:m)
  {
    # sample theta from full conditional distribution
    theta = sample_theta(delta)
    THETA = c(THETA, theta)
    delta = sample_delta(theta)
    DELTA = c(DELTA, delta)
  }
  return(THETA)
}
res = gibbs(1000)
plot(res)
hist(res, breaks = 20)
```

For 1000 MCMC samples, we have

![](/files/-LNYSZXkbs5FJkxJPXhE)

Actually, re-run the above code, we can get much different figure. Let's try 10000 MCMC samples,

![](/files/-LNYSZXmOjy5UI06r3Wz)

It turns out to be much stable when you re-run the above code.

## sample autocorrelation

Use R-function `acf`. If a Markov chain with high autocorrelation, then it will move around the parameter space slowly, taking a long time to achieve the correct balance among the different regions of the parameter space.

## effective sample size

Use R command `effectiveSize` in the `coda` package, which can be interpreted as the number of independent Monte Carlo samples necessary to give the same precision as the MCMC samples.


---

# Agent Instructions: Querying This Documentation

If you need additional information that is not directly available in this page, you can query the documentation dynamically by asking a question.

Perform an HTTP GET request on the current page URL with the `ask` query parameter:

```
GET https://mc.hohoweiya.xyz/mcmc/mcmcdiagnostics.md?ask=<question>
```

The question should be specific, self-contained, and written in natural language.
The response will contain a direct answer to the question and relevant excerpts and sources from the documentation.

Use this mechanism when the answer is not explicitly present in the current page, you need clarification or additional context, or you want to retrieve related documentation sections.