Field Dynamics Inference for Local and Causal Interactions

Inference of fields defined in space and time from observational data is a core discipline in many scientific areas. This work approaches the problem in a Bayesian framework. The proposed method is based on statistically homogeneous random fields defined in space and time and demonstrates how to reconstruct the field together with its prior correlation structure from data. The prior model of the correlation structure is described in a non‐parametric fashion and solely builds on fundamental physical assumptions such as space‐time homogeneity, locality, and causality. These assumptions are sufficient to successfully infer the field and its prior correlation structure from noisy and incomplete data of a single realization of the process as demonstrated via multiple numerical examples.


I. INTRODUCTION
Modeling the dynamic evolution of systems has always been at the very core of many scientific areas.While fundamental classical evolution is fully determined via ordinary and partial differential equations, many real-world models only adapt an effective probabilistic approach.Stochastic differential equations (SDE) provide an excellent way to construct such models for many dynamical systems.Generally SDEs consist of two constituents, a part describing the deterministic evolution of the system and another part describing stochastic evolution.The presence of a stochastic term in the equation can have a variety of reasons.It may be that the prescribed system is influenced by an external quantity which can only be described up to its statistical properties.Or the fundamental deterministic evolution is too complex to be fully modeled, for example due to non-linear self interactions or unknown inhomogeneous couplings, and should be approximated via a simpler (e.g.linear and homogeneous) deterministic evolution.All other influences can be approximated in terms of a stochastic constituent.
As a consequence, SDEs provide a great flexibility in order to decide up to which level of detail an effective description should capture deterministic evolution, and at which level a probabilistic description is sufficient.In the past, SDEs have been extensively used to simplify or approximate dynamic evolution, e.g. when modeling Brownian motion [1], in the Langevin dynamics [2] or the Black Scholes model [3].On the other hand SDEs can also be used in order to extract information about the dynamic evolution of a system from observations.The possibility to add stochastic constituents to SDEs is useful to decide a priori up to which level dynamic evolution is modeled.However, inference of SDEs is complicated due to a degeneracy between the two constituents, the deterministic and the stochastic term.As an illustrative example, consider a case of a simple random walk in time with an additional unknown linear drift.From observations we know which values the system obtained at several moments in time.Given only this incomplete information about a single realization of the stochastic process the reconstruction problem would be completely degenerate as it is impossible to decide whether observed temporal changes are due to the drift or due to the random walk.Only due to additional prior information such as stationarity of the drift it is possible to arrive at a non-degenerate interpretation.Degeneracies can appear on all levels of the inference and therefore prior modeling of the considered problem at hand is inevitable.In many applications this problem is solved by restricting a priori possible dynamics via parametrization with only a few free parameters which have to be inferred [4].However under certain conditions it is also possible to perform non-parametric inference.In [5] it was shown that for a linear autonomous system the spectrum of the corresponding Green's function, or propagator, can be reconstructed non-parametrically given a mild regularization on the spectrum as well as the statistics of the stochastic constituent.We will call these constituents excitations in the following.
The spectrum of the Green's function determines how the field reacts to excitations.The delay of the reaction of a system to an excitation is however still ambiguous.This potential time lag can be encoded in the phase of the Green's function, which therefore needs to be inferred as well.This appears to be problematic since in many realistic settings the statistical properties of the excitations can only be constrained up to the expected statistical amplitude, but it is not possible to a priori constrain the phase configuration of excitations.In this work we show which conditions still allow to reconstruct the full dynamical process in one-dimensional (only temporal) as well as multi-dimensional (spatio-temporal) settings.In particular we show that two fundamental concepts, namely locality and causality, appear to be sufficient to construct a prior for the dynamic Green's function that ultimately allows for a non-degenerate inference of it.
This paper structures as following: In section II we outline some basic properties of linear autonomous SDEs which are needed for the remaining discussion and section III introduces the prior assumptions that enter the inference problem.In section IV we outline a resulting inference algorithm given a generic measurement problem which is applied to several mock data examples in section V. Ultimately we conclude the paper in section VI with a brief summary of the results and an outlook to possible applications and extensions.

II. FIELD AND DYNAMICS
We are interested in reconstructing a field s, defined over space and time, and its dynamical law, governing its evolution in space-time.We assume that s solves a linear stochastic differential equation (SDE) of the form where L is a linear differential operator, and ξ denotes an unknown excitation field with assumed to be known statistics P (ξ).If we restrict the analysis to homogeneous and autonomous processes then L is a homogeneous operator, i.e.
Here x = (t, x) is a space-time vector, q(∂ xµ ) is a differential operator encoding function, and µ is an index of the coordinates.Fourier transforming eq. 1, by applying F with to it, and using eq. 2 yields the Fourier space representation of the differential equation which can be solved trivially where we also introduced the complex valued differential operator encoding field f k := q(ik).Note that if we want the model (eq. 1) to be real, both fields ξ k and f k have to be hermitian.Our goal is to infer s and its governing differential equation from observations of s alone.These observations may be noisy, incomplete and / or are retrieved via a complicated (possibly non-linear) instrument response.Consequently we rely on Bayesian inference in order to reconstruct the desired quantities.We note that the differential equation is fully characterized by f and ξ and we can also calculate the corresponding solution from them using eq. 5.This indicates that we can reformulate the problem such that we only aim to reconstruct f and ξ from observations.Before we can solve the corresponding measurement equation, we have to construct a meaningful prior for f .Up to this point, the statistics of ξ remain known.

III. DIFFERENTIAL OPERATOR ENCODING FIELD PRIORS
Physical processes that follow an autonomous differential equation usually are described in terms of a low order polynomial in the derivatives.While we do not want to restrict the analysis to such structures, we note that this property serves as an ingredient for the prior of f .In particular we aim to impose a certain degree of smoothness for f .What we can do is to construct a prior for f such that it punishes second (or any other) order derivatives of the field.One way to do so is to propose a Gaussian prior where the exponential factor reads For details see [6].

A. Symmetries
Due to the fact that we aim to reconstruct ξ and f from observations on s alone we notice that a corresponding likelihood will be invariant under simultaneous manipulations of these quantities.I.e.
where we split the manipulations r k into its amplitude a k and its phase φ k as we will treat them separately.
If we aim to use the effective model (eq. 1) only to construct a prior for the observed quantity s, there is no need to break this symmetry since by definition, it leaves the resulting s invariant.In this case the relation between ξ and s is only artificial, it is only used to construct a prior for s and we do not care if such an interaction is actually realizable in reality.
This picture completely changes once we aim to find a physical interpretation of f and ξ from observations of s.Although implemented as an effective model, there exist possible physical interpretations of the excitation field as an external force applied to the system, or a nonlinear self-interaction of the system.Therefore, it is of great interest whether or not it is possible to a posteriori analyze the resulting reconstruction regarding these questions.This is obviously not possible if there is a degeneracy between f and ξ as the reconstructed configuration is not unique.Therefore, in applications where the effective model we aim to infer should permit a specific physical interpretation, we need to break this degeneracy via priors which are not invariant under such transformations.Fortunately, the physical interpretation of ξ and f k provide the necessary prior constraints.

B. Amplitude degeneracy
In many applications the expected amplitude of the excitations can be assumed to be known (i.e. if the spectrum of ξ is known).In addition, we require f k to be a smooth function of k which already reduces the amplitude degeneracy.Throughout this work we restrict ourselves to such cases.

C. Phase degeneracy
The degeneracy resulting from phase manipulations is somewhat more problematic, since in many applications we cannot set the phases of the excitations.We notice that if the prior of ξ is assumed to be Gaussian distributed with a homogeneous covariance, such phase rotations wont affect it in any way.I.e.
In order to resolve this we need to provide additional information to the reconstruction problem.One trivial way would be to provide direct, additional measurements of the excitation, which we can use to fix the phase of ξ.However, in many realistic settings such measurements are difficult or even impossible to receive.
In cases where the statistics of ξ can be characterized via a distribution which is not invariant under local phase transformations, this prior knowledge also breaks the invariance of the posterior distribution, allowing for a phase sensitive reconstruction.For example if the excitations are assumed to be sparse, or if we deal with inhomogeneous distributed excitations.
Depending on the strength of such priors, or the quality of possible direct measurements, these approaches (partially) resolve the phase degeneracy.However, there exist also many applications in which we have to rely on less informative priors of the excitations, e.g. if we want to model multiple unknown external interactions or even unknown non-linear self-interactions via excitations.
Consequently, if we cannot break the phase degeneracy via manipulating ξ, we seek to find a way to constrain f such that the phase configuration is unique.To some extend this is already achieved via the smoothness prior.As we can see in Figure 1, the locality of a response of the system is related to the order in the derivatives of the operator.Since lower order derivatives are a priori favoured, the most local configuration that is consistent with the data is the one selected during reconstruction.This is perfectly motivated physically, since we know that for direct (classical) interactions, locality is a fundamental concept.Consequently apparent non-localities (realized via higher order derivatives) indicate additional hidden layers within the system, since the excitation does not influence the observed quantity s directly, but couples to an additional unobserved quantity which then couples to s.Such a configuration is of course possible and perfectly realizable, but if the system can also be described via a more direct interaction, the corresponding model should be favoured.This is achieved via our prior.
Although the smoothness prior determines localization, it does not set the temporal direction of propagation, as can be seen in Figure 2. Therefore the phase degeneracy corresponding to the direction of propagation remains and needs to be fixed.This is a well known degeneracy when constructing a response (Greens function) of a physical system.In general there exist both solutions, the advanced and the retarded Greens function as a response of a system, which could appear as a superposition in the reconstructed response.This issue can be addressed once we add causality as another physical principle besides locality.While it remains a topic to be discussed whether the excitations caused a response of the system or vice versa, a superposition of causal (retarded) and anti-causal (advanced) responses usually does not make any sense.
In order to fix this degeneracy we follow the standard way of obtaining retarded or advanced responses, by convolving with a step function Θ in time.This reads where Θ(±t) denotes a diagonal operator in space and time with Θ(±t) on its diagonal.
In principle we are free to choose if ξ should be the cause or the effect, but we notice that the way we constructed our model clearly indicates that ξ should cause the system to respond and not vice versa.We also want to emphasize at this point that this is a choice we make and the direction of causality is not part of our inference.I.e.we do not aim to perform causal inference, we rather define that ξ causes s.
Consequently we obtain a modified relation for our solution of s: where we chose G ≡ G ret. .Retardation trivially ensures that no anti-causal response of the system is allowed.
Ultimately we are left with one global phase degeneracy (φ k ≡ φ ).Noticing the requirement that the observed process (eq. 1) should be real, the global phase also has to be real.We get In other words, the only remaining degeneracy is an overall multiplication of the r.h.s. of eq. 1 with a minus sign.

D. Relation to minimum-phase
There are some structural similarities of our prior to a well known signal processing approach to phase configurations, namely the minimum-phase filter [7].
Excitation ξ and responses of various systems defined in terms of powers of f , with ft = ∂t + γ and γ small.We see that higher order derivatives result in apparent non-local responses.Since time derivatives are the generators of temporal translation, the case where the response is the exponential of −∆t f leads to translations by ∆t.
In signal processing, one often faces the problem of constructing a meaningful LTI (linear translation invariant) filter response, given a corresponding spectrum.A minimum phase filter can be constructed by relating the cepstrum (Fourier transformed logarithmic spectrum), with the phase via where denotes the cepstrum.A minimum phase filter is defined to be the causal and stable LTI filter with the least amount of group delay, corresponding to a given spectrum.This indicates that the resulting filter response is maximally localized in time.Here we notice the similarities to our prior: The smoothness prior favors more local and stable responses and retardation ensures causality of the response.However, in our approach these restrictions are formulated as a priori assumptions and therefore the posterior solution can deviate from these assumptions if there is enough evidence in the data to support this.This is not possible in a minimum-phase setup as the phase configuration is uniquely defined via the spectrum.In addition the inverse of the response of a minimum phase configuration is also causal and stable.While this is a desirable property in control systems, where the minimum phase filter is commonly used, we do not aim to make a statement on stability of the inverse of the response since there is not necessarily a meaningful physical interpretation of this response.This inverse should not be confused with a back-reaction of the system on excitations.These should also be causal and stable but have to be distinguished since they might have their own Green's function.
It is possible to incorporate the minimum phase condition as a forward model to construct corresponding propagators.This yields (16) with f being a priori distributed according to the smoothness prior given by eq. 6.

E. Light cone prior
So far, we restricted the discussion of causality to the temporal domain.As we generally aim to reconstruct the space-time evolution of a field, we notice that the physical concept of causality also implies a maximal finite propagation speed of interactions.This leads to propagation within a light cone, as depicted in figure 2, t.This information can be used to further restrict the reconstruction to propagators that are within the light cone, given a maximal propagation speed c.We can implement this constraint by rewriting eq.11 as where This ensures that the propagator is only non-zero within the light cone.In general we might expect different maximal propagation speeds in different directions, and consequently C becomes a symmetric tensor.In case propagation is isotropic in space with speed c we have C ij = c 2 δ ij .Note that even in cases where we do not know C, for example if the medium in which the interaction is realized is unknown, such a constraint can also be useful if we elevate C to be an additional unknown parameter that has to be inferred.Since C is the same for all scales, an inference algorithm can use large scale information, where the signal to noise ration (SNR) is usually higher, to determine C, which then effectively increases the SNR for smaller scales.A detailed description how C enters the reconstruction can be found in appendix A.

IV. INFERENCE
In general we aim to reconstruct ξ and the dynamical operator encoding field f from observations of s alone.Throughout this work we assume a generic measurement equation of the form where d denotes the measured data, R is an (in general non-linear) instrument response and n denotes the measurement noise.Furthermore we assume signal independent Gaussian distributed noise, with N being the known covariance structure of n.

A. Field inference only
Let us assume that the dynamics encoding field f is known to be f .Our goal is now to infer ξ given d and f .Note that s is also inferred, it can be calculated from f and ξ using eq. 5.The posterior reads (21) If we assume R to be a linear operator, and also assume that ξ is a priori Gaussian distributed with a known covariance structure Ξ, the posterior distribution is also a Gaussian where and This resembles the solution of the famous Wiener filter equations [8], applied to a field theoretical setting [9,10].

B. Inference of a field and its dynamics
Now we want to infer both ξ and f from the data.The posterior reads where s(ξ, f ) is defined according to eq. 12.Note that even for a linear measurement response R, the inference problem becomes intrinsically non-linear, due to the nonlinear dependency of s on f .We rely on a variational approach to solve the inference problem, where we approximate the posterior with a Gaussian distribution by minimizing the corresponding forward Kullback-Leibler divergence (KL).Here we follow the approach introduced via [11] and [12].Introducing the vector θ = (ξ, f, C) into which we include all parameters of the inference, we approximate the posterior via a Gaussian in θ where we also reconstruct the correlation structure of θ including all cross-correlations.This reads KL( P |P ) = H P P − H P P , where and Minimizing the KL w.r.t. the mean m and the covariance D leads to the desired solution.Following [11], we sample from the approximate posterior which enables us to calculate approximate posterior expectation values for any desired quantity.This is important since calculating s from θ involves non-linear operations and consequently s(θ) P = s(m).

V. APPLICATION
In order to demonstrate the applicability of our method we apply it to several mock data examples.
These include examples consistent with our prior assumptions as well as an inconsistent one.First, we demonstrate the performance of the algorithm in a one dimensional setting, where we only aim to infer the Green's function of the system, given the excitations and data.In the second example we perform a reconstruction of the excitations as well as the Green's function from data alone, in a 1+1 dimensional setting.Next, we add another level of complexity via non-Gaussian statistics in the excitations.Hereby we perform source detection, the task of inferring sparse excitations at various locations in space-time, in a case where also the Green's function is unknown.Ultimately we apply the reconstruction algorithm to data generated from a non-linearly evolving system, to demonstrate the performance in an inconsistent setup.
All examples were implemented using the NIFTy software package [13,14] in its version 5[15].This package allows to infer field quantities and to represent the remaining uncertainties via posterior samples.We included an implementation of the here introduced prior structure for the Green's function in this publicly available package.

A. Temporal evolution
In our first application, shown in fig.3, we aim to demonstrate the inference of the dynamics encoding field alone in a one dimensional setting where there is only a temporal evolution to be reconstructed.The excitation field is known during inference.We generate mock data according to eq. 19, with R being the identity.The excitations were drawn from an inverse gamma distribution to model known, sparse excitations of the system (i.e. in a laboratory setting where the unknown system is driven via sparse excitations).The mock signal φ as well as corresponding data d is shown in the top-left panel of fig. 3. The dynamical system used to generate signal and data is of the form: with ξ t being the inverse gamma distributed excitation field and (m, γ, m, γ) = (0.6, 0.23, 0.2, 0.03).The goal is now to reconstruct the differential operator encoding field f given d and ξ.The results are shown in fig. 3 as well.
We see that the reconstruction of the Green's function behaves as expected and is in agreement with the ground truth in the temporal domain, within uncertainties.Due to the fact that the reconstructed dynamics is uncertain, the recovered signal also has uncertainty although the excitations are known.The reconstructed Green's function indicates that it is indeed possible to reconstruct an apparent non-local response of the system (due to higher order derivatives in this setup) since the true propagator as well as its reconstruction show oscillations that grow in the beginning of the propagator before decaying exponentially.We also notice that there is relatively high uncertainty in the beginning of the propagator.This is caused by the low initial response to excitations of the true propagator.The initial part of the reconstructed propagator is purely dominated by noise and thus only constrained up to the noise level (the noise variance σ n is set to be σ n = 20 and the temporal domain is discretized via 512 equidistant pixels).
We also notice that the posterior solution for the spectrum levels out for small modes (large ω) below the signal to noise ratio, while the true spectrum shows a polynomial divergence (ultimately also the true spectrum levels out in the numerical example due to the finite size of the considered space).The uncertainty increases in this region, but not enough to capture the true solution.
Here we notice the limitations of the variational inference used by NIFTy5 [16], which provides a local approximation of the posterior with a Gaussian.Consequently the true uncertainty might be underestimated, as in this case.However this deviation occurs orders of magnitudes below the peak of the spectrum and therefore has only a barely visible effect on the reconstruction of the signal.One way to allow for a better extrapolation to higher frequencies would be to provide a more restrictive prior for the dynamics encoding field e.g. by defining it on a polynomial basis which is a more suitable basis for this particular setup as the true dynamics is also described in terms of a polynomial.However we aim to provide a general and less restrictive approach here capable of also reconstructing non-polynomial dynamics encoding functions and consequently being less able to extrapolate to regions where we have no information provided by data.

B. Space-time evolution
In our second example, shown in fig.4, we aim to reconstruct the dynamics as well as the excitations in a spatio-temporal (1+1 dimensional) setting from incomplete and noisy observations.We aim to infer the dynamical field f , the propagation speed c, and the excitations ξ from measurements d of φ alone.The dynamical system used for the generation of mock data is a product of a damped harmonic oscillator and an advection-diffusion generating term.This reads The signal φ, the data d and the reconstruction of φ are shown in fig. 4 for a case with (c, m, α, β, m, γ, β) = (0.4,0.16, −0.19, −0.05, 0.1, 0.5, 0.2).We generate the data according to eq. 19 with a linear measurement response which partially (≈ 25%) masks the observed region, and Gaussian distributed noise with σ n = 10.The space-time is discretized via a regular grid with 256×200 pixels, respectively.
The reconstruction algorithm is capable of reconstructing the signal in regions where we have observations thereof, while being relatively blind in unobserved regions.Consequently the posterior uncertainty is higher there.In addition, we notice that unlike the posterior mean, posterior samples consistently fill unobserved regions.Although in these regions the samples deviate strongly from the true signal, they obey the same statistical properties, which is important for posterior analysis.The results indicate that the statistical properties (and ultimately the dynamics encoding field) are correctly reconstructed.
As a validation, in figures 5 and 6, we compare the  reconstructed propagator as well as the corresponding spectrum with the underlying ground truth.The spectrum is comparable to the ground truth in regions with sufficient SNR while it levels out in regions where we have no information given via data.In addition, the reconstructed propagator also shows oscillations consistent with the true propagator.However we notice that modes that propagate "downwards" are reconstructed well while the weaker "upwards" propagating modes are not reconstructed due to the fact that they are below the noise level.In addition, we see that deviations in posterior samples of the propagator only occur within a "cone" and remain zero outside.This is due to the fact that we also reconstruct the maximal propagation speed of the process, which is reconstructed to be c ≈ 0.45 (±0.08) enclosing the correct value of c = 0.4.

C. Source detection
Due to the fact that the prior breaks the phase degeneracy of the problem, the algorithm is also able to perform source detection in the excitation field, even in cases where the dynamics is unknown.To demonstrate this scenario we again generate a 1+1-dimensional mock example where in this case we assume that we are only able to measure the temporal evolution of the system at several locations.In particular we measure the temporal evolution at 50 randomly selected locations of the space under consideration.This results in ≈ 80% of the discrete space-time being unobserved, as the resolution is the same as in the previous example.As before, we assume that the measurements are subject to additive Gaussian noise with σ n = 0.3 and also assume the system to be at rest at t = 0.The resulting data is shown in the top-left panel of fig. 7. The unknown excitations are inverse gamma distributed to model strong but sparse excitations.We infer those from measurements of the system at multiple locations together with the dynamic response.The system used to generate the data exhibits damped traveling waves described by with (c, m, α, β) = (1.2,0.04, −0.013, 0.005).We again aim to reconstruct the dynamics encoding field, the propagation speed and the excitations from the data.From an information theoretical point of view this setup is very similar to the previous one since we only changed the measurement response to describe measurements of the temporal evolution at several locations, as well as the prior for excitations to be an inverse gamma prior.Consequently also the inference can be performed in the same way as before.
The setup as well as the reconstruction of the field evolving in space-time are shown in fig. 7. We see that the reconstruction recovers many sources and the corresponding propagation.The algorithm uses information from the response of strong excitations to reconstruct the Green's function of the system.Due to the assumed homogeneity in space-time, this information helps to improve the overall reconstruction in other regions.The quality of the reconstruction of single excitations depends on the surrounding measurement scenario.
In figures 8 and 9 we depict the dynamic propagator as well as the spectrum and reconstructions.We can validate that the dynamical system was indeed reconstructed correctly, within uncertainties.We conclude that the task of source detection is possible even in cases where the underlying dynamics is unknown, as long as the assumptions of homogeneity, causality, and locality hold.

D. Non-linear evolution
In the last example we aim to study the performance of the proposed method in an inconsistent setup.In particular we want to study the results of a reconstruction with data generated from a system which evolves non-linearly.Our inference algorithm is designed for linear autonomous systems, however we want to see what the resulting reconstruction extracts from data obtained from such systems.

Setup
To do so, we generate a solution of the diffusive Burgers equation, from which we generate data.The differential equation reads where µ and α are constants.We expect that the linear part of the dynamics will be captured by the differential operator we aim to infer, while non-linear interactions consequently have to be explained via consistent excitations.If we rewrite φ as where V is the space-time volume of the total space under consideration, we can rewrite eq.33 as where we separated the linear from the non-linear evolution and identified the non-linear term with an excitation field ξ.
Here we already notice a major difference between the expected result in this case compared to the previous examples.We see that both, the excitations as well as the dynamical operator depend on φ, indicating that we cannot expect the reconstructed dynamics of one particular φ to generalize to a different solution, φ , of the equation.By definition, the Green's function of an autonomous system is the same for all solutions.Therefore, treating a Green's function recovered for one particular φ as a complete description leads to a wrong description of a system described via Burgers' equation.
This is an important observation, since in many applications, the goal of the inference of the dynamical system is to use the posterior solution for the Green's function as an abstract description of the system.Therefore, in the next section we will discuss how to identify whether or not it is possible to generalize the resulting dynamics to consistently produce different solutions of the observed dynamical system.
But first we have to validate if our expectations are met in an actual application.In figures 10 and 11 we show the mock solution as well as the reconstruction for (µ, α) = (0.06, 1).The spectrum of the linear dynamics matches the effective linear evolution from eq. 36 and the excitation field recovers the corresponding non-linear term.In this sense the reconstruction of φ is consistent with the data, but the reconstructed dynamical system does not fully describe the dynamical process.By definition, the method is only capable to infer the dynamic properties up to linear order.As we discussed above, the resulting linearized description is only valid for this particular solution and is no longer valid for a different solution of the Burgers equation.

Generalizing the reconstructed dynamics
The possibility of generalizing a reconstructed Green's function to different realizations of the observed system, strongly depends on the completeness of the description.To be precise, if the reconstruction captures all aspects of the evolution of the system, it can actually be regarded  as a fundamental and complete description of the system rather than an effective one.
However, in an actual application we usually do not know a priori whether a linear autonomous description of the system is sufficient.Therefore we aim to use posterior information in order to make a statement about completeness.The posterior solution is fully characterized in terms of two constituents: the reconstructed Green's function and the excitation field.As we have seen in the previous section, a non-linear interaction term results in a modified Green's function and an excitation field that depends on the solution itself.The modification of the Green's function is still a perfectly valid Green's function in terms of our prior assumptions.Therefore it cannot provide information about the completeness of the description.
Consequently we have to rely on the recovered excitations.A priori the statistical properties of the excitation field are assumed to be known.Therefore, if the posterior solution of the excitations is a typical realization of the prior statistics, we may conclude that the recovered dynamical system is indeed a complete description of the data.
On the other hand, if this is not the case, the description might be incomplete because it is unlikely to reproduce the observed system via sampling from the prior of the excitations.At this point, it is not possible to decide whether the dynamic evolution of the system is not linear or autonomous, or if the prior assumptions on the excitations are wrong.Therefore, we cannot make a statement about the origin of the inconsistency, only that the description is likely to be incomplete, given the data.
As an illustrative example we can take a look at the posterior power spectra of the excitations to test for typicality.We notice that the power spectrum of the excitation field for the inconsistent setup (left panel of fig.12), indicates that the posterior solution is correlated in space-time.In contrast, the power spectrum (right panel of fig.12) of the posterior excitations from the consistent example of section V B is almost flat.This indicates that the posterior excitations are uncorrelated, as it was assumed a priori.
Besides power spectra, there are a variety of additional properties that can be used in order to test the posterior result on consistency.For example the assumption of homogeneity can be tested, only to name one.The Bayesian approach to this question would be to calculate the evidence for this model and compare it to the evidence of a different model.However, calculating the evidence is usually very difficult even approximately for fields.
Nevertheless, questions regarding the origin of deviations from prior assumptions remain a field of speculation without further information about the system.A more detailed discussion of this topic is beyond the scope of this work.FIG.12. Double logarithmic power spectrum of the posterior solutions for excitations in the inconsistent example of section V D (left) and for the consistent example of section V B (right).

VI. CONCLUSION
In this work we considered the general problem of reconstructing the Green's function as well as the excitation field of a linear autonomous SDE, in a blind setting, where the only external information given is the statistical properties of the excitation field as well as noisy and incomplete measurements of one corresponding solution.We show that the principles of locality and causality are sufficient to break the degeneracy in the likelihood between excitations and the Green's function.Furthermore we construct a prior for the dynamic Green's function which makes use of this information and derive the posterior distribution for a generic measurement scenario with Gaussian additive noise.As the posterior distribution of the problem at hand is intractable, we rely on a variational approximation in order to solve the inference problem.
To demonstrate its applicability, the proposed method is then applied to several consistent as well as inconsistent mock data examples.The results indicate that the method behaves as expected.The prior structure breaks the degeneracy of the likelihood, allowing for a unique reconstruction.In addition, the prior assumptions appear to be reasonable for physical applications.
In general, the proposed method aims to provide a solution to the task of inferring linear autonomous SDEs for systems that evolve in space and time.Hereby we rely on minimal prior assumptions, locality and causality, necessary to break the degeneracy between excitations and the Green's function, which are physically motivated.This renders the method to be applicable in a wide range of problems.One particular strength is the non-parametric formulation of the Green's function.This becomes important in scenarios where physical models cannot provide a simple parametric description of evolution so far, to describe the Green's function.In addition, a probabilistic description of excitations is sufficient for inference.Consequently the method is still applicable in cases where external influences cannot be described completely.

FIG. 2 .
FIG.2.Left: Causal response of a system defined as a combination of two damped harmonic oscillators with different masses.Right: Anti-causal response of the system where one of the oscillations travels backwards in time.

FIG. 3 .
FIG. 3. Top:On the left we depict the signal (red line), the data (blue dots) as well as 100 over-plotted posterior samples (gray).The right panel shows the mock propagator (Greens function) in the temporal domain (red) as well as corresponding posterior samples.Bottom: Left: Excitation field used to generate the signal.Right: Logarithmic spectrum of the mock propagator (red) and posterior samples.

4 FIG. 4 .
FIG. 4. Top: The left panel shows the spatio-temporal masked and noisy data, drawn from the mock signal (middle panel) and the corresponding signal reconstruction (right panel).Bottom: Logarithmic one-sigma posterior uncertainty (left panel), logarithmic residual between the true signal and the posterior mean (middle panel), as well as an approximate posterior sample (right panel).

FIG. 5 . 7 FIG. 6 .
FIG. 5. Top: True Green's function (left) and corresponding posterior mean (right).Bottom: Residual of the true Green's function and the reconstruction (left) and a posterior sample for the Green's function (right).

FIG. 7 .
FIG. 7. Top: The left panel shows the sparse and noisy measurement data, drawn from the mock signal (middle panel) and the corresponding reconstruction (right panel).Bottom: Logarithmic one-sigma posterior uncertainty (left panel), logarithmic residual of the true signal and the corresponding posterior mean (middle panel), as well as an approximate posterior sample (right panel).

FIG. 8 . 10 FIG. 9 .
FIG. 8. Top: True Green's function of the process defined in eq.32 (left) and corresponding posterior mean (right).Bottom: Residual of the true Green's function and the reconstruction (left) and a posterior sample for the Green's function (right).

FIG. 10 .FIG. 11 .
FIG. 10.Top: Space-time evolution of the solution of Burgers' equation with an initial Gaussian profile at t = 0 (left) and the corresponding non-linear term as defined in eq.36 (right).Bottom: Reconstruction of the solution (left) and the corresponding excitation field (right).