Estimating the Coefficients of a System of Ordinary Differential Equations Based on Inaccurate Observations

Abstract

In this paper, we solve the problem of estimating the parameters of a system of ordinary differential equations from observations on a short interval of argument values. By analogy with linear regression analysis, a sufficiently large number of observations are selected on this segment and the values of the functions on the right side of the system and the values of the derivatives are estimated. According to the obtained estimates, unknown parameters are determined, using the differential equations system. The consistency of the estimates obtained in this way is proved with an increase in the number of observations over a short period of argument values. Here, an algorithm for estimating parameters acts as a system. The error of the obtained estimate is an indicator of its quality. A sequence of inaccurate measurements is a random process. The method of linear regression analysis applied to an almost linear regression function is used as an optimization procedure.

1. Introduction

The problem of estimating the parameters of a system of nonlinear ordinary differential equations, based on inaccurate deterministic observations, using known optimization algorithms, is solved in the papers [1,2,3]. An alternative approach for estimating the parameters of a deterministic recurrent sequence, observed with random additive and multiplicative errors, based on the relationships between the trajectory averages and their approximation from inaccurate observations, is proposed in [4,5].

The advantage of the first approach is the possibility of using known optimization algorithms, and the disadvantage of it is the lack of analytical estimates of the convergence rate to the estimated parameters. The advantage of the second approach is the availability of theoretical estimates of the convergence rate to the estimated parameters, and the disadvantage of it is the need to establish limit cycles or limit distributions for recurrent sequences.

Despite all the differences in these approaches, the common fact is that by increasing in the length of the observation segment, the accuracy of estimates increases and, under certain conditions, may tend to zero. At the same time, the problem of estimating parameters over a small observation interval is interesting, which is closely related to discrete optimization methods of experiment planning (see, for example, [6,7]).

In this paper, this problem is solved for a system of non-linear ordinary differential equations. At the same time, the estimation of the parameters of this system, based on inaccurate observations, is solved under the assumption that a large number of observations may be carried out over a relatively short segment. To estimate the parameters, the method of linear regression analysis is used in relation to a regression function that slightly deviates from the original function in a small neighbourhood of some time moment [8,9,10,11,12,13].

This method is based on minimizing the standard deviation of a sequence of observations from a linear regression function. In this case, such a relationship is selected between the number of observations and the interval between neighbouring observations so that the resulting error in determining the parameters tends to zero when the number of observations tends to infinity.

The final stage of the parameter estimation algorithm is the substitution of estimates of the values of functions and the values of their derivatives into the original system of equations at the selected point. Further, by analogy with the method of moments, unknown parameters of the system of equations are estimated and the consistency of the estimates obtained is proved. This paper also uses the implicit function theorem, which allows us to establish that the obtained parameter estimates are consistent depending on the number of observations. Based on the results obtained, computational experiments were carried out.

Thus, elements of system analysis have been introduced into the solution of the task. Here, an algorithm for estimating parameters acts as a system. The error of the obtained estimate is an indicator of its quality. A sequence of inaccurate measurements is a random process. Furthermore, the process and the method of linear regression analysis applied to an almost linear regression function is used as an optimization procedure. It is evaluated using the theorem on the existence and uniqueness of the solution of a system of ordinary differential equations and with the help of the implicit function theorem. Additionally, known error estimates in the linear regression analysis method are used.

2. Estimating the Coefficients of a System of Ordinary Differential Equations by Inaccurate Observations

2.1. Preliminaries

Consider a system of ordinary differential equations with fixed values of parameters ${\beta }_i={\beta }_i^0,i=1,\dots ,m,$

$$\frac{d{x}_i}{dt}=F_i(x_1,\dots ,x_m,{\beta }_1^0,\dots ,{\beta }_m^0),i=1,\dots ,m,$$

(1)

where $x_1=x_1\left(t\right),\dots ,x_m=x_m\left(t\right)$ are unknown functions. In well-known monographs on the theory of ordinary differential equations (see, for example, [14,15]), the theorem of the existence and uniqueness of the solution of this system in a small neighbourhood of a certain point is formulated and proved in Theorem 1.

Theorem 1.

Assume that functions $F_i=F_i(x_1,\dots ,x_m,{\beta }_1^0,\dots ,{\beta }_m^0)$ are continuous in a rectangular parallelepiped $Q=\{(x_1,...,x_m)\in R^m:x_i^0-a_i\le x_i\le x_i^0+a_i,i=1,..,m\}$ together with their partial derivatives ${\displaystyle \frac{\partial F_i}{\partial x_i},i=1,...,m.}$ Then there is a segment $t_0-r\le t\le t_0+r,$ on which the system of Equation (1) has a unique solution satisfying the initial conditions $x_i\left(t_0\right)=x_i^0,i=1,\dots ,m.$

Remark 1.

From the Weierstrass theorem for continuous functions on a compact, it follows that the functions $x_i\left(t\right),i=1,...,m,$ on the segment $[t_0-r,t_0+r]$ (continuity follows from differentiability) and function ${\displaystyle \left|F_i·\frac{\partial F_i}{\partial x_i}\right|,i=1,...,m,}$ on a set Q (due to the continuity of the multipliers) reach their highest final values $C_i.$

Denote

$$M_0=(x_1^0,\dots ,x_m^0,{\beta }_1^0,\dots ,{\beta }_m^0),F_i\left(M_0\right)=F_i^0,i=1,\dots ,m,$$

$$M_0^{\prime }=(x_1^0,\dots ,x_m^0,F_1^0,\dots ,F_m^0,{\beta }_1^0,\dots ,{\beta }_m^0),$$

$$G_i(x_1,\dots ,x_m,f_1,\dots ,f_m,{\beta }_1,\dots ,{\beta }_m)=F_i(x_1,\dots ,x_m,{\beta }_1,\dots ,{\beta }_m)-f_i,$$

where $F_i$ are described in Theorem 1, and consider the system of equations

$$G_i(x_1,\dots ,x_m,f_1,\dots ,f_m,{\beta }_1,\dots ,{\beta }_m)=0,i=1,\dots ,m.$$

(2)

In monographs on mathematical analysis (see, for example, [16,17]), conditions are formulated, under which the system (2) may be resolved with respect to variables ${\beta }_1,\dots ,{\beta }_m$ (see for example Theorem 2).

Theorem 2.

If the functions $G_i,i=1,\dots ,m$ are continuously differentiable in the neighbourhood of the point $M_0^{\prime }$ and the Jacobian

$$\frac{\partial (G_1,\dots ,G_m)}{\partial ({\beta }_1,\dots ,{\beta }_m)}{|}_{M_0^{\prime }}\ne 0,$$

(3)

then there are neighbourhoods $U,V,W$ of points $(x_1^0,\dots ,x_m^0),$$(F_1^0,\dots ,F_m^0),$$({\beta }_1^0,\dots ,{\beta }_m^0)$, respectively, such that the system of Equation (2) is uniquely solvable in the neighbourhood of $U\times V\times W$ of the point $M_0^{\prime }$ relative to the variables ${\beta }_1,\dots ,{\beta }_m.$ Moreover, if ${\beta }_i=g_i(x_1,\dots ,x_m,f_1,\dots ,$$f_m),i=1,\dots ,m,$ is the specified solution, then all functions $g_i$ are continuously differentiable in the neighbourhood $U\times V$ and ${\beta }_i^0=g_i(x_1^0,\dots ,x_m^0,F_1^0,\dots ,F_m^0).$

Remark 2.

When the conditions of Theorem 2 are met, the functions $g_i,i=1,...,m,$ are continuous at the point $(x_1^0,\dots ,x_m^0,F_1^0,\dots ,F_m^0).$

2.2. Ordinary Differential Equation

Consider the differential equation for a fixed value of the parameter $\beta ={\beta }_0$

$$\frac{dx}{dt}=F(x,{\beta }_0)$$

(4)

with the initial condition $x\left(0\right)=x_0,$ assuming that the function $F(x,\beta )$ is continuously differentiable in the neighbourhood of a point $M_0=(x_0,{\beta }_0)$ and ${\displaystyle \frac{\partial F}{\partial \beta }{|}_{M_0}\ne 0}$. Let the inaccurate observations $y\left(t\right)=x\left(t\right)+\epsilon \left(t\right)$ are known for the state of $x\left(t\right)$ at the moments $t=kh,k=0,\pm 1,\dots ,\pm n,hn\le r.$ Denote

$${\epsilon }_k=\epsilon \left(kh\right),x_k=x\left(kh\right),y_k=y\left(kh\right)=x_k+{\epsilon }_k,F_0=F(x_0,{\beta }_0)$$

and suppose that ${\epsilon }_k,k=0,\pm 1,\dots ,\pm n,$ is a set of independent and identically distributed random variables with zero mean and variance ${\sigma }^2.$ The problem of estimating the parameter ${\beta }_0$ of the differential Equation (4) from these observations is posed.

The solution of this problem is carried out in two stages. First, they are constructed using a modification of the least squares estimation method ${\widehat{x}}_0,{\widehat{F}}_0$ and their convergence to the estimated parameters $x_0,F_0$ is investigated. Then, by analogy with the method of moments, an estimate of ${\widehat{\beta }}_0$ is constructed and its convergence to the estimated parameter ${\beta }_0$ is investigated.

Evaluation of values $x_0,F_0$. Let us introduce the notations, outlining the method for defining ${\widehat{x}}_0,{\widehat{F}}_0$

$${\widehat{x}}_0=\frac{{\sum }_{k=-n}^n{y}_k}{2n+1},{\widehat{F}}_0=\frac{{\sum }_{k=-n}^n{y}_{k}kh}{{\sum }_{k=-n}^n{\left(kh\right)}^2}.$$

(5)

Theorem 3.

If ${\sigma }^2<\infty$ and $h=n^{-\alpha },$ then, for $\alpha >1$, the estimate of ${\widehat{x}}_0$ is an asymptotically unbiased and consistent estimate of the parameter $x_0.$ The estimate ${\widehat{F}}_0$ is an asymptotically unbiased estimate of the parameter $F_0.$ At $1<\alpha <3/2$; the estimate ${\widehat{F}}_0$ is a consistent estimate of $F_0.$

Proof of Theorem 3. 

Denote ${\tilde{y}}_k=x_0+F_{0}kh+{\epsilon }_k$ and put

$${\tilde{x}}_0=\frac{{\sum }_{k=-n}^n{\tilde{y}}_k}{2n+1},{\tilde{F}}_0=\frac{{\sum }_{k=-n}^n{\tilde{y}}_{k}kh}{{\sum }_{k=-n}^n{\left(kh\right)}^2}.$$

Estimates of ${\tilde{x}}_0,{\tilde{F}}_0$ are obtained by the least squares method for coefficients $x_0,F_0$ of linear regression [9] and satisfy the following relations

$$E{\tilde{x}}_0=x_0,E{\tilde{F}}_0=F_0,Var{\tilde{x}}_0=\frac{{\sigma }^2}{2n+1},Var{\tilde{F}}_0=\frac{{\sigma }^2}{{\sum }_{k=-n}^n{\left(kh\right)}^2}.$$

(6)

Here, $Ex$ is mathematical expectation of arbitrary random variable x and $Varx=E{(x-Ex)}^2$ is its variance. In turn, the following equalities are almost certainly fulfilled

$${\widehat{x}}_0-{\tilde{x}}_0=\frac{{\sum }_{k=-n}^n({\widehat{y}}_k-{\tilde{y}}_k)}{2n+1},{\widehat{F}}_0-{\tilde{F}}_0=\frac{{\sum }_{k=-n}^n({\widehat{y}}_k-{\tilde{y}}_k)kh}{{\sum }_{k=-n}^n{\left(kh\right)}^2}.$$

(7)

Moreover, the differences ${\widehat{y}}_k-{\tilde{y}}_k=x_k-x_0-F_{0}kh,k=0,\pm 1,\dots ,\pm n$ are deterministic quantities.

The Remark 1 implies the existence of a number C satisfying the inequality

$$\underset{\left|t\right|\le nh}{sup}\left|x^{″}\left(t\right)\right|=\underset{\left|t\right|\le nh}{sup}\left|\frac{\partial F(x\left(t\right),{\beta }_0)}{\partial x}F(x\left(t\right),{\beta }_0)\right| =2C<\infty .$$

Then, from the Taylor formula with a residual term in the Lagrange form,

$$x\left(kh\right)=x\left(0\right)+F_{0}kh+\frac{{\left(kh\right)}^2}{2}{x}^{″}\left(kh{\tau }_k\right),\left|{\tau }_k\right|\le 1,k=0,\pm 1,\dots ,\pm n,$$

inequalities follow

$$|x_k-x_0-F_0{kh|\le C\left(kh\right)}^2,k=0,\pm 1,\dots ,\pm n.$$

(8)

From the Formulas (7) and (8) for $n\to \infty ,$ the relations follow

$$|{\widehat{x}}_0-{\tilde{x}}_0| \le \frac{{\sum }_{k=-n}^n|x_k-x_0-F_{0}kh|}{2n+1}\le \frac{2C{h}^2{\sum }_{k=1}^n{k}^2}{2n+1}\sim \frac{C{h}^2{n}^2}{3},$$

(9)

$$|{\widehat{F}}_0-{\tilde{F}}_0| \le \frac{{\sum }_{k=-n}^n\left|(x_k-x_0-F_{0}kh)kh\right|}{{\sum }_{k=-n}^n{\left(kh\right)}^2}\le \frac{C{h}^3{\sum }_{k=1}^n{k}^3}{{\sum }_{k=1}^n{h}^2{k}^2}\sim \frac{Chn}{4}.$$

(10)

The Formulas (6), (9) and (10) lead to the relations

$$|E{\widehat{x}}_0-x_0| = |E{\widehat{x}}_0-E{\tilde{x}}_0| ⪯\frac{C{h}^2{n}^2}{2},Var{\widehat{x}}_0=Var{\tilde{x}}_0,$$

(11)

$$|E{\widehat{F}}_0-F_0| = |E{\widehat{F}}_0-E{\tilde{F}}_0| ⪯\frac{3Chn}{4},Var{\widehat{F}}_0=Var{\tilde{F}}_0.$$

(12)

Here $a_n⪯b_n$ means that ${\displaystyle \underset{n\to \infty }{lim\; sup}{a}_n/b_n\le 1.}$ Then from the condition $h=n^{-\alpha },\alpha >1,$ and the Relations (11) and (12) we have

$$|E{\widehat{x}}_0-x_0| \to 0,|E{\widehat{F}}_0-F_0| \to 0,n\to \infty ,$$

(13)

that ${\widehat{x}}_0,{\widehat{F}}_0$ are asymptotic unbiased estimates of $x_0,F_0.$

From the Bieneme–Chebyshev inequality, the Relations (9) and (11) and the conditions $h=n^{-\alpha },\alpha >1,$ we get for any $\delta >0$

$$\begin{array}{c}P\left(\right|{\widehat{x}}_0-x_0| >\delta )\le P\left(\right|{\widehat{x}}_0-{\tilde{x}}_0| + |{\tilde{x}}_0-x_0\left|\right)>\delta )=P\left(\right|{\tilde{x}}_0-x_0| \ge \delta - |{\widehat{x}}_0-{\tilde{x}}_0\left|\right)\le \\ \le \frac{{\sigma }^2}{(2n + 1)(\delta - |{\widehat{x}}_0 - {\tilde{x}}_0{\left|\right)}^2}\to 0,n\to \infty .\end{array}$$

(14)

Thus, for $h=n^{-\alpha },\alpha >1,$ estimate ${\widehat{x}}_0$ is a consistent estimate of $x_0.$

At the same time, from the Relations (10), (12) and (13) for $h=n^{-\alpha },1<\alpha <3/2,$ we get for any $\delta >0$

$$\begin{array}{c}P\left(\right|{\widehat{F}}_0-F_0| >\delta )\le P\left(\right|{\widehat{F}}_0-{\tilde{F}}_0| + |{\tilde{F}}_0-F_0\left|\right)>\delta \left)\right)=\\ P\left(\right|{\tilde{F}}_0-F_0| >\delta - |{\widehat{F}}_0-{\tilde{F}}_0\left|\right)\le \frac{3{\sigma }^2}{h^2{n}^3(\delta - |{\widehat{F}}_0 - {\tilde{F}}_0{\left|\right)}^2}\to 0,n\to \infty .\end{array}$$

(15)

Therefore, if the condition $h=n^{-\alpha },1<\alpha <3/2,$ is true, the estimate ${\widehat{F}}_0$ is a consistent estimate of $F_0.$ □

Remark 3.

It is worth noting that Theorem 3 is true for any distribution of random variables ${\epsilon }_k$ with finite variance ${\sigma }^2.$ Indeed it is necessary to prove limit relation ${\displaystyle H_n=\frac{{\sum }_{k=-n}^n{\epsilon }_{k}k}{h{\sum }_{k=-n}^n{k}^2}\to 0,n\to \infty .}$ However, the most reasonable way to solve this question is to consider such distributions of random variables ${\epsilon }_k$ as normal for ${\sigma }^2<\infty /$ or stable for ${\sigma }^2=\infty ,$ because $H_n$ has normal/stable distribution also.

Evaluation of parameter ${\beta }_0.$ Consider the equation

$$F({\widehat{x}}_0,\beta )={\widehat{F}}_0.$$

(16)

Theorem 4.

In conditions of Theorem 3, Equation (16) has a unique solution ${\widehat{\beta }}_0,$ which is a consistent estimate of the parameter ${\beta }_0.$

Proof of Theorem 4. 

Since the function $F(x,\beta )$ is continuously differentiable in the neighbourhood of the point $M_0=(x_0,{\beta }_0)$ and ${\displaystyle \frac{\partial F}{\partial \beta }{|}_{M_0}\ne 0,}$, then the conditions of the theorem for the function $G(x,f,\beta )=F(x,\beta )-f$ are fulfilled. So, in some neighbourhood of the point $M_0^{\prime }=(x_0,F_0,{\beta }_0)$, the equation is solvable with respect to $\beta =g(x,f),$ while ${\beta }_0=g(x_0,F_0).$ Then, from the Remark 2, we get that for any $\epsilon >0$ there exists $\delta \left(\epsilon \right)>0$ such that in the neighbourhood $\{(x,f):|x-x_0|\le \delta \left(\epsilon \right),|f-F_0|\le \delta \left(\epsilon \right)\}$ of the point $(x_0,F_0)$ the inequality $|\beta - {\beta }_0| \le \epsilon$ is executed.

It follows that with the specified choice of $\delta \left(\epsilon \right)$, the relation is fulfilled

$$|\widehat{x_0}-x_0| \le \delta \left(\epsilon \right),|\widehat{F_0}-F_0| \le \delta \left(\epsilon \right)\Rightarrow |{\widehat{\beta }}_0-{\beta }_0| \le \epsilon .$$

(17)

In turn, from the Relations (14) and (15) it follows that for any $\epsilon$ and $\delta \left(\epsilon \right)$ there is such a $n_0(\epsilon ,\delta \left(\epsilon \right)),$ that for any $n>n_0(\epsilon ,\delta \left(\epsilon \right))$ inequalities are fair

$$P\left(\right|{\widehat{x}}_0-x_0| \le \delta \left(\epsilon \right))\ge 1-\frac{\epsilon }{2},P\left(\right|{\widehat{F}}_0-F_0| \le \delta \left(\epsilon \right))\ge 1-\frac{\epsilon }{2}.$$

(18)

Therefore, from the Relations (17) and (18) we have $P(|{\widehat{\beta }}_0-{\beta }_0|\le \epsilon )\ge 1-\epsilon .$ Thus, for any $\epsilon >0$, there exists $n_0\left(\epsilon \right)$ such that for $n>n_0\left(\epsilon \right)$, the inequality holds $P(|{\widehat{\beta }}_0-{\beta }_0|>\epsilon )<\epsilon ,$ which means consistency (convergence in probability at $n\to \infty$) of the constructed estimate. □

2.3. System of Differential Equations

Consider a system (1) with initial conditions $x_i\left(0\right)=x_i^0,i=1,\dots ,m.$ We assume that the functions $F_i(x_1,\dots ,x_m,{\beta }_1,\dots ,{\beta }_m),i=1,\dots ,m,$ are continuously differentiable in the neighbourhood of the point $M_0$ and the Jacobian ${\displaystyle \frac{\partial (F_1,\dots ,F_m)}{\partial ({\beta }_1,\dots ,{\beta }_m)}{|}_{M_0}\ne 0.}$ Inaccurate observations are known $y_i\left(t\right)=x_i\left(t\right)+{\epsilon }_i\left(t\right)$ for the state $x_i\left(t\right),i=1,\dots ,m,$ at moments $t=kh,k=0,\pm 1,\dots ,\pm n,hn\le r.$ Let ${\epsilon }_i\left(kh\right),k=0,\pm 1,\dots ,\pm nh,i=1,\dots ,m,$ is a set of independent and identically distributed random variables with zero mean and variance ${\sigma }^2.$ The task is to estimate the vector of parameters $({\beta }_1^0,\dots ,{\beta }_m^0)$ of a system of differential Equation (1) based on these observations.

Denote

$${\widehat{x}}_i^0=\frac{{\sum }_{k=-n}^n{y}_i\left(kh\right)}{2n+1},{\widehat{F}}_i^0=\frac{{\sum }_{k=-n}^n{y}_i\left(kh\right)}{{\sum }_{k=-n}^n{\left(kh\right)}^2},i=1,\dots ,m.$$

(19)

Theorem 5.

If ${\sigma }^2<\infty$ and $h=n^{-\alpha },$ then, for $\alpha >1$, the estimate ${\widehat{x}}_i^0$ is an asymptotically unbiased and consistent estimate of the parameter $x_i^0.$ The estimate ${\widehat{F}}_i^0$ is an asymptotically unbiased estimate of the parameter $F_i^0.$ For $1<\alpha <3/2$, the estimate ${\widehat{F}}_i^0$ is a consistent estimate of the value $F_i^0,i=1,\dots ,m.$

Consider a system of equations

$$F_i({\widehat{x}}_1,\dots ,{\widehat{x}}_m,{\beta }_1,\dots ,{\beta }_m)=F_i^0,i=1,\dots ,m.$$

(20)

Theorem 6.

In conditions of Theorem 5 the system of Equation (20) has a unique solution $({\widehat{\beta }}_1^0,\dots {\widehat{\beta }}_m^0),$ which is a consistent estimate of the vector of parameters $({\beta }_1^0,\dots ,{\beta }_m^0).$

The proofs of the Theorems 5 and 6 almost verbatim repeat the proofs of the Theorems 3 and 4.

Remark 4.

Theorems 3–6 are devoted to ordinary differential equations of the first order and their systems. However, it is possible to spread them to ordinary differential equations and their systems of arbitrary order. For this purpose it is possible to use for examples results of [12,13].

2.4. Computational Experiment

Example 1.

The computational experiment was conducted first for the Cauchy problem

$$\frac{dx}{dt}=F(x,{\beta }_0)={\beta }_{0}x,x\left(0\right)=1,{\beta }_0=0.5.$$

The solution of this equation has the form$x=e^{b_{0}t}.$We assumed that by observing the process described by this equation,$\pm kh,k=0,1,...,n,h=n^{-5/4},n=10,000,$inaccurate observations were obtained at time points$y_{\pm k}=e^{\pm b_{0}hk}+{\epsilon }_{\pm k},k=0,1,...,n.$

Here, independent random variables ${\epsilon }_{\pm k},k=0,1,...,n,$ are distributed uniformly on the segment $[-1/2,1/2]$ left and on the segment $[-1/4,1/4]$ right. According to the Formula (5), the parameters $x_0,F_0=F(x_0,{\beta }_0)$ in our notation ${\widehat{x}}_0,{\widehat{F}}_0,$ were evaluated first, then the formula for evaluating the parameter ${\beta }_0$ was found from the equation ${\widehat{F}}_0={\widehat{x}}_0{\widehat{\beta }}_0.$ 

Table 1. Interval distribution of estimate ${\widehat{\beta }}_0$ when ${\epsilon }_{\pm k}$ has variance ${\sigma }^2=1/12$ left and variance ${\sigma }^2=1/48$ right.

Distribution Intervals	Relative Frequencies	Distribution Intervals	Relative Frequencies
0.387413–0.432612	0.027	0.437796–0.460889	0.011
0.432612–0.477811	0.238	0.460889–0.483983	0.171
0.477811–0.52301	0.477	0.483983–0.507076	0.494
0.52301–0.568209	0.229	0.507076–0.53017	0.29
0.568209–0.613408	0.029	0.53017–0.553263	0.034

shows the results of a computational experiment conducted 1000 times, namely, the interval distribution (5 intervals) of relative frequencies ${\widehat{\beta }}_0.$

Consequently, a decrease in variance ${\sigma }^2$ improves the quality of the obtained estimates sufficiently clearly.

Now, consider the case in which independent random variables ${\epsilon }_{\pm k},k=0,1,...,n,$ are distributed normally with mean 0 and variance ${\sigma }^2.$ 

Table 2. Interval distribution of estimate ${\widehat{\beta }}_0$ when ${\epsilon }_{\pm k}$ has variance ${\sigma }^2=1/12$ left and variance ${\sigma }^2=1/48$ right.

Distribution Intervals	Relative Frequencies	Distribution Intervals	Relative Frequencies
0.393462–0.442011	0.062	0.435015–0.458419	0.014
0.442011–0.49056	0.328	0.458419–0.481824	0.144
0.49056–0.53911	0.458	0.481824–0.505228	0.47
0.53911–0.587659	0.143	0.505228–0.528632	0.322
0.587659–0.636208	0.009	0.528632–0.552036	0.05

shows the results of a computational experiment conducted 1000 times, namely, the interval distribution (five intervals) of relative frequencies ${\widehat{\beta }}_0.$

Consequently, the quality of obtained results for disturbances distributed normally behaves like in a case of uniform distribution.

Example 2. A computational experiment was also carried out for the system of Lorentz equations

$$\left\{\begin{array}{c}{\displaystyle \frac{dx}{dt}=F_1(x,y,z,{\sigma }_0,r_0,b_0)={\sigma }_0(y-x),}\hfill \\ {\displaystyle \frac{dy}{dt}=F_2(x,y,z,{\sigma }_0,r_0,b_0)=x(r_0-z)-y,}\hfill \\ {\displaystyle \frac{dz}{dt}=F_3(x,y,z,{\sigma }_0,r_0,b_0)=xy-b_{0}z,}\hfill \end{array}\right.$$

(21)

with the given initial conditions $x\left(0\right)=1,y\left(0\right)=2,z\left(0\right)=1,$ in the case of ${\sigma }_0=1,r_0=2,b_0=3$. The solution of this system is not written out explicitly, but it is solved by the finite difference method. We write out the corresponding equations for the grid $\{\pm kh,k=0,1,...,n\}$ in increments of $h=n^{-5/4},n=\mathrm{10,000}:$

$$\left\{\begin{array}{c}{x}_{\pm (k+1)}=x_{\pm k}\pm {\sigma }_{0}h(y_{\pm k}-x_{\pm k}),\hfill \\ y_{\pm (k+1)}=y_{\pm k}\pm h(x_{\pm k}(r_0-z_{\pm k})-y_{\pm k}),\hfill \\ z_{\pm (k+1)}=z_{\pm k}\pm h(x_{\pm k}{y}_{\pm k}-b_0{z}_{\pm k}),\hfill \end{array}\right.$$

(22)

$x_0=1,y_0=2,z_0=1.$ We assumed that by observing the process described by these equations, inaccurate observations were obtained

$$X_{\pm k}=x_{\pm k}+{\epsilon }_1(\pm hk),Y_{\pm k}=y_{\pm k}+{\epsilon }_2(\pm hk),Z_{\pm k}=z_{\pm k}+{\epsilon }_3(\pm hk),k=0,1,...,n,$$

where ${\epsilon }_i(\pm hk),i=1,2,3,k=0,1,\dots ,n,$ are independent random variables, distributed uniformly over a segment $[-1/2,1/2].$ According to the Formula (19), the parameters were evaluated first $x_0,y_0,z_0,F_i^0=F_i(x_0,y_0,z_0,{\sigma }_0,r_0,b_0),i=1,2,3,$ in our notation ${\widehat{x}}_0,{\widehat{y}}_0,{\widehat{z}}_0,{\widehat{F}}_i^0,i=1,2,3.$ Further, the estimates of the parameters ${\sigma }_0,r_0,b_0$ were found from the relations

$$\left\{\begin{array}{c}{\widehat{F}}_1^0=\sigma ({\widehat{y}}_0-{\widehat{x}}_0),\hfill \\ {\widehat{F}}_2^0={\widehat{x}}_0(r-{\widehat{z}}_0)-{\widehat{y}}_0,\hfill \\ {\widehat{F}}_3^0={\widehat{x}}_0{\widehat{y}}_0-b{\widehat{z}}_0.\hfill \end{array}\right.$$

(23)

Table 3. Interval distribution of estimates ${\widehat{\sigma }}_0,{\widehat{r}}_0,{\widehat{b}}_0$.

distribution intervals ${\widehat{\mathit{\sigma }}}_0$	relative frequencies ${\widehat{\mathit{\sigma }}}_0$
0.883607–0.931473	0.035
0.931473–0.979339	0.275
0.979339–1.0272	0.486
1.0272–1.07507	0.192
1.07507–1.12294	0.015
distribution intervals ${\widehat{r}}_0$	relative frequencies ${\widehat{r}}_0$
1.89817–1.94253	0.038
1.94253–1.98689	0.242
1.98689–2.031262	0.471
2.03126–2.07562	0.224
2.07562–2.11998	0.025
distribution intervals ${\widehat{b}}_0$	relative frequencies ${\widehat{b}}_0$
2.87579–2.92301	0.021
2.92301–2.97022	0.267
2.97022–3.01744	0.457
3.01744–3.06466	0.231
3.06466–3.11188	0.024

shows the results of a computational experiment conducted 1000 times, namely, the interval distribution of relative frequencies ${\widehat{\sigma }}_0,{\widehat{r}}_0,{\widehat{\beta }}_0.$

3. Conclusions

Remarks 3 and 4 indicate the following possible generalizations of the results obtained in Theorems 3–6. First, we should consider the case when the variance of random perturbations ${\sigma }^2$ decreases and so quality of obtained estimates improves. However, if the variance ${\sigma }^2=\infty$ like in a case of heavy tails of disturbances distributions, then it is necessary to consider stable distribution of random variables ${\epsilon }_k.$ Secondly, we should consider the case of ordinary differential Equations (and systems) of higher than the unit order.

Furthermore, at last, along with systems of ordinary differential equations, the proposed method for estimating parameters may be applied to equations or systems of partial differential equations. As a basis for the development of this method of parameter estimation, the theorem of the existence of a solution of partial differential equations system in the vicinity of a certain point may be taken (see, for example, [18]).

4. Discussion

The solution of the considered problem involves the choice of an experimental plan, the use of the theorem of existence and uniqueness for a system of ordinary differential equations, the implicit function theorem and the method of linear regression analysis. Linear regression analysis is based on minimizing of the root-mean-square deviation of the sequence of observations from the linear regression function.

Practically, all the considered generalizations of the results obtained in the paper arise at the junction of several scientific directions. These include probability theory and mathematical statistics, ordinary differential and partial differential equations and their systems, and mathematical analysis. Such tasks arising at the junction of several research directions are usually considered in the system analysis, management and information processing. This circumstance determines the choice of this research topic and the ways to solve the task and an application of optimization procedures.