Let $\mu_0 (dx)$ be the standard one-dimensional Cauchy distribution, i.e.

\begin{align} \mu_0 (dx) = \frac{1}{\pi} \frac{1}{1+x^2} dx. \end{align}

Suppose I fix $h \in [0, 1]$, and form a Markov chain $\{X_n\}_{n \geqslant 0}$ as follows:

  1. At step $n$, I sample $Y_n \sim \mu_0$.
  2. I then set $X_{n+1} = (1 - h) X_n + h Y_n$

It is not so hard to show that this chain admits $\mu_0$ as a stationary measure, as this essentially comes from the fact that Cauchy distribution is a stable distribution.

What I'm interested in is the reversal of this Markov chain. More precisely, if the chain I describe above uses the Markov kernel $q (x \to dy)$, I want to understand the Markov kernel $r$ such that

\begin{align} \mu_0 (dx) q (x \to dy) = \mu_0 (dy) r (y \to dx). \end{align}

Fortunately, all of the quantities involved have densities with respect to Lebesgue measure, and as such, I can write down what $r (y \to dx)$ is:

\begin{align} r (y \to dx) = \frac{1 + y^2}{\pi} \frac{1}{1 + x^2} \frac{h}{ h^2 + (y - (1 - h) x)^2} dx. \end{align}

My question is then: is there a simple, elegant way to draw exact samples from $r$?

I would highlight that this is not a purely algorithmic question; I'd really like to understand what this reversal kernel $r$ is doing. A nice byproduct of that would then be that I could simulate from it easily.

For completeness, some of the `purely algorithmic' solutions I had considered were the following.

  • I could try rejection sampling, and in principle this would work, but it wouldn't really give me insight into the nature of the Markov chain.
  • I could try something like the inverse CDF method, but it seems to me that the CDF of $r$ is not particularly nice to work with. As such, I'd have to use e.g. Newton iterations to use this method, and I'd prefer to not have to do this.
  • 1
    $\begingroup$ The general Cauchy distribution $c_w(dx)=\frac{b dx}{\pi(b^2+(x-a)^2)}$ where $b>0$ and $w=a+ib$ has Fourier transform $e^{iwt}$ for $t>0.$ Therefore your Markov chain defined by $X_{n+1}=(1-h)X_n+hY_n$ where $Y_n\sim c_i=\mu_0$ is such that $X_n\sim c_{w_n}$ where $w_{n+1}=(1-h)w_n+ih$, hence $w_n=(1-h)^nw_0+ih$. I will try to understand your question in terms of these $w$ and $c_w.$ $\endgroup$ Commented Aug 1, 2023 at 12:05


You must log in to answer this question.

Browse other questions tagged .