Process with parameters that are themselves statistical I'd like to work with a pair of statistical processes such that the random variable from one process is the parameter of the second process. The simplest case I can imagine (and which is still appealing to me) is as such:
$$p(x,y|a,c) = p(x|a)p(y|c,x)$$
where $a, c$ are parameters and $p(x|a) = ae^{-ax}$ is the exponential distribution, and $y$ is exponential with parameter $\frac{c}{x}$,  $p(y|c, x) = (\frac{c}{x})e^{\frac{-c*y}{x}}$.
The difficulty, here, is that I can't seem to integrate over all $[0,\infty]$ and find the marginal distribution for $y$. 
Q: Has this sort of formulation (process with parameters that are themselves statistical) been pursued and are there results out there that are published?
 A: To elaborate on the discussion in the comments, your described situation can be viewed as a hierarchical model, like so:
$$X \sim Exp(a)$$
$$Y \sim Exp(\frac{c}{x})$$
There are certainly Bayesian models that provide analytic solutions, e.g. treating the $p$ parameter of a binomial distribution as a beta-distributed random variable. (These are also attractive for the interpretability of their parameters.)
But, when the relationships among variables get complicated, these methods tend to fall short. Where analytical methods are impossible, one can use a MCMC trace of any given parameter as an approximation of its posterior distribution. (Assuming one has properly assessed the chain for convergence.) Because the conjugate prior models that produce analytic solutions tend to be simple ones, their assumptions may be unrealistic for more elaborate processes. Even a conjugate normal model assumes mean depends on precision (p. 5). 
MCMC allows us to work around the absence of analytic solutions. This is a commonly cited reason for the growth in popularity of Bayesian methods: Cheap computing power made them possible.
A: Further to my comments, here's a worked example where we try to infer $a$ and $c$ on the basis of a sample of observed $x$ and $y$ pairs and the belief that the data is generated in the way you describe.
I'll use R and JAGS in this example.  First make the data in R:
set.seed(1234)
c <- 2 # true unknown values of c and a
a <- 1

Sample some $x$s according to an exponential dist with rate $a$
x <- rexp(100, a)

and for each x, sample a y from an exponential dist with rate $c/x$
y <- sapply(x, function(x){ rexp(1, rate=c/x) })

Bundle it into a data.frame
dd <- data.frame(x=x, y=y)

and define a Bayesian model by specifying a data generating process for $x$ and $y$ and priors for unknown $a$ and $c$ in the form of an R function for R2jags to work with
mod <- function(){
  # how we think the data was made
  for (i in 1:100) {
    x[i] ~ dexp(a)
    y[i] ~ dexp(c / x[i])
  }
  # prior dists for unknown a and c
  a ~ dgamma(1, 2) 
  c ~ dunif(0, 10)
}

Now load the interface to JAGS
library(R2jags) # requires JAGS to be installed first

and compile a sampler that will sample from the posterior for this model
jmod <- jags(dd, 
             parameters.to.save=c('c', 'a'), 
             model.file=mod)

An initial round of sampling is actually sufficient (convergence statistic is here Rhat) so we look what the sampler has come up with
jmod

3 chains, each with 2000 iterations (first 1000 discarded)
n.sims = 3000 iterations saved
     mu.vect sd.vect    2.5%     25%     50%     75%   97.5%  Rhat n.eff
a          1.011   0.101   0.824   0.943   1.010   1.076   1.219 1.001  3000
c          1.954   0.194   1.592   1.819   1.948   2.079   2.351 1.001  3000

The marginal posterior means of $a$ and $b$ are 1.011 and 1.954 respectively with quantiles as shown above, etc.  In this case not too far from their true values of 1 and 2.  
You can also summarise the joint posterior in any way you like from the samples themselves, but this shows the general approach.
