Inferring the Dynamics of Underdamped Stochastic Systems

Many complex systems, ranging from migrating cells to animal groups, exhibit stochastic dynamics described by the underdamped Langevin equation. Inferring such an equation of motion from experimental data can provide profound insight into the physical laws governing the system. Here, we derive a principled framework to infer the dynamics of underdamped stochastic systems from realistic experimental trajectories, sampled at discrete times and subject to measurement errors. This framework yields an operational method, Underdamped Langevin Inference, which performs well on experimental trajectories of single migrating cells and in complex high-dimensional systems, including flocks with Viscek-like alignment interactions. Our method is robust to experimental measurement errors, and includes a self-consistent estimate of the inference error.

Figure 1

Inference from discrete time series subject to measurement error. (a) Trajectory $x(t)$ of a stochastic damped harmonic oscillator, $F(x,v)=-\gamma v-kx$. (b) The same trajectory represented in $xv$-phase space. Color coding indicates time. (c) Force field in $xv$-space inferred from the trajectory in (a) using ULI with basis functions $b=\{1,x,v\}$ (blue arrows), compared to the exact force field (black arrows). Inset: inferred components of the force along the trajectory versus the exact values. (d) Convergence of the mean squared error of the inferred force field, obtained using ULI (circles) and with the previous standard approach [13, 14, 25, 27] (squares). Dashed lines indicate the predicted error $\delta {\widehat{F}}^2/{\widehat{F}}^2\sim N_b/2{\widehat{I}}_b$. (e) Inferred friction coefficient $\gamma$ divided by the exact one as a function of the sampling time interval $\mathrm{\Delta }t$, comparing the previous standard approach to ULI. (f) Trajectory $y(t)=x(t)+\eta (t)$ (blue) corresponding to the same realization $x(t)$ in (a), with additional time-uncorrelated measurement error $\eta (t)$ (orange) with small amplitude $|\eta |=0.02$. (g),(h) Force field inferred from $y(t)$ using estimators without and with measurement error corrections, respectively. (i) Inference convergence for data subject to measurement error using estimators without (circles) and with (diamonds) measurement error corrections. (j) Dependence of the inference error on the noise amplitude $|\eta |$ [same symbols as in (i)].

Figure 2

Inferring nonlinear dynamics and multiplicative noise. (a) $xv$ trajectory of the stochastic Van der Pol oscillator, $F(x,v)=\kappa (1-x^2)v-x$ with measurement error. (b) Partial information of the 28 basis functions of a sixth-order polynomial basis in natural information units (1 nat $=1/\log 2\mathrm{bits}$), inferred from the trajectory in (a). Inset: Corresponding force field reconstruction. (c) Convergence of the inference error for the $d$-dimensional Van der Pol oscillator $F_{\mu }(\mathbf{x},\mathbf{v})={\kappa }_{\mu }(1-x_{\mu }^2)v_{\mu }-x_{\mu }$ (no summation, $1\le \mu \le d$) with $d=1,\dots ,6$, using a third-order polynomial basis. (d) Microscopy image of a migrating human breast cancer cell (MDA-MB-231) confined in a two-state micropattern (scale bar: $20\mu \mathrm{m}$). Experimental trajectory of the cell nucleus position, recorded at a time interval $\mathrm{\Delta }t=10\min$ (blue), and simulated trajectory using the inferred model (red). (e) Partial information for the experimental trajectory in (d), projected onto a third-order polynomial basis. (f) Deterministic flow field inferred from the experimental trajectory in (d). (g) Trajectory of a Van der Pol oscillator with multiplicative noise ${\sigma }^2(x,v)={\sigma }_0+{\sigma }_x{x}^2+{\sigma }_v{v}^2$ (colormap). (h),(i) Inferred versus exact components of the force and noise term, respectively, for the trajectory in (g). (j) Inference convergence of the multiplicative noise amplitude, using Eq. (3) without measurement error (circles), with measurement error (squares), and using the error-corrected estimator (diamonds). The error saturation at large $\tau$ is due to the finite time step. Dashed line: predicted error $\delta \widehat{{\sigma }^2}/\widehat{{\sigma }^2}\sim \sqrt{N_b\mathrm{\Delta }t/\tau }$ [8].

Figure 3

Interacting flocks. (a) Trajectory (green) of $N=27$ Viscek-like particles [Eq. (5)] in the flocking regime (1000 frames). We perform ULI on this trajectory using a translation-invariant basis of pair interaction and alignment terms, both fitted with $n=8$ exponential kernels. (b) Exact (blue) and inferred (orange) cohesion $rf(r)$. Exact form includes short-range repulsion and long-range attraction, $f(r)={\epsilon }_0[1-(r/r_0{)}^3]/[(r/r_0{)}^6+1]$. Dotted inference dependence indicates distances not sampled by the initial data. (c) Exact and inferred alignment kernel $g(r)$. Exact form: $g(r)={\epsilon }_1\exp (-r/r_1)$. (d) Inferred versus exact components of the force field. (e) Convergence of the inferred force as a function of trajectory length. Dashed line is the predicted error $\delta {\widehat{F}}^2/{\widehat{F}}^2\sim N_b/2{\widehat{I}}_b$. (f) Simulated trajectory (red) employing the inferred force and noise, showing qualitatively similar flocking behavior.