What is the distribution of a nested Laplace? We have three random variables $X_1, X_2, X_3$ having the following conditional distributions:
\begin{align}
p(X_2 \mid X_1=x_1) &\sim \mathrm{Laplace}(X_2;x_1, b) = \frac{1}{2b} \exp(-\frac{|X_2 - x_1|}{b}) \\
p(X_3 \mid X_2=x_2) &\sim \mathrm{Laplace}(X_3;x_2, b) = \frac{1}{2b} \exp(-\frac{|X_3 - x_2|}{b}) 
\end{align}
We want to find the conditional distribution $p(X_3 \mid X_1 = x_1)$. We can write the following:
$$p(X_3 \mid X_1) = \sum\limits_{X_2} p(X_3 \mid X_2) p(X_2 \mid X_1)$$
However, I don't know how to compute this because of the absolute values. I also thought I can write $X_2$ and $X_3$ as
$$X_3 = X_2 \pm \sigma_2~, \qquad X_2 = X_1 \pm \sigma_1$$
Thus we have $X_3 = X_1 \pm \sigma_1 \pm \sigma_2$, where $\sigma_1$ and $\sigma_2$ are i.i.d. exponential random variables and $\pm$ has a $50\%$ chance of being positive and a $50\%$ chance of being negative. 
Does $p(X_3 \mid X_1=x_1)$ have a closed from solution? Perhaps, $\mathrm{Laplace}(X3\,;~x_1, 2b)$? 
Update:
I did the following experiment:
Chose $X_1 = 5$, got 100,000 samples from $p(X_2|X_1 = 5, b = 1)$ with Laplace distribution, got one sample from $p(X_3|X_2, b = 1)$ for all 100,000 values of $X_2$ and plotted the histograms. I also fitted a Laplace distribution to $p(X_3|X_1)$ where $b \approx 1.5$.

It seems that although Laplace is a good fit for the tails of the data, it doesn't fit the data well around zero. Maybe suggesting that the nested Laplace distribution is not Laplace? 
 A: $$\begin{align*}
f_{X_3 \mid X_1}(x_3 \mid x_1) 
&= \int_{x_2=-\infty}^\infty f_{X_3 \mid X_2}(x_3 \mid x_2) f_{X_2 \mid X_1}(x_2 \mid x_1) \, dx_2 \\
&= \frac{1}{4b^2} \int_{x_2 = -\infty}^\infty e^{-|x_3 - x_2|/b} e^{-|x_2 - x_1|/b} \, dx_2 \\
&= \frac{1}{4b^2} \int_{x_2 = -\infty}^\infty e^{-(|x_3 - x_2|+|x_2 - x_1|)/b} \, dx_2.
\end{align*}$$
At this point, it becomes clear that the integration will depend on whether $x_3 > x_1$ or $x_3 < x_1$.  Without loss of generality, suppose the former.  Then note $$|x_3 - x_2| + |x_2 - x_1| = \begin{cases}
x_1 + x_3 - 2x_2, & x_2 \in (-\infty, x_1) \\ 
x_3 - x_1, & x_2 \in [x_1, x_3] \\
2x_2 - x_1 - x_3, & x_2 \in (x_3, \infty).
\end{cases}$$  So when $x_3 > x_1$, 
$$\begin{align*}
f_{X_3 \mid X_1}(x_3 \mid x_1) 
&= \frac{1}{4b^2} \left( \int_{x_2 = -\infty}^{x_1} \!\!\!\! e^{-(x_1 + x_3)/b} e^{2x_2/b} \, dx_2 + \int_{x_2 = x_1}^{x_3} \!\!\!\! e^{(x_1 - x_3)/b} \, dx_2 + \int_{x_2 = x_3}^\infty \!\!\!\! e^{(x_1 + x_3)/b} e^{-2x_2/b} \, dx_2 \right) \\
&= \frac{1}{4b^2} \left( \frac{b}{2}e^{-(x_1 + x_3)/b} e^{2x_1/b} + (x_3 - x_1) e^{(x_1 - x_3)/b} + \frac{b}{2} e^{(x_1 + x_3)/b} e^{-2x_3/b} \right) \\
&= \frac{1}{8b} \left( e^{-(x_3 - x_1)/b} + \frac{2}{b} (x_3 - x_1)e^{-(x_3 - x_1)/b} + e^{-(x_3 - x_1)/b}  \right) \\
&= \frac{1}{4b} \left( 1 + \frac{x_3 - x_1}{b} \right)e^{-(x_3 - x_1)/b}.
\end{align*}$$
When $x_3 < x_1$, then the roles of $x_1$ and $x_3$ are switched, thus we can write the full conditional density as
$$f_{X_3 \mid X_1}(x_3 \mid x_1) = \frac{1}{4b} \left( 1 + \frac{|x_3 - x_1|}{b}\right) e^{-|x_3 - x_1|/b}, \quad x_3 \in (-\infty, \infty).$$
Clearly this is not Laplace but it is in a sense analogous to how a gamma distribution is a generalization of the exponential distribution.  With this expression, I think now you should be able to use simulation to confirm it is correct.
