  We propose a multilevel Monte Carlo-FEM algorithm to solve elliptic Bayesian
inverse problems with "Besov random tree prior". These priors are given by a
wavelet series with stochastic coefficients, and certain terms in the expansion
vanishing at random, according to the law of so-called Galton-Watson trees.
This allows to incorporate random fractal structures and large deviations in
the log-diffusion, which occur naturally in many applications from geophysics
or medical imaging. This framework entails two main difficulties: First, the
associated diffusion coefficient does not satisfy a uniform ellipticity
condition, which leads to non-integrable terms and thus divergence of standard
multilevel estimators. Secondly, the associated space of parameters is Polish,
but not a normed linear space. We address the first point by introducing
cut-off functions in the estimator to compensate for the non-integrable terms,
while the second issue is resolved by employing an independence
Metropolis-Hastings sampler. The resulting algorithm converges in the
mean-square sense with essentially optimal asymptotic complexity, and
dimension-independent acceptance probabilities.
