## 4 Statistical Optimization and Model Validation
### Bayesian Inference
Similarly, the posterior predictive distribution \(P_{pred}^{\prime\prime}(y|\mathcal{Y})\) of \(\mathbf{Y}\) can be computed by averaging the conditional distribution \(P(y|\theta)\) over the posterior parameter distribution \(P(\theta|\mathcal{Y})\), defined in Equation (15):  
\[P_{pred}^{\prime\prime}(y):=\int_{\mathcal{D}_{\Theta}}p(y|\mathbf{\theta})P(\mathbf{ \theta}|\mathcal{Y})d\mathbf{\theta}=\frac{1}{P(\mathcal{Y})}\int_{\mathcal{D}_{ \Theta}}p(y|\mathbf{\theta})p(\mathcal{Y}|\mathbf{\theta})P(\mathbf{\theta})d\mathbf{\theta}. \tag{19}\]  
As mentioned earlier, the computation of the posterior distribution is nothing but trivial for arbitrary choices of the prior distributions. MCMC methods are widely used in the statistical community to approximate the posterior distribution in Equation (15). One common approach to solve such inference problems is to employ _Markov Chain Monte Carlo_ (MCMC) simulations [43; 44]. The main goal of MCMC methods is to draw \(M\) samples \(\{\mathbf{\theta}_{i}\}\) from the posterior probability density \(P(\mathbf{\theta}|\mathcal{Y})\) in Equation (15). It is worth mentioning that there is no need to compute the normalization factor \(P(\mathcal{Y})\) in Equation (15), as long as one intends to perform parameter inference task. Generating the samples \(\{\mathbf{\theta}_{i}\}\) is a non-trivial process unless \(p(\mathbf{\theta}|\mathcal{Y})\) consists of a highly specific analytic distribution, for example, a Gaussian distribution. MCMC is a procedure for generating a random walk in the parameter space that, over time, draws a representative set of samples from the distribution. Stated differently, it constructs a Markov Chain \((\mathbf{\theta}^{(1)},\mathbf{\theta}^{(2)},...)\) in the parameter space \(\mathcal{D}_{\Theta}\) with an invariant distribution that is equivalent to the target posterior distribution. For more information, the reader is referred to Reference [45]. The simplest and most widely used MCMC method is the Metropolis-Hastings (MH) algorithm [46; 47]. To efficiently sample from the posterior distribution, this algorithm pro poses and subsequently accepts or rejects candidate points. However, when the posterior distribution of parameters shows high correlation, most MCMC algorithms deliver unsatisfactory performances. To circumvent this problem, [PERSON] & Weare [48] proposed the _Affine Invariant Ensemble Sampler_ (AIES). This algorithm runs an ensemble of \(K\) Markov chains, a.k.a _walkers_\(\{\theta_{1},\ldots,\theta_{K}\}\). The walker's locations \(\theta_{i}^{(t)}\) are updated by randomly selecting a conjugate walker \(\theta_{j}^{(t)}\) from the walker sets excluding the current \(i\)-th walker (\(i\
eq j\)). A so-called _stretch move_ generates proposal that guarantee the affine invariant property of the algorithm. The candidate proposal is made by sampling from distribution \(p(z)\) defined by the tuning parameter \(a>1\):  
\[\theta_{i}^{(\star)}=\theta_{i}^{(t)}+\mathcal{Z}\Big{(}\theta_{j}^{(t)}- \theta_{i}^{(t)}\Big{)},\quad\mathcal{Z}\sim p(z)=\left\{ \begin{array}{ll} \frac{1}{\sqrt{2}\left(2\sqrt{a}-\frac{2}{\sqrt{a}}\right)}&\text{if }z\in[1/a,a]\\ 0&\text{otherwise}\end{array} \right. \tag{20}\]  
For the \(i\)-th walker in the ensemble, the candidate \(\theta_{i}^{(\star)}\) is accepted with the following probability as the new location:  
\[\alpha\left(\theta_{i}^{(\star)},\theta_{i}^{(t)},z\right)=\min\left\{1,z^{M -1}\frac{p\Big{(}\theta_{i}^{(\star)}\mid\mathcal{Y}\Big{)}}{p\Big{(}\theta_ {i}^{(t)}\mid\mathcal{Y}\Big{)}}\right\}. \tag{21}\]  
Contrary to other MCMC algorithms, AIES has only one scalar tuning parameter, namely \(a\). This parameter can be set to 2 [48; 49; 50]. In this study, we have used the Python implementation of AIES, _emcee_[51]. Judging the convergence and performance of an algorithm is nothing but trivial and there is a huge associated literature (see, for example Reference [52] for a comprehensive review). In this study, we monitor the convergence of the sampler using the integrated autocorrelation time.