Numerical analytic continuation: answers to well-posed questions
Olga Goulko
Department of Physics, University of Massachusetts, Amherst, MA 01003, USA
Andrey S. Mishchenko
RIKEN Center for Emergent Matter Science (CEMS),
2-1 Hirosawa, Wako, Saitama, 351-0198, Japan and
National Research Center “Kurchatov Institute,” 123182 Moscow, Russia
Lode Pollet
Department of Physics, Arnold Sommerfeld Center for Theoretical Physics,
arXiv:1609.01260v1 [cond-mat.other] 5 Sep 2016
University of Munich, Theresienstrasse 37, 80333 Munich, Germany
Nikolay Prokof’ev
Department of Physics, University of Massachusetts, Amherst, MA 01003, USA
Department of Physics, Arnold Sommerfeld Center for Theoretical Physics,
University of Munich, Theresienstrasse 37, 80333 Munich, Germany and
National Research Center “Kurchatov Institute,” 123182 Moscow, Russia
Boris Svistunov
Department of Physics, University of Massachusetts, Amherst, MA 01003, USA
National Research Center “Kurchatov Institute,” 123182 Moscow, Russia and
Wilczek Quantum Center, Zhejiang University of Technology, Hangzhou 310014, China
(Dated: September 6, 2016)
We formulate the problem of numerical analytic continuation in a way that lets us draw meaning-
ful conclusions about properties of the spectral function based solely on the input data. Apart from
ensuring consistency with the input data (within their error bars) and the a priori and a posteriori
(conditional) constraints, it is crucial to reliably characterize the accuracy—or even ambiguity—of
the output. We explain how these challenges can be met with two approaches: stochastic optimiza-
tion with consistent constraints and the modified maximum entropy method. We perform illustrative
tests for spectra with a double-peak structure, where we critically examine which spectral properties
are accessible and which ones are lost. For an important practical example, we apply our protocol
to the Fermi polaron problem.
PACS numbers: 02.70.-c, 71.15.Dx, 71.28.+d, 71.10.Fd
I. INTRODUCTION The NAC problem is often characterized as ill-posed.
Mathematically, the near-degeneracy of the kernel im-
Numerous problems in science, from spectral analysis plies two closely related circumstances: (i) the absence of
to image processing, require that we restore properties of the resolvent, and (ii) a continuum of solutions satisfying
a function A(z) from a set of integrals the input data within their error bars (even when inte-
Z ∞ grals over z are replaced with finite sums containing less
or equal to N terms). Nowadays, the first circumstance
gn = G[n, A] ≡ dzK(n, z)A(z), n = 1, . . . , N, (1)
−∞ is merely a minor technical problem, as there exists a
number of methods allowing one to find solutions to (1)
where K(n, z) is a known kernel and {gn } is a finite set of without compromising the error bars of gn .
experimental or numerical input data with error bars. An The second circumstance—the ambiguity of the
important class of such problems—known as numerical solution—is a more essential problem. It is clear that
analytic continuation (NAC)—deals with “pathological” if one formulates the goal as to restore A(z) as a con-
kernels featuring numerous eigenfunctions with anoma- tinuous curve, or to determine its value on a given grid
lously small eigenvalues. An archetypal NAC problem of points, then the goal cannot be reached as stated, ir-
is the numerical spectral analysis at zero temperature, respective of the properties of the kernel K(n, z). The
where the challenge is to restore the non-negative spec- input data set is finite and noisy, thereby introducing a
tral function A(z ≥ 0) satisfying the equation natural limit on the resolution of fine structures in A(z).
Z ∞ Fortunately, the above-formulated goal has little to do
gn = dze−zτn A(z), (2) with the practical world. In an experiment, all devices
0
are characterized by a finite resolution function and the
from numerical data for gn = g(τn ≥ 0). data they collect always correspond to integrals. The
, 2
R
data are processed by making certain assumptions about dzA(z) is typically known with an accuracy that is or-
the underlying function. This motivates an alternative ders of magnitude better than what would be predicted
formulation of the NAC goal involving integrals of A(z) by the central limit theorem if this integral is represented
that render the problem well-defined. With additional by a finite sum of integrals over nonoverlapping inter-
assumptions about the smoothness and other properties vals. Second, the errors are not necessarily distributed
of A(z) behind these integrals, consistent with both a as a Gaussian. Atypical fluctuations can have a signifi-
priori and a posteriori knowledge, the ambiguity of the cant probability and their analysis should not be avoided
solution can be substantially suppressed. The simplest as the actual physical solution may well be one of them.
setup is as follows: To this end, it is important to explore the minimal and
maximal values that the integral im can take, and check
Given a set of finite intervals {∆m }, determine the inte- that these are not significantly different from the typical
grals of the spectral function over these intervals: value of im . In certain cases this criterion cannot be met
Z without increasing the intervals ∆m to an extent when
im = ∆−1 dzA(z), m = 1, . . . , M, (3) the assumption of linearity of A(z) becomes uncontrolled,
m
z∈∆m implying that an important piece of information about
the shape of A(z) in this interval is missing. A charac-
along with the corresponding dispersions of fluctuations teristic example that plays a key role in the subsequent
{σm } (straightforwardly extendable to the dispersion discussion is presented in Fig. 1, where the challenge is
correlation matrix {σmm′ }). to extract the shape of the second peak.
Naively one might think that nothing is achieved by going
from the integrals in (1) to the integrals in (3) because A(z)
2.5
the latter have exactly the same form with the kernel
K̄(m, z) = ∆−1m for z ∈ ∆m and zero otherwise (other
forms of the “resolution function” K̄(m, z) are discussed 2.0
in Sec. II):
Z ∞
1.5
im = I[m, A] = dz K̄(m, z)A(z), m = 1, . . . , M.
−∞
1.0
(4)
This impression, however, is false because kernel proper- 0.5
ties are at the heart of the problem. If for appropriately
small intervals (sufficiently small to resolve the variations
0.0
of A(z)), the uncertainties for im remain small, then one 0 1 2 3 Z
can draw reliable conclusions for the underlying behav-
ior of A(z) itself. The difference between “good” (e.g. FIG. 1. Challenging example of a spectrum A(z) for the NAC
as in Fourier transforms) and pathological kernels is that problem (2). As shown in the main part of the text, the signif-
for the latter, due to the notorious saw-tooth instability, icant width of the first peak makes it essentially impossible to
the uncertainties for im quickly become too large for a controllably restore the width of the second peak, even with
meaningful analysis of fine structures in A(z). small relative error bars (∼ 10−5 ) on gn . On the other hand,
To obtain a solution from the integrals (3), one has the first two moments of the second peak, characterizing its
to invoke the notion of conditional knowledge. The most weight and position, can be extracted reliably.
straightforward approach is to set the spectral function
values at the middle points zm of the intervals ∆m to In the more sophisticated approach used in this work,
Afin (zm ) = im . This is only possible if the intervals can the values of Afin (zm ) can be further optimized (without
be made appropriately narrow without losing accuracy compromising the accuracy of the solution with respect
for the integrals. With this approach we assume that to gn ) to produce a smooth curve. This protocol has
the function is nearly linear over the intervals in ques- the additional advantage of eliminating minima, maxima,
tion. This is a typical procedure for experimental data. gaps, and sharp features that are not guaranteed to exist
Quantifying the error bar on Afin (zm ) necessarily involves by the quality of input data. The nature of the problem is
two numbers: the “vertical” dispersion σm is directly in- such that very narrow peaks (or gaps) with tiny spectral
herited from im , and the “horizontal” error bar ∆m /2 weight can always be imagined to be present (for narrow
represents the interval half-width. intervals they will certainly emerge due to the saw-tooth
The reader should be aware of two issues regard- instability). Our philosophy with respect to these fea-
ing such error bars. First, the error bars for different tures is to erase them within the established error bounds
points are not independent but contain significant multi- and obtain a solution that is insensitive to the interval
point correlations. For example, an unrestricted integral parameters.
Olga Goulko
Department of Physics, University of Massachusetts, Amherst, MA 01003, USA
Andrey S. Mishchenko
RIKEN Center for Emergent Matter Science (CEMS),
2-1 Hirosawa, Wako, Saitama, 351-0198, Japan and
National Research Center “Kurchatov Institute,” 123182 Moscow, Russia
Lode Pollet
Department of Physics, Arnold Sommerfeld Center for Theoretical Physics,
arXiv:1609.01260v1 [cond-mat.other] 5 Sep 2016
University of Munich, Theresienstrasse 37, 80333 Munich, Germany
Nikolay Prokof’ev
Department of Physics, University of Massachusetts, Amherst, MA 01003, USA
Department of Physics, Arnold Sommerfeld Center for Theoretical Physics,
University of Munich, Theresienstrasse 37, 80333 Munich, Germany and
National Research Center “Kurchatov Institute,” 123182 Moscow, Russia
Boris Svistunov
Department of Physics, University of Massachusetts, Amherst, MA 01003, USA
National Research Center “Kurchatov Institute,” 123182 Moscow, Russia and
Wilczek Quantum Center, Zhejiang University of Technology, Hangzhou 310014, China
(Dated: September 6, 2016)
We formulate the problem of numerical analytic continuation in a way that lets us draw meaning-
ful conclusions about properties of the spectral function based solely on the input data. Apart from
ensuring consistency with the input data (within their error bars) and the a priori and a posteriori
(conditional) constraints, it is crucial to reliably characterize the accuracy—or even ambiguity—of
the output. We explain how these challenges can be met with two approaches: stochastic optimiza-
tion with consistent constraints and the modified maximum entropy method. We perform illustrative
tests for spectra with a double-peak structure, where we critically examine which spectral properties
are accessible and which ones are lost. For an important practical example, we apply our protocol
to the Fermi polaron problem.
PACS numbers: 02.70.-c, 71.15.Dx, 71.28.+d, 71.10.Fd
I. INTRODUCTION The NAC problem is often characterized as ill-posed.
Mathematically, the near-degeneracy of the kernel im-
Numerous problems in science, from spectral analysis plies two closely related circumstances: (i) the absence of
to image processing, require that we restore properties of the resolvent, and (ii) a continuum of solutions satisfying
a function A(z) from a set of integrals the input data within their error bars (even when inte-
Z ∞ grals over z are replaced with finite sums containing less
or equal to N terms). Nowadays, the first circumstance
gn = G[n, A] ≡ dzK(n, z)A(z), n = 1, . . . , N, (1)
−∞ is merely a minor technical problem, as there exists a
number of methods allowing one to find solutions to (1)
where K(n, z) is a known kernel and {gn } is a finite set of without compromising the error bars of gn .
experimental or numerical input data with error bars. An The second circumstance—the ambiguity of the
important class of such problems—known as numerical solution—is a more essential problem. It is clear that
analytic continuation (NAC)—deals with “pathological” if one formulates the goal as to restore A(z) as a con-
kernels featuring numerous eigenfunctions with anoma- tinuous curve, or to determine its value on a given grid
lously small eigenvalues. An archetypal NAC problem of points, then the goal cannot be reached as stated, ir-
is the numerical spectral analysis at zero temperature, respective of the properties of the kernel K(n, z). The
where the challenge is to restore the non-negative spec- input data set is finite and noisy, thereby introducing a
tral function A(z ≥ 0) satisfying the equation natural limit on the resolution of fine structures in A(z).
Z ∞ Fortunately, the above-formulated goal has little to do
gn = dze−zτn A(z), (2) with the practical world. In an experiment, all devices
0
are characterized by a finite resolution function and the
from numerical data for gn = g(τn ≥ 0). data they collect always correspond to integrals. The
, 2
R
data are processed by making certain assumptions about dzA(z) is typically known with an accuracy that is or-
the underlying function. This motivates an alternative ders of magnitude better than what would be predicted
formulation of the NAC goal involving integrals of A(z) by the central limit theorem if this integral is represented
that render the problem well-defined. With additional by a finite sum of integrals over nonoverlapping inter-
assumptions about the smoothness and other properties vals. Second, the errors are not necessarily distributed
of A(z) behind these integrals, consistent with both a as a Gaussian. Atypical fluctuations can have a signifi-
priori and a posteriori knowledge, the ambiguity of the cant probability and their analysis should not be avoided
solution can be substantially suppressed. The simplest as the actual physical solution may well be one of them.
setup is as follows: To this end, it is important to explore the minimal and
maximal values that the integral im can take, and check
Given a set of finite intervals {∆m }, determine the inte- that these are not significantly different from the typical
grals of the spectral function over these intervals: value of im . In certain cases this criterion cannot be met
Z without increasing the intervals ∆m to an extent when
im = ∆−1 dzA(z), m = 1, . . . , M, (3) the assumption of linearity of A(z) becomes uncontrolled,
m
z∈∆m implying that an important piece of information about
the shape of A(z) in this interval is missing. A charac-
along with the corresponding dispersions of fluctuations teristic example that plays a key role in the subsequent
{σm } (straightforwardly extendable to the dispersion discussion is presented in Fig. 1, where the challenge is
correlation matrix {σmm′ }). to extract the shape of the second peak.
Naively one might think that nothing is achieved by going
from the integrals in (1) to the integrals in (3) because A(z)
2.5
the latter have exactly the same form with the kernel
K̄(m, z) = ∆−1m for z ∈ ∆m and zero otherwise (other
forms of the “resolution function” K̄(m, z) are discussed 2.0
in Sec. II):
Z ∞
1.5
im = I[m, A] = dz K̄(m, z)A(z), m = 1, . . . , M.
−∞
1.0
(4)
This impression, however, is false because kernel proper- 0.5
ties are at the heart of the problem. If for appropriately
small intervals (sufficiently small to resolve the variations
0.0
of A(z)), the uncertainties for im remain small, then one 0 1 2 3 Z
can draw reliable conclusions for the underlying behav-
ior of A(z) itself. The difference between “good” (e.g. FIG. 1. Challenging example of a spectrum A(z) for the NAC
as in Fourier transforms) and pathological kernels is that problem (2). As shown in the main part of the text, the signif-
for the latter, due to the notorious saw-tooth instability, icant width of the first peak makes it essentially impossible to
the uncertainties for im quickly become too large for a controllably restore the width of the second peak, even with
meaningful analysis of fine structures in A(z). small relative error bars (∼ 10−5 ) on gn . On the other hand,
To obtain a solution from the integrals (3), one has the first two moments of the second peak, characterizing its
to invoke the notion of conditional knowledge. The most weight and position, can be extracted reliably.
straightforward approach is to set the spectral function
values at the middle points zm of the intervals ∆m to In the more sophisticated approach used in this work,
Afin (zm ) = im . This is only possible if the intervals can the values of Afin (zm ) can be further optimized (without
be made appropriately narrow without losing accuracy compromising the accuracy of the solution with respect
for the integrals. With this approach we assume that to gn ) to produce a smooth curve. This protocol has
the function is nearly linear over the intervals in ques- the additional advantage of eliminating minima, maxima,
tion. This is a typical procedure for experimental data. gaps, and sharp features that are not guaranteed to exist
Quantifying the error bar on Afin (zm ) necessarily involves by the quality of input data. The nature of the problem is
two numbers: the “vertical” dispersion σm is directly in- such that very narrow peaks (or gaps) with tiny spectral
herited from im , and the “horizontal” error bar ∆m /2 weight can always be imagined to be present (for narrow
represents the interval half-width. intervals they will certainly emerge due to the saw-tooth
The reader should be aware of two issues regard- instability). Our philosophy with respect to these fea-
ing such error bars. First, the error bars for different tures is to erase them within the established error bounds
points are not independent but contain significant multi- and obtain a solution that is insensitive to the interval
point correlations. For example, an unrestricted integral parameters.
