Decomposing Laplace prior for a hierarchical representation of a model I have a conditional Laplace prior:
$$
\pi(\boldsymbol{\beta}|\sigma^2) = \prod\limits_{j=1}^{p}\frac{\lambda}{2\sqrt{\sigma^2}}e^{-\lambda|\beta_j|/\sqrt{\sigma^2}}
$$
and a marginal prior on $\sigma^2$, $\pi(\sigma^2)$.
I want to decompose this Laplace prior for a hierarchical representation of a model.
I know that that the following relations will be necessary:
\begin{equation}
f(x|z)=\int_{Y}f(x|y)f(y|z)dy
\end{equation}
and 
$$
\frac{a}{2}e^{-a|z|} = \int_0^\infty \frac{1}{\sqrt{2\pi s}}e^{-z^2/(2s)}\frac{a^2}{2}e^{-a^2s/2} ds, \ \ \ \ \ \ a>0.
$$
However, I have problems with some details. According to the paper Casella, Park (2008) I should get that, this
\begin{align*}
\boldsymbol{\beta} | \sigma^2, \tau_1^2, \dots, \tau_p^2 &\sim N_p(\mathbf{0}_p, \sigma^2\mathbf{D}_{\boldsymbol{\tau}}), \\
\mathbf{D}_{\boldsymbol{\tau}} &= \text{diag}(\tau_1^2, \dots, \tau_p^2), \\
\sigma^2, \tau_1^2, \dots, \tau_p^2 &\sim \pi(\sigma^2) d\sigma^2 \prod\limits_{j=1}^{p}\frac{\lambda^2}{2}e^{-\lambda^2\tau_j^2/2}d\tau_j^2
\end{align*}
represents the conditional Laplace prior, however I still end up with something slightly different. Can someone tell me how to get there (e.g. by telling me, plug $z = \beta_j$ etc.)? Thank you for your help!
 A: There are a few errors here.  Basically you put the $\sigma^2$ in the wrong place.
Note that technically speaking you should really have $f(x|yz)$ and not $f(x|y)$ in the integral formula.  Although we will see that for your specific representation the result as you have written it is still true using $f(x|y)$.
The distribution $p(\beta|\sigma^2,\tau_1^2,\dots,\tau_p^2)$ and $p(\sigma^2,\tau_1^2,\dots,\tau_p^2)$ are both inconsistent with the integral formula you provide.  The first mistake is that the latent variable distribution should also be conditional on $\sigma^2$, so you should have $p(\tau_1^2,\dots,\tau_p^2|\sigma^2)$.  In relating this to the general integral forumla you have $z=\sigma^2,y=\tau_j^2,x=\beta_j$ for $j=1,...,p$.
In your definition of $p(\beta|\sigma^2)$, the value of $a$ should be $a=\frac{\lambda}{\sqrt{\sigma^2}}=\frac{\lambda}{\sigma}$ in order for you to use the integral result.  Now note that $a$ does not appear in the term $\frac{z^2}{2s}$ inside your integral result.  This is the normal distribution part so you can remove $\sigma^2$ from the distribution $p(\beta|\sigma^2,\tau_1^2,\dots,\tau_p^2)$.  So you have:
$$
\begin{align*}
\boldsymbol{\beta} | \sigma^2, \tau_1^2, \dots, \tau_p^2 &\sim N_p(\mathbf{0}_p, \mathbf{D}_{\boldsymbol{\tau}}), \\
\mathbf{D}_{\boldsymbol{\tau}} &= \text{diag}(\tau_1^2, \dots, \tau_p^2)
\end{align*}
$$
This confirms the previous note that $f(x|z)=\int_{Y}f(x|y)f(y|z)$ holds for your specific representation.  For a different representation where this isn't true, make the change of variables $s=ra$ in the integral.  The other term is straigh-forward, and an exponential density with rate parameter $\frac{a^2}{2}=\frac{\lambda^2}{2\sigma^2}$.  So the product should be:
$$f(\tau_1^2,\dots,\tau_p^2|\sigma^2)=\prod\limits_{j=1}^{p}\frac{\lambda^2}{2\sigma^2}\exp\left(-\frac{\lambda^2}{2\sigma^2}\tau_j^2\right)$$
UPDATE
On having a bit more think about it, the end result you have is actually correct, because if $\tau_j^2\sim Expo\left(\frac{\lambda^2}{2\sigma^2}\right)$ then we also have that $\sigma^2\tau_j^2\sim Expo\left(\frac{\lambda^2}{2}\right)$.  This shows that the result is correct.  You can also make a change of variables in your integral definition of $s=r\sigma^2$ and you have:
$$\frac{a}{2}e^{-a|z|} = \int_0^\infty \frac{1}{\sqrt{2\pi r\sigma^2}}e^{-z^2/(2r\sigma^2)}\frac{a^2\sigma^2}{2}e^{-a^2r\sigma^2/2} dr, \ \ \ \ \ \ a>0$$
This is consistent with your answer, by setting $r=\tau_j^2$
