Section: Scientific Foundations

Fast parametric estimation and its applications

Parametric estimation may often be formalized as follows:



  • the measured signal y is a functional F of the "true" signal x, which depends on a set Θ of parameters,

  • n is a noise corrupting the observation.

Finding a "good" approximation of the components of Θ has been the subject of a huge literature in various fields of applied mathematics. Most of those researches have been done in a probabilistic setting, which necessitates a good knowledge of the statistical properties of n. Our project is devoted to a new standpoint which does not require this knowledge and which is based on the following tools, which are of algebraic flavor:

  • differential algebra(Differential algebra was introduced in nonlinear control theory by one of us almost twenty years ago for understanding some specific questions like input-output inversion. It allowed to recast the whole of nonlinear control into a more realistic light. The best example is of course the discovery of flat systems which are now quite popular in industry.), which plays with respect to differential equations a similar role to commutative algebra with respect to algebraic equations;

  • module theory, i.e., linear algebra over rings which are not necessarily commutative;

  • operational calculus which was the most classical tool among control and mechanical engineers(Operational calculus is often formalized via the Laplace transform whereas the Fourier transform is today the cornerstone in estimation. Note that the one-sided Laplace transform is causal, but the Fourier transform over R is not.).

Linear identifiability

In most problems appearing in linear control as well as in signal processing, the unknown parameters are linearly identifiable: standard elimination procedures are yielding the following matrix equation

Pθ 1 θ r =Q,(2)


  • θ i , 1ir, represents unknown parameter,

  • P is a r×r square matrix and Q is a r×1 column matrix,

  • the entries of P and Q are finite linear combinations of terms of the form t ν d μ ξ dt μ , μ,ν0, where ξ is an input or output signal,

  • the matrix P is generically invertible, i.e.,  det(P)0.

How to deal with perturbations and noises?

With noisy measurements equation (2 ) becomes:

Pθ 1 θ r =Q+R,(3)

where R is a r×1 column matrix, whose entries are finite linear combination of terms of the form t ν d μ η dt μ ,μ,ν0, where η is a perturbation or a noise.

Structured perturbations

A perturbation π is said to be structured if, and only if, it is annihilated by a linear differential operator of the form  finite a k (t)d k dt k , where a k (t) is a rational function of t, i.e.,   finite a k (t)d k dt k π=0. Note that many classical perturbations like a constant bias are annihilated by such an operator. An unstructured noise cannot be annihilated by a non-zero differential operator.

By well known properties of the non-commutative ring of differential operators, we can multiply both sides of equation (3 ) by a suitable differential operator Δ such that equation (3 ) becomes:

ΔPθ 1 θ r =ΔQ+R ' ,(4)

where the entries of the r×1 column matrix R ' are unstructured noises.

Attenuating unstructured noises

Unstructured noises are usually dealt with stochastic processes like white Gaussian noises. They are considered here as highly fluctuating phenomena, which may therefore be attenuated via low pass filters. Note that no precise knowledge of the statistical properties of the noises is required.


Although the previous noise attenuation(It is reminiscent to what most practitioners in electronics are doing.) may be fully explained via formula (4 ), its theoretical comparison(Let us stress again that many computer simulations and several laboratory experiments have been already successfully achieved and can be quite favorably compared with the existing techniques.) with today's literature(Especially in signal processing.) has yet to be done. It will require a complete resetting of the notions of noises and perturbations. Besides some connections with physics, it might lead to quite new "epistemological" issues  [80] .

Some hints on the calculations

The time derivatives of the input and output signals appearing in equations (2 ), (3 ), (4 ) can be suppressed in the two following ways which might be combined:

  • integrate both sides of the equation a sufficient number of times,

  • take the convolution product of both sides by a suitable low pass filter.

The numerical values of the unknown parameters Θ=(θ 1 ,,θ r ) can be obtained by integrating both sides of the modified equation (4 ) during a very short time interval.

A first, very simple example

Let us illustrate on a very basic example, the grounding ideas of the algebraic approach, based on algebra. For this, consider the first order, linear system:

y ˙(t)=ay(t)+u(t)+γ 0 ,(5)

where a is an unknown parameter to be identified and γ 0 is an unknown, constant perturbation. With the notations of operational calculus and y 0 =y(0), equation (5 ) reads:

sy ^(s)=ay ^(s)+u ^(s)+y 0 +γ 0 s(6)

where y ^(s) represents Laplace transform.

In order to eliminate the term γ 0 , multiply first the two hand-sides of this equation by s and, then, take their derivatives with respect to s:

d dsssy ^(s)=ay ^(s)+u ^(s)+y 0 +γ 0 s(7)

2sy ^(s)+s 2 y ^ ' (s)=asy ^ ' (s)+y ^(s)+su ^ ' (s)+u ^(s)+y 0 .(8)

Recall that y ^ ' (s)dy ^(s) ds corresponds to -ty(t). Assume y 0 =0 for simplicity's sake(If y 0 0 one has to take above derivatives of order 2 with respect to s, in order to eliminate the initial condition.). Then, for any ν>0,

s -ν 2sy ^(s)+s 2 y ^ ' (s)=s -ν a(sy ^ ' (s)+y ^(s))+su ^ ' (s)+u ^(s).(9)

For ν=3, we obtained the estimated value a:

a=2 0 T dλ 0 λ y(t)dt- 0 T ty(t)dt+ 0 T dλ 0 λ tu(t)dt- 0 T dλ 0 λ dσ 0 σ u(t)dt 0 T dλ 0 λ dσ 0 σ y(t)dt- 0 T dλ 0 λ ty(t)dt(10)

Since T>0 can be very small, estimation via (10 ) is very fast.

Note that equation (10 ) represents an on-line algorithm that only involves two kinds of operations on u and y: (1) multiplications by t, and (2) integrations over a pre-selected time interval.

If we now consider an additional noise, of zero mean, in (5 ), say:

y ˙(t)=ay(t)+u(t)+γ 0 +n(t),(11)

it will be considered as fast fluctuating signal. The order ν in (9 ) determines the order of iterations in the integrals (3 integrals in (10 )). Those iterated integrals are low-pass filters which are attenuating the fluctuations.

This example, even simple, clearly demonstrates how algebraic's techniques proceed:

  • they are algebraic: operations on s-functions;

  • they are non-asymptotic: parameter a is obtained from (10 ) in finite time;

  • they are deterministic: no knowledge of the statistical properties of the noise n is required.

A second simple example, with delay

Consider the first order, linear system with constant input delay(This example is taken from  [69] . For further details, we suggest the reader to refer to it.):

y ˙(t)+ay(t)=y(0)δ+γ 0 H+bu(t-τ).(12)

Here we use a distributional-like notation where δ denotes the Dirac impulse and H is its integral, i.e., the Heaviside function (unit step)(In this document, for the sake of simplicity, we make an abuse of the language since we merge in a single notation the Heaviside function H and the integration operator. To be rigorous, the iterated integration (k times) corresponds, in the operational domain, to a division by s k , whereas the convolution with H (k times) corresponds to a division by s k /(k-1)!. For k=0, there is no difference and H*y realizes the integration of y. More generally, since we will always apply these operations to complete equations (left- and right-hand sides), the factor (k-1)! makes no difference.). Still for simplicity, we suppose that the parameter a is known. The parameter to be identified is now the delay τ. As previously, γ 0 is a constant perturbation, a, b, and τ are constant parameters. Consider also a step input u=u 0 H. A first order derivation yields:

y ¨+ay ˙=ϕ 0 +γ 0 δ+bu 0 δ τ ,(13)

where δ τ denotes the delayed Dirac impulse and ϕ 0 =(y ˙(0)+ay(0))δ+y(0)δ (1) , of order 1 and support {0}, contains the contributions of the initial conditions. According to Schwartz theorem, multiplication by a function α such that α(0)=α ' (0)=0, α(τ)=0 yields interesting simplifications. For instance, choosing α(t)=t 3 -τt 2 leads to the following equalities (to be understood in the distributional framework):

t 3 [y ¨+ay ˙]=τt 2 [y ¨+ay ˙],bu 0 t 3 δ τ =bu 0 τt 2 δ τ .(14)

The delay τ becomes available from k1 successive integrations (represented by the operator H), as follows:

τ=H k (w 0 +aw 3 ) H k (w 1 +aw 2 ),t>τ,(15)

where the w i are defined, using the notation z i =t i y, by:

w 0 =t 3 y (2) =-6z 1 +6z 2 (1) -z 3 (2) ,w 1 =t 2 y (2) =-2z 0 +4z 1 (1) -z 2 (2) ,w 2 =t 2 y (1) =2z 1 -z 2 (1) ,w 3 =t 3 y (1) =3z 2 -z 3 (1) .

These coefficients show that k2 integrations are avoiding any derivation in the delay identification.

Figure 1. Delay τ identification from algorithm (15 )

Figure 1 gives a numerical simulation with k=2 integrations and a=2,b=1,τ=0.6, y(0)=0.3,γ 0 =2,u 0 =1. Due to the non identifiability over (0,τ), the delay τ is set to zero until the numerator or the denominator in the right hand side of (15 ) reaches a significant nonzero value.

Again, note the realization algorithm (15 ) involves two kinds of operators: (1) integrations and (2) multiplications by t.

It relies on the measurement of y and on the knowledge of a. If a is also unknown, the same approach can be utilized for a simultaneous identification of a and τ. The following relation is derived from (14 ):

τ(H k w 1 )+aτ(H k w 2 )-a(H k w 3 )=H k w 0 ,(16)

and a linear system with unknown parameters (τ,aτ,a) is obtained by using different integration orders:

H 2 w 1 H 2 w 2 H 2 w 3 H 3 w 1 H 3 w 2 H 3 w 3 H 4 w 1 H 4 w 2 H 4 w 3 τ ^aτ ^-a ^=H 2 w 0 H 3 w 0 H 4 w 0 .

The resulting numerical simulations are shown in Figure 2 . For identifiability reasons, the obtained linear system may be not consistent for t<τ.

Figure 2. Simultaneous identification of a and τ from algorithm (16 )