Analytical and Numerical Investigations of Probabilistic Monochromatic Problem

Stefan Schmuck^1*, Jakob Svensson²

¹Istituto di Fisica del Plasma, CNR, Milan, Italy.
²Max-Planck-Institut für Plasmaphysik Teilinstitut Greifswald, Greifswald, Germany.

DOI: 10.4236/jamp.2019.74054                

Abstract

A probabilistic formalism, relying on Bayes’ theorem and linear Gaussian inversion, is adapted, so that a monochromatic problem can be investigated. The formalism enables an objective test in probabilistic terms of the quantities and model concepts involved in the problem at hand. With this formalism, an amplitude (linear parameter), a frequency (non-linear parameter) and a hyperparameter of the Gaussian amplitude prior are inferred jointly given simulated data sets with Gaussian noise contributions. For the amplitude, an analytical normal posterior follows which is conditional on the frequency and the hyperparameter. The remaining posterior estimates the frequency with an uncertainty of MHz, while the convolution of a standard approach would achieve an uncertainty of some GHz. This improvement in the estimation is investigated analytically and numerically, revealing for instance the positive effect of a high signal-to-noise ratio and/or a large number of data points. As a fixed choice of the hyperparameter imposes certain results on the amplitude and frequency, this parameter is estimated and, thus, tested for plausibility as well. From abstract point of view, the model posterior is investigated as well.

Keywords

Bayes’ Theorem, Linear Inversion, Uncertainty, Monochromatic, Spectral Limitation

1. Introduction

Fourier transformation tools are used to obtain information about spectra for a given data set. As any data has an uncertainty, Fourier transformation techniques can be supported by probabilistic theory captured by Bayes’ theorem [1] to improve scientific results and conclusions. The works [2] [3] demonstrate some advantages of the probabilistic ansatz over a conventional approach.

A common and misleading assumption in the field is that the Nyquist theorem determines the spectral band limitation. However, as the Nyquist frequency follows from finite data sampling, it can only be an upper limit but cannot be an estimator for the limit caused actually by a source and/or diagnostic throughput. A similar reasoning applies to the lower spectral limit. After a basic formalism has been derived, it could be shown in Ref. [4] that the band limitation can be well inferred from experimental data. Thereby, a linear Gaussian inversion technique was used to infer spectral amplitudes. Furthermore, a settings posterior was introduced which estimates non-linear parameters and hyperparameters of a problem. For instance, the spectral band limits and the uncertainty on the derived quasi-continuous spectrum, originating in non-probed Fourier coefficients, have been inferred jointly.

For many applied analysis schemes, certain model assumptions and their implications are not tested but assumed to be valid with infinite trust. From scientific point of view, any analysis scheme should be tested given simulated noisy data for which all model assumptions are clear. Then, analysis results and model assumptions can be investigated which is achieved objectively by a probabilistic ansatz. If this can be carried out analytically, valuable information is available when only actual measured data are given for a scientific problem.

The problem of a monochromatic source is a good example to show the powerfulness of Bayesian inference and to test model assumptions. In this work, the formalism derived in Ref. [4] is adapted to a monochromatic problem in Section 2 and applied to simulated data with different noise levels in Section 3. In addition, analytical and numerical investigations are carried out, so that the probabilistic findings can be understood in more detail. For example, a better signal-to-noise ratio improves the estimation of the frequency, while more data points compensate a worsened signal-to-noise ratio. The conclusion section can be found at the end.

2. Adapted Probabilistic Formalism

The formalism derived in Ref. [4] is adapted below to a monochromatic problem, involving one frequency parameter $f_1$ and one amplitude parameter $S_1$. In addition, a hyperparameter is at hand, entering in the normal prior for $S_1$. Abstractly, Bayes’ theorem reads then

$p\left(S_1\mathrm{,}{f}_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)=\frac{p\left(D\mathrm{<mo>\; |\; </mo>}{S}_1\mathrm{,}{f}_1\right)}{p\left(D\right)}p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr}\right)p\left({\sigma }_{S_1\mathrm{,}Pr}\right)p\left(f_1\right)$ (1)

with the joint posterior on the left-hand, the likelihood $p\left(D\mathrm{<mo>\; |\; </mo>}{S}_1\mathrm{,}{f}_1\right)$ devided by the evidence, the amplitude prior $p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr}\right),$ the prior for the hyperparameter, and the prior for the frequency $p\left(f_1\right)$.

The amplitude $S_1$ maps linearly to the data domain via the vector $M_1=M\left(f_1\right)$ dependent on the frequency in a non-linear manner. The data $D$ with $N_D$ entries are assumed to be acquired independently, and the measurement uncertainty of each data point follows a normal distribution with standard deviation ${\sigma }_D$. With these assumptions, the Gaussian likelihood becomes

$p\left(D\mathrm{<mo>\; |\; </mo>}{S}_1\mathrm{,}{f}_1\right)=\frac{\exp \left[-\frac{1}{2}{\left(D-M_1{S}_1\right)}^{\text{T}}{\Sigma }_D^{-1}\left(D-M_1{S}_1\right)\right]}{{\left(2\text{π}\right)}^{N_D/2}{\left|{\Sigma }_D\right|}^{1/2}}$ (2)

with the covariance matrix ${\Sigma }_D={\sigma }_D^2{\delta }_{ij}$. The amplitude prior takes the form

$p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr}\right)=\frac{\exp \left[-\frac{1}{2}{S}_1^2/{\sigma }_{S_1\mathrm{,}Pr}^2\right]}{{\left(2\text{π}\right)}^{1/2}{\sigma }_{S_1\mathrm{,}Pr}}$ (3)

for vanishing prior mean and variance ${\sigma }_{S_1\mathrm{,}Pr}^2$. After some analytical operations, one obtaines the Gaussian amplitude posterior

$p\left(S_1\mathrm{<mo>\; |\; </mo>}{f}_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{,}D\right)=\frac{\exp \left[-\frac{1}{2}\left(S_1-S_{\mathrm{1,}Po}\right){\sigma }_{S_1\mathrm{,}Po}^{-2}\left(S_1-S_{\mathrm{1,}Po}\right)\right]}{{\left(2\text{π}\right)}^{1/2}{\sigma }_{S_1\mathrm{,}Po}}$ (4)

conditional on $f_1$ and ${\sigma }_{S_1\mathrm{,}Pr}$. In the above equation, the posterior variance is given by ${\sigma }_{S_1\mathrm{,}Po}^2={\left(M_1^{\text{T}}{M}_1{\sigma }_D^{-2}+{\sigma }_{S_1\mathrm{,}Pr}^{-2}\right)}^{-1}$ , and the posterior mean reads $S_{\mathrm{1,}Po}={\sigma }_{S_1\mathrm{,}Po}^2{M}_1^{\text{T}}D/{\sigma }_D^2$.

From the remaining terms, the so-called settings posterior

$p\left(f_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)=\frac{1}{K^{\prime }}\frac{{\sigma }_{S_1\mathrm{,}Po}}{{\sigma }_{S_1\mathrm{,}Pr}}\exp \left[\frac{1}{2}{S}_{\mathrm{1,}Po}^2/{\sigma }_{S_1\mathrm{,}Po}^2\right]p\left({\sigma }_{S_1\mathrm{,}Pr}\right)p\left(f_1\right)$ (5)

is derived which can be interpreted as Ockham’s razor with respect to $S_1$. This settings posterior has no general analytical solution, and, thus, the normalisation constant $K^{\prime }$ cannot be stated further. For the remainder of the paper, $p\left({\sigma }_{S_1\mathrm{,}Pr}\right)=1/{\Delta }_{\sigma }$ and $p\left(f_1\right)=1/{\Delta }_f$ are taken as uniform distributions. By this choice, $K^{\prime }=Kp\left({\sigma }_{S_1\mathrm{,}Pr}\right)p\left(f_1\right)$ can be used to investigate analytically and numerically the quantity

$p_{\times K}=\frac{{\sigma }_{S_1\mathrm{,}Po}}{{\sigma }_{S_1\mathrm{,}Pr}}\exp \left[\frac{1}{2}{S}_{\mathrm{1,}Po}^2/{\sigma }_{S_1\mathrm{,}Po}^2\right]\mathrm{.}$ (6)

Finally, the evidence can be identified as

$p\left(D\right)=\frac{K}{{\Delta }_{\sigma }{\Delta }_f}\frac{\exp \left[-\frac{1}{2}{D}^{\text{T}}{\Sigma }_D^{-1}D\right]}{{\left(2\text{π}\right)}^{N_D/2}{\left|{\Sigma }_D\right|}^{1/2}}\mathrm{.}$ (7)

The constant $K/\left({\Delta }_{\sigma }{\Delta }_f\right)$ varies with the chosen model H, including likelihood, priors, and the prior domains for $f_1$ and ${\sigma }_{S_1\mathrm{,}Pr}$. Hence, this constant is linked to the model posterior $p\left(H\mathrm{<mo>\; |\; </mo>}D\right)$ which is of importance, when the model is further abstracted, or a comparison with an alternative model is carried out.

3. Monochromatic Problem

For a monochromatic even source, data sets with different noise levels are simulated, and the formalism derived in the previous section is applied. To explain the results found, analytical and numerical investigations are carried out. While this problem is treated abstractly in the following, the main results will be presented for two examples with low and high noise contribution to ease the presentation.

(A) Simulated Data: Two Examples

Data sets are modelled for the real-world interferometer found in Ref. [5] which achieves an optical path difference to obtain an interferogram. Formally, the interferogram $V_i=2S\cos \left(2\text{π}f{x}_i/c\right)$ is simulated with the amplitude $S=0.5$ and the frequency f = 50 GHz. The optical path difference grid $x_i$ is chosen as $x_i=\left(i-129\right)\Delta x$ with $i=1,\cdots ,N_D=789$ with constant increment $\Delta x=40\mu \text{m}$. Accordingly, the total length L = 31.56 mm and the length $L_{SS}=21.28\text{mm}$ of the single-sided region follow. Two noisy data sets $D_i=V_i+{\mathcal{N}}_i\left(\mathrm{0,}{\sigma }_D^2\right)$ are generated with the noise levels ${\sigma }_D$ at 0.05 and 1.0 (see Figure 1).

(B) Application of Formalism

The linear mapping of the amplitude parameter is given by $M_{\mathrm{1,}i}=2\cos \left(2\pi f_1{x}_i/c\right)$ , yielding $M_1$. For the two simulated data sets, Figure 2 shows $p\left(S_1\mathrm{<mo>\; |\; </mo>}{f}_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{,}D\right)$ , varying in $f_1$ and for three values of ${\sigma }_{S_1\mathrm{,}Pr}$. Even for the high noise level case, $S_{\mathrm{1,}Po}$ peaks close to 50 GHz, and the uncertainty ${\sigma }_{S_1\mathrm{,}Po}$ remains reasonably small. The peak has the full width at half maximum (FWHM) of about 8 GHz which is roughly twice as large as the classical spectral resolution stated by $c/\left(2L\right)\approx 4.75\text{GHz}$. High and low values of ${\sigma }_{S_1\mathrm{,}Pr}$ seem to have little effect on the results.

Figure 1. Simulated interferograms $V$ (without noise, black) and $D$ (with noise, red). Formally, $V=2S\cos \left(2\pi fx/c\right)$ was used with amplitude $2S=1$ , frequency f = 50 GHz, and optical path difference grid x with increment $\Delta x=40\mu \text{m}$ (number of data points $N_D=789$ ). The data-like quantity $D=V+\mathcal{N}\left(\mathrm{0,}{\sigma }_D^2\right)$ resembles a sample of $V$ , assuming a Gaussian noise component with variance ${\sigma }_D^2$. Furthermore, the noise is chosen to be independent of x. (a) ${\sigma }_D=0.05$ (low noise); (b) ${\sigma }_D=1.0$ (high noise).

Figure 2. Conditional Gaussian amplitude posterior $p\left(S_1\mathrm{<mo>\; |\; </mo>}{f}_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{,}D\right)$ with mean $S_{\mathrm{1,}Po}$ and standard deviation ${\sigma }_{S_1\mathrm{,}Po}$. This posterior depends on the frequency $f_1$ , prior standard deviation ${\sigma }_{S_1\mathrm{,}Pr}$ for $S_1$ , and noisy interferogram data $D$. The conditional posterior is obtained for two sets $D$ (see Figure 1) modelled by a monochromatic source with the frequency f = 50 GHz, the amplitude S = 0.5 and the noise levels ${\sigma }_D=0.05$ (left column) and 1 (right column). The posterior mean (cyan) peaks at values of about 0.5 close to 50 GHz. The spectral FWHM of this peak accounts for about 8 GHz. The posterior standard deviation (dashed-white for its multiples) increases with the noise level, as one would expect. Only small values of ${\sigma }_{S_1\mathrm{,}Pr}$ affect $S_{\mathrm{1,}Po}$ when the noise level is high.

The settings posterior $p\left(f_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ is proportional to the quantity $p_{\times K}$ which takes very large values even on logarithmic scale and, thus, $\ln p_{\times K}$ is used in the following. For the two cases, $\ln p_{\times K}$ is shown in Figure 3(a) and Figure 3(b) versus the spectral parameter and for three prior values of ${\sigma }_{S_1\mathrm{,}Pr}$. All three cross-sections look similar, but close to 50 GHz the peak is highest for values near ${\sigma }_{S_1\mathrm{,}Pr}=0.5$. Identifying this global maximum by $\max \left(p_{\times K}\right)$ , the ratio $p_{\times K}/\max \left(p_{\times K}\right)$ takes reasonable numbers and is presented for the

Figure 3. $p_{\times K}$ , being proportional to joint posterior $p\left(f_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ , dependent on frequency $f_1$ and prior standard deviation ${\sigma }_{S_1\mathrm{,}Pr}$ given interferogram data set $D$ (see Figure 1) with low and high noise levels ${\sigma }_D=0.05$ (left column) and 1 (right column), respectively. Since $p_{\times K}$ takes large values, $\ln p_{\times K}$ is presented in (a) and (b). Deviding by the global maximum, the ratio $p_{\times K}/\max \left(p_{\times K}\right)$ becomes manageable from numerical point of view as shown in (c) and (d). This ratio peaks in the vicinity of 50 GHz and ${\sigma }_{S_1\mathrm{,}Pr}=0.5$ ( $\log {\sigma }_{S_1\mathrm{,}Pr}=-0.3$ ). The spectral width of $p_{\times K}/\max \left(p_{\times K}\right)$ increases with the noise level from about 10 to 200 MHz but remains well below the classical spectral resolution of about 4.75 GHz. The width with respect to ${\sigma }_{S_1\mathrm{,}Pr}$ remains constant and covers one order of magnitude. (a) $\ln p_{\times K}$ (low noise case) for broad spectral domain at three values of ${\sigma }_{S_1\mathrm{,}Pr}$ ; (b) $\ln p_{\times K}$ (high noise case) for broad spectral domain at three values of ${\sigma }_{S_1\mathrm{,}Pr}$ ; (c) $p_{\times K}/\max \left(p_{\times K}\right)$ (low noise case) for spectral domain 49.9 - 50.1 GHz; (d) $p_{\times K}/\max \left(p_{\times K}\right)$ (high noise case) for spectral domain 49 - 51 GHz.

relevant domains of $f_1$ and ${\sigma }_{S_1\mathrm{,}Pr}$ by Figure 3(c) and Figure 3(d) for both data sets. This ratio is narrow in the vicinity of 50 GHz and in the domain $0.3\le {\sigma }_{S_1\mathrm{,}Pr}\le 3$. Furthermore, the peaking has a spectral width of about 10 MHz and 200 MHz for the low and high noise case, respectively. When interpreted as spectral resolution, this width is smaller by several orders of magnitude with respect to its classical counterpart of about 4.75 GHz.

(C) Analytical and Numerical Investigations

To have further inside in the results, the conditional amplitude posterior and $p_{\times K}$ are analytically investigated. In order to do so, the noise contribution in the data is neglected, meaning that $D_i$ is replaced with the noise-less interferogram $V_i=2S\cos \left(2\text{π}f{x}_i/c\right)$. This allows the tracking of the dependencies of the posteriors on the original amplitude S and frequency f, noise level ${\sigma }_D$ , and number of data points $N_D$. Furthermore, a second factorisation $p\left(f_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)=p\left(f_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{,}D\right)p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ can be justified. Parallel to this investigation, numerics and conditional posteriors are presented.

1) Conditional Amplitude Posterior

Starting point is the variance ${\sigma }_{S_1\mathrm{,}Po}^2$ of the conditional posterior $p\left(S_1\mathrm{<mo>\; |\; </mo>}{f}_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{,}D\right)$ (see Equation (4)). Formally, one finds for the inverse

${\sigma }_{S_1\mathrm{,}Po}^{-2}=\frac{M_1^{\text{T}}{M}_1}{{\sigma }_D^2}+\frac{1}{{\sigma }_{S_1\mathrm{,}Pr}^2}\mathrm{.}$ (8)

As long as the spatial sampling is well enough, one can use the approximation

$\begin{array}{c}{M}_1^{\text{T}}{M}_1=4{\displaystyle \underset{i=1}{\overset{N_D}{\sum }}}\cos \left(2\text{π}\frac{f_1}{c}{x}_i\right)\cos \left(2\text{π}\frac{f_1}{c}{x}_i\right)\\ \approx 2N_D\left[1+\cos \left(2\text{π}\frac{f_1}{c}{L}_{SS}\right)\text{sinc}\left(2\text{π}\frac{f_1}{c}L\right)\right],\end{array}$ (9)

as shown by Equation (A1), to get the posterior variance

${\sigma }_{S_1\mathrm{,}Po}^2=\frac{{\sigma }_D^2/\left(2N_D\right)}{1+\cos \left(2\text{π}\frac{f_1}{c}{L}_{SS}\right)\text{sinc}\left(2\text{π}\frac{f_1}{c}L\right)+\frac{{\sigma }_D^2}{2N_D{\sigma }_{S_1\mathrm{,}Pr}^2}}\mathrm{.}$ (10)

This variance increases quadratically with the noise level and reduces with the number of data samples. Thus, more data points per spatial unit can compensate the noise contribution, at least partly. However, this assumes that the noise is independent of the data sampling.

The modulated sinc function in the denominator of Equation (10) depends on $f_1$ , L and $L_{SS}$ for which the interferogram is available. The sinc function becomes unity and vanishes for small and large frequencies, respectively. For instance, since both data examples share the same spatial domain, the sinc function is close to 0 above about 20 GHz. The amplitude prior influences ${\sigma }_{S_1\mathrm{,}Po}^2$ , when the ratio ${\sigma }_D^2/\left(2N_D{\sigma }_{S_1\mathrm{,}Pr}^2\right)$ becomes larger than 1 which is obtained for a large noise level, a small value of ${\sigma }_{S_1\mathrm{,}Pr}$ and/or few data points.

For the posterior mean,

$S_{\mathrm{1,}Po}\left(f_1\right)=\frac{{\sigma }_{S_1\mathrm{,}Po}^2}{{\sigma }_D^2}{M}_1^{\text{T}}V=\frac{2{\displaystyle \underset{i\mathrm{<mo>\; =\; </mo>1}}{\overset{N_D}{\sum }}}{V}_i\cos \left(2\text{π}\frac{f_1}{c}{x}_i\right)/\left(2N_D\right)}{1+\cos \left(2\text{π}\frac{f_1}{c}{L}_{SS}\right)\text{sinc}\left(2\text{π}\frac{f_1}{c}L\right)+\frac{{\sigma }_D^2}{2N_D{\sigma }_{S_1\mathrm{,}Pr}^2}}$ (11)

follows, leaving aside the noise contribution but the uncertainty on the mean is still captured by ${\sigma }_{S_1\mathrm{,}Po}$. With the approximation (A1) and the trigonometric identity (B1), one can resolve

$S_{\mathrm{1,}Po}\left(f_1\right)=S\frac{\cos \left(\text{π}\frac{f-f_1}{c}{L}_{SS}\right)\text{sinc}\left(\text{π}\frac{f-f_1}{c}L\right)+\cos \left(\text{π}\frac{f+f_1}{c}{L}_{SS}\right)\text{sinc}\left(\text{π}\frac{f+f_1}{c}L\right)}{1+\cos \left(2\text{π}\frac{f_1}{c}{L}_{SS}\right)\text{sinc}\left(2\text{π}\frac{f_1}{c}L\right)+\frac{{\sigma }_D^2}{2N_D{\sigma }_{S_1\mathrm{,}Pr}^2}}\mathrm{.}$ (12)

Indeed, one finds the original amplitude $S_{1,Po}=S$ for $f_1=f$ and prior values with the condition ${\sigma }_D^2/\left(2N_D\right)\ll {\sigma }_{S_1\mathrm{,}Pr}^2$. Furthermore, due to $\cos \left(\text{π}\frac{f-f_1}{c}L\right)\text{sinc}\left(\text{π}\frac{f-f_1}{c}{L}_{SS}\right)$ there is a peaking of $S_{\mathrm{1,}Po}\left(f_1\right)$ in the vicinity of $f_1\approx f$ with the FWHM determined by both trigonometric arguments. Away from f, the oscillating amplitude drops like $1/\left(f-f_1\right)$ due to the sinc function which can be seen in Figure 2.

2) Conditional Posterior for f₁

Now the information is available to investigate $p\left(f_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$. Neglecting the noise, the triple product in the exponent of the exponential in $p_{\times K}$ reads

$\begin{array}{l}\frac{1}{2}{S}_{\mathrm{1,}Po}{\sigma }_{S_1\mathrm{,}Po}^{-2}{S}_{\mathrm{1,}Po}\approx \frac{1}{2}\frac{M_1^{\text{T}}V}{{\sigma }_D^2}{\sigma }_{S_1\mathrm{,}Po}^2\frac{M_1^{\text{T}}V}{{\sigma }_D^2}\\ =\frac{1}{2}{\sigma }_{S_1\mathrm{,}Po}^2\frac{4N_D^2{S}^2}{{\sigma }_D^4}[\cos \left(\text{π}\frac{f-f_1}{c}{L}_{SS}\right)\text{sinc}\left(\text{π}\frac{f-f_1}{c}L\right)\\ +{\cos \left(\text{π}\frac{f+f_1}{c}{L}_{SS}\right)\text{sinc}\left(\text{π}\frac{f+f_1}{c}L\right)]}^2\mathrm{.}\end{array}$ (13)

The global maximum is close to $f_1=f$ , and, assuming a sufficiently large frequency f, one gets the approximation

$\begin{array}{l}{\left[\cos \left(\text{π}\frac{f-f_1}{c}{L}_{SS}\right)\text{sinc}\left(\text{π}\frac{f-f_1}{c}L\right)\right]}^2\\ \approx \left[1-{\left(\text{π}\frac{f-f_1}{c}{L}_{SS}\right)}^2-\frac{1}{3}{\left(\text{π}\frac{f-f_1}{c}L\right)}^2\right]\end{array}$ (14)

when the Taylor series expansions $\text{sinc}\left(z\right)\approx 1-z^2/\mathrm{3\; <mo>\; !\; </mo>}$ and $\cos \left(z\right)\approx 1-z^2/2$ are used. Then,

$\frac{1}{2}{S}_{\mathrm{1,}Po}{\sigma }_{S_1\mathrm{,}Po}^{-2}{S}_{\mathrm{1,}Po}\approx \frac{1}{2}{\sigma }_{S_1\mathrm{,}Po}^2\frac{4N_D^2{S}^2}{{\sigma }_D^4}\left[1-{\left(f-f_1\right)}^2{\left(\frac{\text{π}}{c}\right)}^2\left(L_{SS}^2+\frac{1}{3}{L}^2\right)\right]$ (15)

remains, and the term which is independent on $f-f_1$ takes very large numerical values and is treated in more detail in the next subsection. The remaining term can be rewritten by a quadratic exponent of a Gaussian with posterior mean $f_{\mathrm{1,}Po}=f$ and square root of the variance

$\begin{array}{c}{\sigma }_{f_1\mathrm{,}Po}=\frac{1}{{\sigma }_{S_1\mathrm{,}Po}}\frac{{\sigma }_D^2}{2N_{D}S}\frac{1}{\text{π}}\frac{c}{L}{\left(\frac{3}{3L_{SS}^2/L^2+1}\right)}^{1/2}\\ \approx \frac{{\sigma }_D}{S{N}_D^{1/2}}\frac{1}{\text{π}}\frac{c}{L}{\left(\frac{3}{2}\frac{1}{3L_{SS}^2/L^2+1}\right)}^{1/2}\mathrm{,}\end{array}$ (16)

using Equation (10) for large frequencies and ${\sigma }_{S_1\mathrm{,}Pr}^2\gg {\sigma }_D^2/\left(2N_D\right)$. This uncertainty increases with the noise but decreases with the square root of the number of data points per spatial unit, the signal level and the spatial domain covered. For the data set examples, one inserts $L=31.56$ mm, $L_{SS}=21.28$ mm, $N_D=789$ , ${\sigma }_D/S=0.1$ and 2 to find the astonishing numbers ${\sigma }_{f_1\mathrm{,}Po}=8.6$ and 171.9 MHz, respectively. This explains the narrowness of $p\left(f_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ with respect to $f_1$ as shown in Figure 3(c) and Figure 3(d).

In fact, ${\sigma }_{S_1\mathrm{,}Po}$ changes with $f_1$ on a GHz, but the peaking in $f_1$ is of the order of MHz. Hence, the factorisation $p\left(f_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)\approx p\left(f_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{,}D\right)p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ with the conditional Gaussian posterior

$p\left(f_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr}\mathrm{,}D\right)=\frac{\exp \left[-\frac{1}{2}\frac{{\left(f_1-f_{\mathrm{1,}Po}\right)}^2}{{\sigma }_{f_1\mathrm{,}Po}^2}\right]}{{\left(2\text{π}\right)}^{1/2}{\sigma }_{f_1\mathrm{,}Po}}$ (17)

and

$p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)=\frac{{\sigma }_{S_1\mathrm{,}Po}}{{\sigma }_{S_1\mathrm{,}Pr}}\frac{\exp \left[\frac{1}{2}{\sigma }_{S_1\mathrm{,}Po}^2\frac{4N_D^2{S}^2}{{\sigma }_D^4}\right]}{K}{\left(2\text{π}\right)}^{1/2}{\sigma }_{f_1\mathrm{,}Po}$ (18)

is justified. Thereby, ${\sigma }_{S_1\mathrm{,}Po}$ needs to be taken at $f_1=f_{\mathrm{1,}Po}$. The remaining probability $p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ is investigated in the next subsection.

Numerics

Because the variance ${\sigma }_{S_1\mathrm{,}Po}^2$ changes in $f_1$ on GHz scale but

$\exp \left[\frac{1}{2}{S}_{\mathrm{1,}Po}{\sigma }_{S_1\mathrm{,}Po}^{-2}{S}_{\mathrm{1,}Po}\right]$ on MHz scale, only the exponential needs to be taken into account for the determination of $f_{\mathrm{1,}Po}$ and ${\sigma }_{f_1\mathrm{,}Po}$. Furthermore, one faces large numbers, and, hence, the exponent is dealt with directly.

At a given ${\sigma }_{S_1\mathrm{,}Pr}$ , the frequency at which the maximum of the exponent occurs is identified as the posterior mean $f_{\mathrm{1,}Po}\left({\sigma }_{S_1\mathrm{,}Pr}\right)$. The posterior standard deviation ${\sigma }_{f_1\mathrm{,}Po}\left({\sigma }_{S_1\mathrm{,}Pr}\right)$ is obtained when the ratio of the exponent to its maximum at $f_{\mathrm{1,}Po}$ reads $\frac{1}{2}{S}_{\mathrm{1,}Po}{\sigma }_{S_1\mathrm{,}Po}^{-2}{S}_{\mathrm{1,}Po}/\max \left(\frac{1}{2}{S}_{\mathrm{1,}Po}{\sigma }_{S_1\mathrm{,}Po}^{-2}{S}_{\mathrm{1,}Po}\right)=-1/2$. For both data sets, the ratios are shown in Figure 4(a) and Figure 4(b) in the vicinity of the peak. While the uncertainty of $f_{\mathrm{1,}Po}$ and, hence, the width of the Gaussian augment with the noise level (see also Figure 5(a) and Figure 5(b)) as described by Equation (16), f remains included in the 2- ${\sigma }_{f_1\mathrm{,}Po}$ band centred at $f_{\mathrm{1,}Po}$. Furthermore, the choice of ${\sigma }_{S_1\mathrm{,}Pr}$ affects $f_{\mathrm{1,}Po}$ and ${\sigma }_{f_1\mathrm{,}Po}$ only for small values (see Figure 5(a) and Figure 5(b)).

3) Posterior for ${\sigma }_{1,Pr}$

After the factorisation, one combines Equations (16) and (18) like

$\begin{array}{c}p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)=\frac{P^{\prime }}{K}=\frac{{\sigma }_{S_1\mathrm{,}Po}}{{\sigma }_{S_1\mathrm{,}Pr}}\frac{\exp \left[\frac{1}{2}{\sigma }_{S_1\mathrm{,}Po}^2\frac{4N_D^2{S}^2}{{\sigma }_D^4}\right]}{K}{\left(2\text{π}\right)}^{1/2}{\sigma }_{f_1\mathrm{,}Po}\\ =\frac{P}{K}\frac{{\sigma }_D^2}{N_{D}S}\frac{c}{L}\sqrt{\frac{3}{\text{2π}}\frac{L^2}{3L_{SS}^2+L^2}}\end{array}$ (19)

with the important term

Figure 4. Ratio of exponent to its maximum dependent on frequency $f_1$ and keeping constant prior standard deviation ${\sigma }_{S_1\mathrm{,}Pr}$. For the data set with low noise level ( ${\sigma }_D=0.05$ , left), the peaking is more narrow than the one for the high noise level ( ${\sigma }_D=1$ , right). The maximum location and the width of the ratio are used to estimate posterior mean $f_{\mathrm{1,}Po}\left({\sigma }_{S_1\mathrm{,}Pr}\right)$ and standard deviation ${\sigma }_{f_1\mathrm{,}Po}\left({\sigma }_{S_1\mathrm{,}Pr}\right)$. (a) Ratio for low noise data set at three values of ${\sigma }_{S_1\mathrm{,}Pr}$. For ${\sigma }_{S_1\mathrm{,}Pr}=0.5$ , $f_{1,Po}=50.016$ GHz and ${\sigma }_{f_1\mathrm{,}Po}=8.5$ MHz are extracted which are consistent with 50 GHz used to simulate the data; (b) Ratio for high noise data set at three values of ${\sigma }_{S_1\mathrm{,}Pr}$. For ${\sigma }_{S_1\mathrm{,}Pr}=0.5$ , $f_{\mathrm{1,}Po}=50.167$ GHz and ${\sigma }_{f_1\mathrm{,}Po}=173$ MHz are estimated. ${\sigma }_{f_1\mathrm{,}Po}$ and the deviation of $f_{\mathrm{1,}Po}$ from 50 GHz increase with the noise level. Furthermore, small prior values for ${\sigma }_{S_1\mathrm{,}Pr}$ affect ${\sigma }_{f_1\mathrm{,}Po}$.

Figure 5. Posterior mean $f_{\mathrm{1,}Po}$ and standard deviation ${\sigma }_{f_1\mathrm{,}Po}$ for frequency $f_1$ dependent on prior standard deviation ${\sigma }_{S_1\mathrm{,}Pr}$ for noise levels ${\sigma }_D=0.05$ (left) and 1 (right). The original frequency f = 50 GHz (red line) is included usually in the band $f_{\mathrm{1,}Po}\pm 2{\sigma }_{f_{\mathrm{1,}Po}}$. However, ${\sigma }_{f_{\mathrm{1,}Po}}$ increases from about 10 MHz to 200 MHz given both noise levels. The posterior $p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ (bottom row) changes little for the different noise levels and has a global maximum close to S = 0.5 (dashed-black). (a) $f_{\mathrm{1,}Po}\pm 2{\sigma }_{f_{\mathrm{1,}Po}}$ for ${\sigma }_D=0.05$ ; (b) $f_{\mathrm{1,}Po}\pm 2{\sigma }_{f_{\mathrm{1,}Po}}$ for ${\sigma }_D=1$ ; (c) $p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ for ${\sigma }_D=0.05$ ; (d) $p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ for ${\sigma }_D=1$.

$P=\frac{\exp \left[\frac{1}{2}{\sigma }_{S_1\mathrm{,}Po}^2\frac{4N_D^2{S}^2}{{\sigma }_D^4}\right]}{{\sigma }_{S_1\mathrm{,}Pr}}\mathrm{.}$ (20)

The exponent of $P$ becomes

$\frac{{\sigma }_{S_1\mathrm{,}Pr}^2\frac{N_D{S}^2}{{\sigma }_D^2}}{{\sigma }_{S_1\mathrm{,}Pr}^2\left[1+\cos \left(2\text{π}\frac{f_{\mathrm{1,}Po}}{c}{L}_{SS}\right)\text{sinc}\left(2\text{π}\frac{f_{\mathrm{1,}Po}}{c}L\right)\right]+\frac{{\sigma }_D^2}{2N_D}}$ (21)

which is monotonically rising in ${\sigma }_{S_1\mathrm{,}Pr}$ from 0 to about $N_D{S}^2/{\sigma }_D^2$. In case ${\sigma }_{S_1\mathrm{,}Pr}^2>{\sigma }_D^2/\left(2N_D\right)$ and $S/{\sigma }_D=1$ , the number of data points dominates. For the two examples with $N_D=789$ , large numbers $\exp \left(78900\right)$ and $\exp \left(197\right)$ are at hand as could be seen in Figure 3(a) and Figure 3(b).

As the exponent vanishes for ${\sigma }_{S_1\mathrm{,}Pr}=0$ , $1/{\sigma }_{S_1\mathrm{,}Pr}$ and, thus, P and $P^{\prime }$ diverge. Hence, this global maximum should be excluded by a proper choice for the prior of ${\sigma }_{S_1\mathrm{,}Pr}$. Doing so, P has the global maximum $P_{GMax}$ at

${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}^2=\frac{S^2}{2}-\frac{{\sigma }_D^2}{2N_D}+\frac{S^2}{2}{\left(1-2\frac{{\sigma }_D^2}{N_D{S}^2}\right)}^{1/2}$ (22)

supposed a large $f_{\mathrm{1,}Po}$. Interestingly, $N_D$ reduces the noise impact drastically, so that even for elevated noise levels ${\sigma }_D\approx S$ , the maximum is still given by S with minor corrections.

The uniform prior $p\left({\sigma }_{S_1\mathrm{,}Pr}\right)$ may be finite for the domain ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}\le {\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}$ , and the maximum of P in this domain is denoted by $P_{Max}$. Then, the peaking can be made more obvious by

$\frac{P}{P_{Max}}=\frac{{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}Max}}{{\sigma }_{S_1\mathrm{,}Pr}}\exp \left[-\frac{1}{2}\frac{{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}Max}^2-{\sigma }_{S_1\mathrm{,}Pr}^2}{{\sigma }_{S_1\mathrm{,}Pr}^2+\frac{{\sigma }_D^2}{2N_D}}\frac{S^2}{{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}Max}^2+\frac{{\sigma }_D^2}{2N_D}}\right]$ (23)

after some algebra for a large $f_{\mathrm{1,}Po}$. With the normalisation of the above equation

${\displaystyle \underset{{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}}{\overset{{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}}{\int }}}\frac{P}{P_{Max}}\text{d}{\sigma }_{S_1\mathrm{,}Pr}=K_{{\sigma }_{S_1\mathrm{,}Pr}}\left({\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}\right)\mathrm{,}$ (24)

the posterior becomes

$p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)=\frac{1}{K_{{\sigma }_{S_1\mathrm{,}Pr}}}\frac{P}{P_{Max}}\mathrm{.}$ (25)

In case, the prior is finite for the domain ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}\le {\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\le {\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}$ , so that $P_{Max}=P_{GMax}$ , one can set ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}Max}^2=S^2-{\sigma }_D^2/\left(2N_D\right)$. This reveals the simple forms

$P_{GMax}=\frac{\exp \left[\frac{N_D{S}^2}{{\sigma }_D^2}-\frac{1}{2}\right]}{\sqrt{S^2-\frac{{\sigma }_D^2}{2N_D}}}$ (26)

and

$\frac{P}{P_{GMax}}=\frac{\sqrt{S^2-\frac{{\sigma }_D^2}{2N_D}}}{{\sigma }_{S_1\mathrm{,}Pr}}\exp \left[\frac{1}{2}\left(1-\frac{S^2}{{\sigma }_{S_1\mathrm{,}Pr}^2+\frac{{\sigma }_D^2}{2N_D}}\right)\right]$ (27)

of Equation (23). With the above expression and Equation (19), the constant

$K={P^{\prime }}_{GMax}{K}_{{\sigma }_{S_1\mathrm{,}Pr}}=\frac{{\sigma }_D^2}{N_{D}S}\frac{c}{L}\sqrt{\frac{3}{2\text{π}}\frac{L^2}{3L_{SS}^2+L^2}}\frac{\exp \left[\frac{N_D{S}^2}{{\sigma }_D^2}-\frac{1}{2}\right]}{\sqrt{S^2-\frac{{\sigma }_D^2}{2N_D}}}{K}_{{\sigma }_{S_1\mathrm{,}Pr}}$ (28)

follows with proper unit Hz. K depends in a complicated way on the noise level, the signal, the number of data points, the spatial domain for which the data is acquired, and the limits ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}$ and ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}$. For the two examples, one gets $\ln \left({P^{\prime }}_{GMax}/\text{Hz}\right)=78910$ ( ${\sigma }_D=0.05$ ) and $\ln \left({P^{\prime }}_{GMax}/\text{Hz}\right)=214$ ( ${\sigma }_D=1$ ).

Numerics

Starting point is Equation (19) which is rewritten like

$\begin{array}{l}\ln p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)+\ln K=\ln P^{\prime }\\ =\ln \left(\frac{{\sigma }_{S_1\mathrm{,}Po}}{{\sigma }_{S_1\mathrm{,}Pr}}\max \left(\exp \left[\frac{1}{2}{S}_{\mathrm{1,}Po}{\sigma }_{S_1\mathrm{,}Po}^{-2}{S}_{\mathrm{1,}Po}\right]\right){\left(2\text{π}\right)}^{1/2}{\sigma }_{f_1\mathrm{,}Po}\right)\end{array}$ (29)

to handle large numerical values. Thereby, ${\sigma }_{S_1\mathrm{,}Po}$ is evaluated at $f_1=f_{1,Po}$. The far right-hand side of the above equation is numerically available, and the maximum $\ln {P^{\prime }}_{Max}$ at ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}Max}$ and the normalisation $K_{{\sigma }_{S_1\mathrm{,}Pr}}$ can be determined for the domain from ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}$ to ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}$.

Choosing for the moment ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}={10}^{-2}$ and ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}={10}^3$ , so that the global maximum ${P^{\prime }}_{GMax}$ is included, the above formalism yields $\ln \left({P^{\prime }}_{GMax}/\text{Hz}\right)=76861$ (213) and $K_{{\sigma }_{S_1\mathrm{,}Pr}}=6.32$ (6.36) for the data set simulated with ${\sigma }_D=0.05$ (1). For both examples, the order of magnitude of $\ln \left({P^{\prime }}_{GMax}/\text{Hz}\right)$ is well described by Equation (26).

For the data sets, both ${P^{\prime }}_{GMax}$ are located at ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}=0.499$ and 0.501 which is in the vicinity of S (see Figure 5(c) and Figure 5(d)) and is in agreement with expression (22). For the different noise levels, $p\left({\sigma }_{S_1\mathrm{,}Pr}\mathrm{<mo>\; |\; </mo>}D\right)$ has a similar shape; the probability reduces drastically below ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}$ but remains finite above up to ${\sigma }_{S_1\mathrm{,}Pr}=10$ , at least.

4) Model Posterior

The posterior of the model with uniform priors for ${\sigma }_{S_1}$ and $f_1$ is given by

$\text{}\text{}\text{}\text{}\text{}\text{}\text{}p\left(H\mathrm{<mo>\; |\; </mo>}D\right)=\frac{{P^{\prime }}_{Max}{K}_{{\sigma }_{S_1\mathrm{,}Pr}}}{\left({\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}-{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}\right)\left(f_{\mathrm{1,}Pr\mathrm{,}U}-f_{\mathrm{1,}Pr\mathrm{,}L}\right)}\mathrm{.}$ (30)

This posterior weighs the constants obtained with respect to the spectral and amplitude domains. Hence, $p\left(H\mathrm{<mo>\; |\; </mo>}D\right)$ increases with the constants in the numerator and smaller domains in the denominator. When one investigates Figure 3(a) and Figure 3(b), it is obvious that $p\left(H\mathrm{<mo>\; |\; </mo>}D\right)$ increases, when $f_{\mathrm{1,}Pr\mathrm{,}L}$ and $f_{\mathrm{1,}Pr\mathrm{,}U}$ are chosen to include $f_{\mathrm{1,}Po}$ and are separated by few GHz or even some 100s MHz. Similarly, choosing ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}$ and ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}$ to include the global maximum ${P^{\prime }}_{GMax}$ at ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}$ (see Figure 5(c) and Figure 5(d)) elevates $p\left(H\mathrm{<mo>\; |\; </mo>}D\right)$. Furthermore, $p\left(H\mathrm{<mo>\; |\; </mo>}D\right)$ has a global maximum when ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}$ and ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}$ approach ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\approx S$ from either side. This is numerically confirmed by scanning ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}$ and ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}$ over several orders of magnitude and evaluating the model posterior for the two data sets (see Figure 6(a) and Figure 6(b)).

(D) Joint Posterior and Marginal Posteriors for S₁ and f₁

Because the model posterior peaks at ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}$ , the marginal posterior for $f_1$ is simply given by the Gaussian $p\left(f_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr}={\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ with posterior mean $f_{\mathrm{1,}Po}$ and standard deviation ${\sigma }_{f_1\mathrm{,}Po}$ (see Section III C2 and Figure 5(a) and Figure 5(b)). Furthermore, the joint posterior

$p\left(S_1\mathrm{,}{f}_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)=p\left(S_1\mathrm{<mo>\; |\; </mo>}{f}_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)p\left(f_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ (31)

can be evaluated and is shown for both data sets in Figure 7(a) and Figure 7(b). This joint posterior peaks close to S and f used to simulate the data, and its width increases with the noise level. However, S and f are usually located inside the peak, and, hence, one can trust the results within the posterior uncertainties.

The marginal posterior for $S_1$ is obtained by the marginalisation

$p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)={\displaystyle \int }p\left(S_1\mathrm{,}{f}_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)\text{d}{f}_1\mathrm{.}$ (32)

The above integration is performed numerically, and the results are shown in Figure 7(c) and Figure 7(d) for both data sets. Another way to achieve the marginalisation is to draw samples for $f_1$ from $p\left(f_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ , and for each of those samples draw samples for $S_1$ from $p\left(S_1\mathrm{<mo>\; |\; </mo>}{f}_1\mathrm{,}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$. Histograms of these posterior samples are given in Figure 7(c) and Figure 7(d).

Figure 6. Model posterior $p\left(H\mathrm{<mo>\; |\; </mo>}D\right)$. The stated model uses a uniform prior for the prior standard deviation ${\sigma }_{S_1\mathrm{,}Pr}$ in the domain from ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}L}$ to ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}U}$. As these limits approach ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\approx 0.499$ for ${\sigma }_D=0.05$ and 0.501 for ${\sigma }_D=1$ (dashed-white), the posterior increases. Thereby, ${\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\approx \sqrt{S^2-{\sigma }_D^2/\left(2N_D\right)}$ is closely linked with the original amplitude S = 0.5. (a) Low noise case ${\sigma }_D=0.05$ ; (b) High noise case ${\sigma }_D=1$.

Figure 7. Joint posterior $p\left(S_1\mathrm{,}{f}_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ and marginal posterior $p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ of frequency $f_1$ and amplitude $S_1$ for low noise ( ${\sigma }_D=0.05$ , left) and high noise data set ( ${\sigma }_D=1$ , right). The joint posterior peaks close to the original values S = 0.5 and f = 50 GHz used to simulate the data sets. For the higher noise level, the posterior width increases, and the location of the global maximum separates further from S and f. By numerical integration, the marginal $p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ is evaluated, and its shape is close to a normal distribution. This is confirmed by a second approach which uses 10⁵ samples from the joint posterior to estimate the marginal posterior mean and standard deviation for $S_1$. The histogram of the samples is given as well. (a) $p\left(S_1\mathrm{,}{f}_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ for low noise level; (b) $p\left(S_1\mathrm{,}{f}_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ for high noise level; (c) $p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ for low noise level. The marginal mean and standard deviation read 0.4998 and 0.0013, respectively; (d) $p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ for high noise level. One finds that the marginal mean and standard deviation account for 0.505 and 0.026, respectively.

Evaluating the mean and standard deviation from these samples approximates $p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ by a Gaussian. In doing so for the data set with ${\sigma }_D=0.05$ (1), the marginal posterior mean $S_{\mathrm{1,}Po\mathrm{,}Ma}=0.4998$ (0.505) and standard deviation 0.0013 (0.026) are obtained. Both normal distributions are good approximations for $p\left(S_1\mathrm{<mo>\; |\; </mo>}{\sigma }_{S_1\mathrm{,}Pr\mathrm{,}GMax}\mathrm{,}D\right)$ (see Figure 7(c) and Figure 7(d)).

4. Conclusions

The presented work investigates a monochromatic problem with an adapted version of Bayes’ theorem to infer an amplitude parameter, a frequency parameter and a hyperparameter which acts on the amplitude prior only. The amplitude at a given frequency is a parameter which has a linear mapping to the data domain. If the measurements have a normal noise contribution (Gaussian likelihood) and the amplitude prior is chosen to be Gaussian, an analytical linear inversion technique is enabled. The amplitude posterior, being conditional on the frequency and prior information, has similiarities with the result of a conventional Fourier transformation analysis approach. However, the estimation of the frequency reveals a well localised domain which is in accordance with the measured data. The corresponding posterior for the frequency is orders of magnitude narrower than the width of the convolution function used by conventional approaches to estimate the spectral resolution. The prior information captured by one hyperparameter is tested as well, and its posterior peaks over a one order of magnitude and is very robust for significantly different noise levels on the data.

The findings for the monochromatic problem are investigated analytically. This reveals, for instance, the influence of the signal-to-noise ratio, the amount of data points and the spatial domains covered by the measurements on the results. The presented approach can be followed to examine more complex problems which involve more than one frequency and one amplitude and other diagnostic imperfections like a variable offset. Since this offset has a linear correspondence to the data domain, it would enter in the conditional amplitude posterior. Any non-linear parameter and additional hyperparameter would be estimated by the joint settings posterior specific for a certain problem. In doing so, a profound understanding of model implications could be established in analytical terms for applications relying on Fourier transformations. This is essential for designing a diagnostic for a given problem with somewhat known signal-to-noise ratio and hardware limitations, or for comparing results of different models in an objective way when a data set is given and, for example, the number of contributing frequencies is unknown.

Acknowledgements

This work has been carried out within the framework of the EURO fusion Consortium and has received funding from the Euratom research and training programme 2014-2018 and 2019-2020 under grant agreement No. 633053. The views and opinions expressed herein do not necessarily reflect those of the European Commission.

Appendix A: Integral Approximation

The spatial grid $x_i$ with $i\in \left[\mathrm{1,}N\right]$ may be given by the set $x_i=i\Delta x+x_0$ of uniformly distributed locations separated by the constant increment $\Delta x$ , so that $x_1=-L_{DS}$ and $x_N=-L_{DS}+L$. Then, the sum

$\begin{array}{l}{\displaystyle \underset{i=1}{\overset{N}{\sum }}}\cos \left(2\text{π}\frac{f_1}{c}{x}_i\right)\cos \left(2\text{π}\frac{f_2}{c}{x}_i\right)\\ \approx \frac{1}{2\Delta x}{\displaystyle \underset{-L_{DS}}{\overset{-L_{DS}+L}{\int }}}\left[\cos \left(2\text{π}\frac{f_{-}}{c}x\right)+\cos \left(2\text{π}\frac{f_{+}}{c}x\right)\right]\text{d}x\\ =\frac{1}{2\Delta x}{\left[\frac{c}{2\text{π}{f}_{-}}\sin \left(2\text{π}\frac{f_{-}}{c}x\right)\right]}_{-L_{DS}}^{-L_{DS}+L}+\frac{1}{2\Delta x}{\left[\frac{c}{2\text{π}{f}_{+}}\sin \left(2\text{π}\frac{f_{+}}{c}x\right)\right]}_{-L_{DS}}^{-L_{DS}+L}\\ =\frac{N}{2}\cos \left(\text{π}\frac{f_{-}}{c}{L}_{SS}\right)\text{sinc}\left(\pi \frac{f_{-}}{c}L\right)+\frac{N}{2}\cos \left(\text{π}\frac{f_{+}}{c}{L}_{SS}\right)\text{sinc}\left(\text{π}\frac{f_{+}}{c}L\right)\end{array}$ (A1)

is approximated by an integral over the spatial domain using the trigonometric identities (B1) and (B2). The approximation reveals two modulated sinc functions with a spatial and a spectral dependence. The spectral dependence is due to the sum and difference frequencies $f_{-}=f_1-f_2$ and $f_{+}=f_1+f_2$ , respectively. Spatially, the modulation depends on the length $L_{SS}=L-2L_{DS}$ of single-sided domain which is present on side only. The spatial dependence of the sinc function is given by the length L of the total domain covered.

Appendix B: Trigonometric Identities

The trigonometric identities are

$\cos \left(x\right)\cos \left(y\right)=\frac{1}{2}\left[\cos \left(x-y\right)+\cos \left(x+y\right)\right],$ (B1)

$\sin \left(x\right)-\sin \left(y\right)=2\cos \left(\frac{x+y}{2}\right)\sin \left(\frac{x-y}{2}\right).$ (B2)

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Sivia, D. and Skilling, J. (2006) Data Analysis: A Bayesian Tutorial. Oxford University Press, Oxford.
[2]	Jaynes, E.T. (1987) Maximum Entropy and Bayesian Spectral Analysis and Estimation Problems. Fundamental Theories of Physics, 21, 1-37.
[3]	Bretthorst, G.L. (1988). Bayesian Spectrum Analysis and Parameter Estimation. Springer-Verlag, Berlin, Heidelberg. https://doi.org/10.1007/978-1-4684-9399-3
[4]	Schmuck, S. and Svensson, J. (2017) Fourier Spectroscopy: A Bayesian Way. International Journal of Spectroscopy, 2017, Article ID: 9265084. https://doi.org/10.1155/2017/9265084
[5]	Schmuck, S., Fessey, J., Boom, J.E., Meneses, L., Abreu, P., Belonohy, E. and Lupelli, I. (2016) Electron Cyclotron Emission Spectra in X- and O-Mode Polarisation at JET: Martin-Puplett Interferometer, Absolute Calibration, Revised Uncertainties, Inboard/Outboard Temperature Profile, and Wall Properties. Review of Scientific Instruments, 87, Article ID: 093506. https://doi.org/10.1063/1.4962809