Analytic continuation of quantum Monte Carlo data by stochastic analytical inference

We present an algorithm for the analytic continuation of imaginary-time quantum Monte Carlo data which is strictly based on principles of Bayesian statistical inference. Within this framework we are able to obtain an explicit expression for the calculation of a weighted average over possible energy spectra, which can be evaluated by standard Monte Carlo simulations, yielding as by-product also the distribution function as function of the regularization parameter. Our algorithm thus avoids the usual ad hoc assumptions introduced in similar algorithms to fix the regularization parameter. We apply the algorithm to imaginary-time quantum Monte Carlo data and compare the resulting energy spectra with those from a standard maximum-entropy calculation.

Figure 1

(Color online) Simulated spectra for a range of regularization parameters $\alpha$ $\left(\beta =14W^{-1}\right)$. For large $\alpha$ the Gaussian shape of the default model is visible. For decreasing $\alpha$ several features begin to appear.

Figure 2

(Color online) The probability distributions $P\left[\alpha \mid \overline{G}\right]$ (a) based on a Wang-Landau simulation of the density of states (b). The different choices for $P\left[\alpha \right]$ only have a weak influence on the position of the peak. The resulting spectra (c) are calculated by either averaging all spectra over $P\left[\alpha \mid \overline{G}\right]$ or by taking the spectrum that maximizes them. The four different spectra are practically identical.

Figure 3

(Color online) The double-logarithmic plot of ${\chi }^2$ (a) shows a kink at $\alpha \sim 0.8$ which is very close to the choice of Syljuåsen $\left(\alpha =1\right)$. It is also the region where ${\chi }^2\sim N$. Thus, all three methods give about the same answer and the resulting spectra (c) are identical. A comparison shows that this solution is quite close to the spectrum shown in Fig. 2c. The entropy (b) exhibits no significant features and gives at least for this data set no indication how to determine $\alpha$.

Figure 4

(Color online) Results of a maximum-entropy calculation for the same QMC data as in Figs. 2, 3. The probability distribution $P\left[\alpha \mid \overline{G}\right]$ (a) shows a noticeable dependence on $P\left[\alpha \right]$, but analogous to SAI, the resulting spectra (c) are identical. The $\alpha$ where ${\chi }^2\sim N$ (b) is larger $\left(\alpha \sim 500\right)$ than the one for which $P\left[\alpha \mid \overline{G}\right]$ is maximal $\left(\alpha \sim 90\right)$. Accordingly the spectrum chosen by the ${\chi }^2\sim N$ rule is more regularized than the one calculated by Bayesian inference.

Figure 5

(Color online) Spectra simulated by stochastic analytical inference compared to maximum-entropy calculations. All calculations are based on a Gaussian default model

Figure 6

(Color online) Spectra for $\beta =34W^{-1}$ based on a flat default model calculated by SAI (a) and by the MEM (b). These spectra are only marginally different compared to those calculated using a Gaussian default model. We conclude that for the QMC data under consideration the calculated spectra are quite independent of the default model.

Figure 7

(Color online) Results for a square-root-shaped spectrum. Only positive frequencies are shown since the spectrum is symmetric around zero. This symmetry was enforced for all simulations. The SAI including the Wang-Landau criterion for the determination of $\alpha$ is the most successful in resolving the sharply peaked input spectrum. Note that the choice $\alpha =1$ gives the same result as ${\chi }^2\sim N$.