From physical linear systems to discrete-time series. A guide for analysis of the sampled experimental data

TL;DR

This paper establishes a theoretical connection between discrete-time ARFIMA models and continuous-time linear stochastic differential equations, enabling better physical interpretation of sampled experimental data.

Contribution

It introduces formulas linking ARFIMA parameters to the coefficients of underlying physical stochastic systems, bridging discrete-time models with continuous-time physical models.

Findings

ARFIMA models can be viewed as sampled trajectories of continuous-time systems

Formulas relate ARFIMA parameters to physical system coefficients

Provides a framework for physically interpretable time series analysis

Abstract

Modelling physical data with linear discrete time series, namely Fractionally Integrated Autoregressive Moving Average (ARFIMA), is a technique which achieved attention in recent years. However, these models are used mainly as a statistical tool only, with weak emphasis on physical background of the model. The main reason for this lack of attention is that ARFIMA model describes discrete-time measurements, whereas physical models are formulated using continuous-time parameter. In order to remove this discrepancy we show that time series of this type can be regarded as sampled trajectories of the coordinates governed by system of linear stochastic differential equations with constant coefficients. The observed correspondence provides formulas linking ARFIMA parameters and the coefficients of the underlying physical stochastic system, thus providing a bridge between continuous-time linear dynamical systems and ARFIMA models.

Equations

$${\Delta }^d\left(X_n-\sum _{k=1}^p{\varphi }_k{X}_{n-k}\right)={\xi }_n+\sum _{j=1}^q{\theta }_j{\xi }_{n-j},$$

$${\rho }_X(k)=\frac{⟨X_n{X}_{n+k}⟩}{\sqrt{⟨X_n^2⟩⟨X_{n+k}^2⟩}}={\mathrm{e}}^{-{\varphi }_{1}k}.$$

$${\Delta }^{d_f}{X}_n=\sum _{k=0}^{\infty }\frac{d_f(d_f-1)\cdots (d_f-k+1)}{k!}(-1{)}^k{X}_{n-k}.$$

$$f_X(\omega )={\int }_{-\infty }^{\infty }d\tau ⟨X(t)X(t+\tau )⟩{\mathrm{e}}^{-\mathrm{i}\omega \tau }.$$

$$f_X^{\Delta t}(\omega )=\frac{1}{\Delta t}\sum _{k=-\infty }^{\infty }⟨X_n{X}_{n+k}⟩{\mathrm{e}}^{-\mathrm{i}\omega \Delta tk},$$

$$\underset{\Delta t\to 0}{\lim }{f}_X^{\Delta t}=f_X;$$

$$\frac{d}{dt}\mathbf{S}(t)=A\mathbf{S}(t)+\mathbf{F}(t),$$

$$m\frac{d^2}{d{t}^2}X(t)=-\kappa X(t)-\beta \frac{d}{dt}X(t)+F(t),$$

$$L\frac{d^2}{d{t}^2}Q(t)+R\frac{d}{dt}Q(t)+\frac{1}{C}Q(t)=\mathcal{E}(t),$$

$$\mathbf{S}(t)={\int }_{-\infty }^{t}ds{\mathrm{e}}^{A(t-s)}\mathbf{F}(s).$$

$$\mathbf{S}(n\Delta t)$$

$$+{\int }_{(n-1)\Delta t}^{n\Delta t}ds{\mathrm{e}}^{A(n\Delta t-s)}\mathbf{F}(s),$$

$${\mathbf{S}}_n=E{\mathbf{S}}_{n-1}+{\mathbf{\Xi }}_n.$$

$${\mathbf{S}}_n=\sum _{k=0}^{\infty }{E}^k{\mathbf{\Xi }}_{n-k}.$$

$$E={\mathrm{e}}^{A\Delta t},{\mathbf{\Xi }}_n={\int }_{(n-1)\Delta t}^{n\Delta t}ds{\mathrm{e}}^{A(n\Delta t-s)}\mathbf{F}(s).$$

$$p(E)=E^N-\sum _{k=1}^N{\varphi }_k{E}^{N-k}=0.$$

$${\varphi }_k=(-1{)}^{k+1}\sum _{D_k}\prod _{i\in D_k}{\lambda }_i,$$

$${\varphi }_k=(-1{)}^{k+1}\sum _{D_k}{\mathrm{e}}^{\Delta t{\sum }_{i\in D_k}{\nu }_i}.$$

$${\mathbf{S}}_{n-N}$$

$${\mathbf{S}}_{n-N+1}$$

$${\mathbf{S}}_{n-N+2}$$

$${\mathbf{S}}_{n-N+3}$$

$$⋮$$

$${\mathbf{S}}_n$$

$${\mathbf{S}}_n-\sum _{k=1}^N{\varphi }_k{\mathbf{S}}_{n-k}={\mathbf{\Xi }}_n+\sum _{k=1}^{N-1}{R}_k{\mathbf{\Xi }}_{n-k}.$$

$${\mathbf{\eta }}_n={\mathbf{\Xi }}_n+\sum _{k=1}^{N-1}{R}_k{\mathbf{\Xi }}_{n-k}$$

$$R_k=E^k-\sum _{j=1}^k{\varphi }_j{E}^{k-j}$$

$$S_n^1-\sum _{k=1}^N{\varphi }_k{S}_{n-k}^1={\eta }_n^1.$$

$$\frac{\Delta S_n^1}{\Delta t}-\sum _{k=1}^N{\varphi }_k^{\prime }{S}_{n-k}^1=\frac{{\eta }_n^1}{\Delta t},$$

$${\eta }_n^1={\xi }_n^1+\sum _{j=1}^{N-1}{\theta }_j^1{\xi }_{n-j}^1,$$

$$X(t)=\frac{1}{\beta }{\int }_{-\infty }^{t}ds{\mathrm{e}}^{-\frac{\kappa }{\beta }(t-s)}\xi (s)$$

$$f_{\text{OU}}(\omega)=\frac{1}{\beta^{2}}\frac{\sigma^{2}}{\big{(}\frac{\kappa}{\beta}\big{)}^{2}+\omega^{2}},$$

$$f_{\text{AR(1)}}^{\Delta t}(\omega )=\left({\varphi }_1^{-2}-1\right)\frac{1}{2\kappa \beta }\frac{{\sigma }^2\Delta t}{1+{\varphi }_1^2-2{\varphi }_1\cos (\omega \Delta t)}.$$

Full text

ARMA, ARFIMA, FARIMA, Langevin equation

I. INTRODUCTION

Discrete-time series methods based on ARMA (Autoregressive Moving Average) model, and more recently, its generalisation ARFIMA (Fractionally Integrated ARMA, also called FARIMA) granger ; BDbook ; beran provide powerful and flexible statistical tools which were successful in analysing data in econometrics (which resulted in 2003 Nobel Prize in Economic Sciences for C. W. J. Granger and R. Engel), finance and engineering crato ; gil-alana ; fou . It is a model which fully describes the behaviour of time series using small number of parameters, which can be estimated from the data using well-established techniques and widely available statistical packages BDbook ; farima . Moreover, these techniques allow for control of the estimation’s quality, checking the correctness of the model or even forecasting future values of the time series. In recent years, new physical solarFlare ; Nose , biological burneckiBiol and medical MEDarma applications of ARFIMA model were found, allowing for empirical description of complex systems with long (power-like), short (exponential) and finite-range dependencies STarma ; farima , see Fig. 1. ARMA processes were also studied as models of physical data governed by discrete-time Langevin equations GLEarma ; MZarma .

The main physical interpretation of this model was based on the fact that ARFIMA approximates processes such as fractional Brownian motion, Lévy stable motions farima ; conv and the corresponding noises, whereas its special case ARMA can model properties of various stationary processes with finite or exponentially decaying memory BDbook ; beran . However, most of these continuous-time processes themselves reflect rather behaviour of the process then its internal physical dynamics. On the other hand, the mathematical theory proposed by A. W. Philips philips , in recent years developed further by P. Brockwell, R. Davies and Y. Yang BD , establish a connection between ARMA model and a class of continuous-time stochastic dynamical systems. Here we show that the physical importance of these results is significant, and after suitable refinement, this connection establishes reliable physical basis for ARMA and ARFIMA models.

II. ARFIMA MODEL

The studied model ARFIMA($p,d,q$) states that the considered time series $X_{n}$ fulfils the recursive relation farima

[TABLE]

where $\phi_{k},\theta_{j},d$ are deterministic coefficients and $\xi_{n}$ is white noise which generates the stochastic dynamics; it is Gaussian or non-Gaussian, e.g. $\alpha$-stable farima ; kokoszka which determines the distribution of $X_{n}$. The above equation is comprised of three basic building blocks: AR part (left parentheses), FI part (operator $\Delta^{d}$) and MA part (right side). Each of these blocks models different type of memory and has distinct interpretation. If no memory is present we deal with ARFIMA(0,0,0), which is a white noise: $X_{n}=\xi_{n}$.

The left parentheses of Eq. 1 is AR($p$) (Auto Regressive) part in which coefficients $\phi_{k}$ determine how the present value of the time series $X_{n}$ depends linearly on the past values $X_{n-k}$; in other words it models internal dynamics of the system. Because this dynamics is linear, it describes the exponential components of the memory. The most basic process from this class is called AR(1) or ARMA(1,0), which is the simplest exponential memory process with correlation function BDbook

[TABLE]

More general AR($p$) models have memory functions which are mixtures of exponential decays BDbook . These processes have great importance for statistics because AR($p$) are maximal entropy processes for fixed first p+1 values of the autocovariance function entr .

The right side of Eq. 1 is MA($q$) process (Moving Average) which determines the external dynamics through modification of the white noise. It models finite-range components of the memory which depend on coefficients $\theta_{j}$. Actual value of MA($q$) process depends only on $q$ last values of the the generating noise $\xi_{n}$, and because of that it does not contain any information about the history of the process older then $q\Delta t$. For example, process called MA(1) or ARMA(0,1) is the simplest type of coloured noise with ultra short memory and correlation function $\rho_{\text{MA(1)}}(1)=\theta_{1}$, $\rho_{\text{MA(1)}}(k)=0$ for $k>1$ BDbook ; BJbook .

The last part of Eq. 1, the operator $\Delta^{d}$ denoting the FI (Fractional Integration) part reflects both non-stationarity and fractional memory. The symbol $\Delta$ denotes the discrete difference operator: $\Delta X_{n}=X_{n+1}-X_{n}$. When $d$ is a natural number, the non-stationary process ARFIMA($p,d,q$), in this case also called ARIMA($p,d,q$), is understood as a process witch after $d$ differentiations is stationary ARMA($p,q$). Basic example is ARFIMA(0,1,0) which is summed white noise, that is sampled trajectory of Brownian motion.

In a situation when $d$ is real, it can be decomposed into natural-number part $d_{n}$ and fractional remainder $d_{f}$, $-1/2\leq d_{f}\leq 1/2$, such that $\Delta^{d}=\Delta^{d_{n}}\Delta^{d_{f}}$. This remainder accounts for the power-law type of memory common for, e.g., anomalous diffusion anom . This operator is understood as a series farima

[TABLE]

Applying $\Delta^{-d_{f}}$ to both sides of Eq. 1, it can be confirmed that ARFIMA($p,d,q$) can be interpreted as modification of ARFIMA($p,d_{n},q$) generated by noise $\Delta^{-d_{f}}\xi_{n}$. This noise, which is denoted as FI($d_{f}$) or ARFIMA($0,d_{f},0$) is a stationary process with power-law memory function witch has a tail $\sim n^{2d_{f}-1}$; when $\xi_{n}$ are Gaussian this time series is very similar to discrete-time fractional Brownian noise conv ; barkai ; kou .

III. CONTINOUS- VERSUS DISCRETE-TIME PROCESSES IN EXPERIMENTS

A continuous-time process $X(t)$ in a natural way contain a lot more information than its discrete-time counterpart $X_{n}=X(n\Delta t)$, see Fig. 2. During sampling we lose information, e.g. about the geometrical properties of trajectories. Memory functions become discrete and do not contain information about the dependence within intervals smaller than sampling rate $\Delta t$. However, only in discrete time we can define some more refined memory functions, like partial autocorrelation function, which is correlation of $X_{n}$ and $X_{n+k}$ with the influence of all in-between $X_{n+j}$ removed BJbook ; BDbook , which in many cases has a simple form which leads to a greater usability.

Some state functions, like power spectral densities (psd) differ considerably for discrete and continuous time. These are memory functions of the considered process in the Fourier space. Continuous-time psd (cpsd) of the process $X$: $f_{X}$, is most easily defined as a Fourier transform of the covariance function

[TABLE]

But, in contrast to psd $f_{X}$, discrete-time power spectral density (dpsd) is a Fourier series of the sampled covariance function

[TABLE]

and is a periodic function witch repeats after $2\pi/\Delta t$. The relation between these two functions is given by Poisson summation formula spec ; spec2 . It allows to calculate dpsd numerically or analytically given cpsd, but for the processes considered in this paper it is not a very practical tool. However, using Poisson summation formula, one can prove that dpsd converges to cpsd as $\Delta t\to 0$ BJbook ; spec

[TABLE]

this fact can be interpreted as convergence of the time series $X_{n}$ to the process $X(t)$ as $\Delta t\to 0$ (this limit is often called infill asymptotics). But, in realistic conditions, we often are far away from this limit and only discrete-time model properly reflect the behaviour of the observed system. Note also that all distortions of the data caused by the measurement equipment are changes of the sampled series $X_{n}$ as this is the object that is actually processed by the hardware. Thus, accounting for unwanted effects like blur or measurement noise must be performed in the discrete time setting ja . In order to perform this procedure, the discrete-time model of the undistorted observations is needed.

IV. LINEAR DYNAMICAL SYSTEMS

To link these two notions using ARFIMA model we will consider a linear stochastic system with $N$-dimensional state vector $\boldsymbol{S}(t)=[S^{1}(t),S^{2}(t),\ldots,S^{N}(t)]^{\mathrm{T}}$, which evolves in continuous time according to stochastic differential equation of first ordersobczyk

[TABLE]

where $A$ is $N\times N$ matrix with constant coefficients and $\boldsymbol{F}(t)$ is some stationary noise, acting as a random force. This general model describes class of systems with time-independent environment and additive stochastic disturbance. Note, that if we study state vector $\boldsymbol{S}(t)$ in which some coordinates $S^{i}(t)$ are described by differential equations of order bigger than one, we can complement the state vector $\boldsymbol{S}(t)$ by auxiliary coordinates $\{\frac{d}{dt}S^{i}(t),\frac{d^{2}}{dt^{2}}S^{i}(t),\ldots\}$ and also reduce the problem to the form 7. One well-known property described by Eq. 7 is a position of particle trapped in the harmonic potential within liquid lang

[TABLE]

where $m$ is the mass of particle, $\kappa$ is stiffness of the harmonic trap, $\beta$ is friction coefficient of the liquid and $F(t)$ is white noise modelling the exchange of momenta with surrounding particles. Phase plot of this equation is shown in Fig. 3. Other examples include evolution of charge $Q(t)$ in a linear RLC circuit disturbed by noise electromotive force $\mathcal{E}(t)$

[TABLE]

as well as other types of linear disturbed circuits RLC , harmonic heat bath models zwanzig , Brownian magnetic particle in constant magnetic field Bmag and many more.

V. TIME DISCRETISATION PROCEDURE (SAMPLING)

In any of these cases the real system evolves in continuous time, but the experimental observations must be discrete and usually have the form of time series $\boldsymbol{S}_{n}$ sampled with constant sampling time $\Delta t$: $\boldsymbol{S}_{n}=\boldsymbol{S}(n\Delta t)$. The necessary condition for $\boldsymbol{S}(t)$ and $\boldsymbol{S}_{n}$ to be stationary is for matrix $A$ to be negative-definite, in other words it needs to have eigenvalues with negative real part. In such case there exists stationary solution of 7 given by convolution of the force $\boldsymbol{F}(t)$ with matrix exponential $\mathrm{e}^{At}$ sobczyk

[TABLE]

From elementary properties of the matrix exponential and integration it follows that

[TABLE]

that is the sampled process $\boldsymbol{S}_{n}$ fulfils the equation

[TABLE]

The obtained recursive formula is a vector counterpart of the process AR(1) and is called VAR(1) (Vector AR) BDbook . The VAR(1) process $\boldsymbol{S}_{n}$ can be explicitly expressed in terms of the generating noise $\boldsymbol{\Xi}_{n}$ as

[TABLE]

The above equation is a discrete counterpart of the convolution formula 10.

From Eq. 13 it fallows that neither matrix $A$ directly, nor all values of $\boldsymbol{F}(t)$ affect the state $\boldsymbol{S}_{n}$, they do this only through their discrete-time counterparts

[TABLE]

The left part of the above equation defines the discretized evolution operator $E$, which is responsible for the deterministic part of the transition from state $\boldsymbol{S}_{n-1}$ to $\boldsymbol{S}_{n}$. The stochastic deviation from the deterministic path is described by the right part of 14, which defines the packing operator $\boldsymbol{F}(t)\mapsto\boldsymbol{\Xi}_{n}$. The whole influence of the process $\boldsymbol{F}(t)$ on $\boldsymbol{S}_{n}$ is fully determined by the packed force $\boldsymbol{\Xi}_{n}$. Each of values $\boldsymbol{\Xi}_{n}$ gathers values of $\boldsymbol{F}(t)$ from interval $((n-1)\Delta t,n\Delta t)$. If $\boldsymbol{F}(t)$ is stationary, also the discrete noise $\boldsymbol{\Xi}_{n}$ is stationary. The packed force process inherits most of the properties of the underlying continuous-time $\boldsymbol{F}(t)$: the same type of distribution (Gaussian $\boldsymbol{F}(t)$ gives Gaussian $\boldsymbol{\Xi}_{n}$, $\alpha$-stable $\boldsymbol{F}(t)$ gives $\alpha$-stable $\boldsymbol{\Xi}_{n}$) and the same type of memory (white noise, finite range, exponential or power law).

Because of the above correspondence, the statistical methods available for VAR(1) model pfaff can be used in physical applications. Using estimators for matrix $E$ and taking matrix logarithm we obtain estimates for the underlying matrix $A=\ln(E)/\Delta t$. The properties of force $\boldsymbol{F}(t)$ can be studied through analysis of $\boldsymbol{\Xi}_{n}$, which can be estimated as $\boldsymbol{\Xi}_{n}=\boldsymbol{S}_{n}-E\boldsymbol{S}_{n-1}$.

VI. SINGLE COMPONENT ANALYSIS

The above approach has a huge disadvantage: it requires that we can observe the whole sampled state vector $\boldsymbol{S}_{n}$, which is unrealistic for more complex systems. To avoid this difficulty we want to analyse the behaviour of one chosen component of $\boldsymbol{S}_{n}$, without loss of generality we will assume it is $S_{n}^{1}$. The evolution of this component is coupled with the evolution of the rest of state vector by action of non-diagonal components of the evolution matrix $E$.

In order to decouple component $S_{n}^{1}$, we use Cayley-Hamilton theorem, which states that if $N\times N$ matrix $E$ has characteristic equation $p(\lambda)=\det(E-\lambda I)=\lambda^{N}-\sum_{k=1}^{N}\phi_{k}\lambda^{N-k}=0$, then the matrix $E$ itself fulfils this equation, that is

[TABLE]

This polynomial equation of order $N$ has coefficients $\phi_{k}$ determined by the eigenvalues $\lambda_{i}$ of the matrix $E$

[TABLE]

where $D_{k}$ denotes family of all $k$-element subsets of the set $\{1,2,\ldots,N\}$. The eigenvalues $\nu_{i}$ of matrix $A$ from 7 and the eigenvalues of matrix $E$ are related by formula $\lambda_{i}=\mathrm{e}^{\nu_{i}\Delta t}$. Therefore $E$ fulfils 15 with coefficients

[TABLE]

As a preparation for using Cayley-Hamilton theorem, we recursively use 12 and express variables $\boldsymbol{S}_{n-k}$ as functions of $\boldsymbol{S}_{n-N}$

[TABLE]

thus coupling them to only this one variable. After multiplying these equations by $\phi_{k}$ in order to obtain terms $\phi_{k}E^{N-k}\boldsymbol{S}_{n-N}$, and subtracting them from the last one, we remove all action of the matrix $E$ on the time series $\boldsymbol{S}_{n}$ using equality $(E^{N}-\sum_{k=1}^{N}\phi_{k}E^{N-k})\boldsymbol{S}_{n-N}=0$ obtained from the Cayley-Hamilton theorem. The cost of this decoupling is performing complicated transformations of the discretized force $\boldsymbol{\Xi}_{n}$ along the way. The equation that we obtain after this procedure has a form

[TABLE]

The left side of the above equation is AR($N$) part described by coefficients $\phi_{k}$, which depend only on the deterministic matrix $A$. The left side acts as an effective noise $\boldsymbol{\eta}_{n}$

[TABLE]

which generates the stochastic dynamics of the vector $\boldsymbol{S}_{n}$. The behaviour of this noise is determined by the matrices $R_{k}$

[TABLE]

composed of mixtures of time-shifting operators $E^{k-j}$. They affect the evolution by mixing different $\boldsymbol{\Xi}_{n-k}$; as a result the first component $S_{n}^{1}$ fulfils the recurrence relation

[TABLE]

There is a deep connection between the above formula and the classical Mori-Zwanzig theory zwanzig ; mori . The above equation can be written in slightly different manner

[TABLE]

where $\phi_{1}^{\prime}=(\phi_{1}+1)/\Delta t,\phi_{k}^{\prime}=\phi_{k}/\Delta t,k>1$. Now it becomes clear that this is the discrete-time analogue of the generalised Langevin equation kou ; zwanzig ; mori . It is no coincidence, as derivations in both cases use the same ideas, decoupling most of the coordinates of freedom, at the same time introducing the memory kernel and effective noise. Despite many similarities, this analogy has its limitations, e.g. there is no clear discrete-time equivalent of the fluctuation-dissipation theorem.

VII. ANALYSIS OF THE STOCHASTIC FORCE

The effective noise $\eta_{n}^{1}$ is a mixture of $N-1$ last values and all components of $\boldsymbol{\Xi}_{n}$. If the underlying force $\boldsymbol{F}(t)$ was white Gaussian noise, then the resulting $\eta_{n}^{1}$ are a mixture of Gaussian white noises which forms a time series with $N-1$ finite-range dependence. When analysing only one component, we can ignore its internal structure and represent it as MA($N-1$) process

[TABLE]

generated by a white noise $\xi^{1}_{n}$, which is in fact orthogonalised series $\eta^{1}_{n}$ MAortog . Such orthogonalisation can be always performed and coefficients $\theta_{j}^{1}$ may be obtained solving system of equations resulting from comparing the covariance function of left side of Eq. 24 and effective noise $\eta_{n}^{1}$ BDbook .

Thus, we have come to conclusion that coefficient $S_{n}^{1}$ is ARMA($N,N-1$) process with coefficients $\phi_{k}$ and $\theta_{j}^{1}$. Similar statements hold for more complex models of the force $\boldsymbol{F}(t)$.

If the force $\boldsymbol{F}(t)$ has finite range of memory smaller than $K\Delta t$, then $\eta_{n}^{1}$ by the same reasoning as above can be regarded as MA($N-1+K$) process. In such case $S_{n}^{1}$ is ARMA($N,N-1+K$).

If the force $\boldsymbol{F}(t)$ has power-law memory tails $\sim t^{2d_{f}-1}$, then $\eta_{n}^{1}$ is a composition of finite-range mixing introduced by operators $R_{k}$ and the power-law behaviour. The best approximation of such time series is FIMA($d_{f},q$), i.e. process similar to 24, but where $\xi_{n}^{1}$ are FI($d_{f}$) time series farima . So, component $S_{n}^{1}$ is ARFIMA($N,d_{f},q$).

If the force $\boldsymbol{F}(t)$ has exponential tails of memory, that is, if it can be represented in form similar to 10, then we may treat it as a time-dependent state of the same class as $\boldsymbol{S}(t)$, which confirms that $\eta_{n}^{1}$ is ARMA($L,L-1+N$) for some $L$; $+N$ results from mixing of $\boldsymbol{\Xi}_{n-k}$. The AR($L$) part, understood as an operator, can be freely moved from acting on $\eta_{n}^{1}$ to acting on $S_{n}^{1}$, leading to ARMA($N+L,N+L-1$) model. The operator AR($N+L$) can be easily calculated as composition AR($L$)AR($N$).

In our considerations we assumed that $\boldsymbol{S}(t)$ was stationary. But, we may observe the possibly non-stationary integral of the stationary coordinate $X(t)=\int_{0}^{t}d\tau\ S^{1}(\tau)$, $X(t)$ being position, charge, etc. In this case the process of differences $\Delta X_{n}=\int_{(n-1)\Delta t}^{n\Delta t}d\tau\ S^{1}(\tau)$ is ARFIMA($p,d_{f},q+1$) which follows from the fact that $S^{1}(\tau+n\Delta t)$ is some ARFIMA($p,d_{f},q$) for any $\tau$, with $p,d_{f},q$ determined by the proper model from the description above. Increase by one order in MA part accounts for additional short-time memory introduced by the integral $\int_{(n-1)\Delta t}^{n\Delta t}ds$. And so, $X_{n}$ is ARFIMA($p,d_{f}+1,q+1$).

Analogical statement holds for all other components $S^{i}_{n}$. The AR coefficients $\phi_{k}$ are identical for all of them, MA coefficients $\theta^{i}_{j}$ and noises $\xi^{i}_{n}$ from Eq. 24 vary. The equations governing evolution of different $S^{i}_{n}$ were decoupled, however these components are dependent, because $\xi_{n}^{i}$ for different $i$ are mixtures of $\boldsymbol{\Xi}_{n}$ components, so they are dependent set of variables.

We stress that all modelling is performed at the level of the discretized stochastic force, the obtained ARFIMA model of the observed coordinate is by no means phenomenological as is often the case of discrete-time models, but derived from the theory of the continuous-time dynamical system. For all cases except the power-law memory the correspondence is exact, in the latter case FIMA approximation must be made for the discretized stochastic force, AR part is still exact. As FI part can reflect any type of power-law long-time memory asymptotics and MA part can account for any finite-range deviations, such model is most often sufficient burneckiBiol ; farima .

In our reasoning we used the fact that Gaussian process is fully determined by its covariance structure at the moment when we orthogonalised the effective noise series $\eta_{n}^{1}$. For non-Gaussian processes it is no longer true, as they can have richer than linear memory structure kokoszka . Therefore, for non-Gaussian forces $\boldsymbol{F}(t)$, the obtained ARFIMA model reflect only linear aspect of the memory. In this case it is approximate, but has the same autocovariance and power spectral density as the original process.

VIII. PARTICLE IN A HARMONIC POTENTIAL

Let us come back to Eq. 8, which is the second order differential equation describing the particle trapped in harmonic potential. An approximation, in which the inertial term $m\,d^{2}X/dt^{2}$ is considered negligible, simplifies analysis. In such conditions the state of particle evolves according to force-balance equation $\beta\,dX/dt=-\kappa X+\xi$. Its stationary solution,

[TABLE]

is called Ornstein-Uhlenbeck process OU ; protter , and has a well-known Lorenzian continuous-time power spectral density lang ; ja

[TABLE]

where $\sigma^{2}$ is the variance of the noise $F(t)$; when fluctuation-dissipation theorem holds $\sigma^{2}=2k_{B}T\beta$ zwanzig . Sampled trajectory of 25 is AR(1) process with a coefficient $\phi_{1}=\mathrm{e}^{-\Delta t\kappa/\beta}$. The discrete-time power spectral density can be calculated from the general formula for all ARFIMA processes spec ; beran , which for AR(1) yields ja

[TABLE]

As we see cpsd and dpsd functions differ when $\Delta t$ is not considerably smaller then $\beta/\kappa$ (see Fig. 4), which is often the case for mesoscopic objects observed in normal conditions ja . As we provide exact formula of the observed spectral density, there is no reason to use approximate Eq. 26.

In a case when the mass is not negligible, the effects of sampling become more complex. The full state vector $\boldsymbol{S}$ is then composed of position and momentum, $\boldsymbol{S}(t)=[X(t),P(t)]^{\mathrm{T}}$. The stochastic force affects the change of momenta $\boldsymbol{F}(t)=[0,F(t)]^{\mathrm{T}}$. Note that because of the identical form of Eq. 8 and Eq. 9, all further results would follow also for the RLC circuit 9 after a simple change of letters. In this case the state vector would consists of charge and electric current.

Calculating the eigenvalues

[TABLE]

of the evolution matrix A=\bigl{[}\begin{smallmatrix}0&1/m\\ -\kappa&-\beta/m\end{smallmatrix}\bigr{]}, we obtain the AR coefficients of the sampled position process $X_{n}$ which, if $F(t)$ is a white noise, is ARMA(2,1) with AR(2) coefficients

[TABLE]

The MA coefficient $\theta^{1}_{1}$ is also determined by the calculated eigenvalues and is given by complicated formula, but can be easily calculated numerically. Estimating the AR(2) coefficients from the data, the ratios $\kappa/m$ and $\beta/m$ of the underlying process can be assessed and the parameter $m$ can be subsequently estimated from variance of the sampled process $X_{n}$. If $F(t)$ is not white noise, the MA part may differ and if $F(t)$ would have the power-law dependence it would be reflected in the FI part of ARFIMA model.

IX. SUMMARY

In our work we tried to construct a bridge between continuous-time linear dynamical systems and discrete-time ARMA or, more generally, ARFIMA models. The studied correspondence for many cases might serve as physical interpretation of the ARFIMA model and justification for its usage. Additionally, we have shown what order physical ARFIMA model should have for given dynamical system and we have given explicit formulas for its AR coefficients, which allows for estimation the dynamical system’s parameters using standard statistical tools. Its MA and FI coefficients can also be calculated, but they depend on the assumed model of the stochastic force. The coefficients of ARFIMA model determine its characteristics, such as power spectral density, linking the basic dynamical system model with functions that can be estimated from the sampled data measured during experiment.

This research was supported by NCN Maestro Grant No. 2012/06/A/ST1/00258.