DATA ANALYSIS - FOURIER SERIES

THEORY PHYSICS MECHANICS MATHEMATICS		OBSERVATIONS EXPERIMENTS IN SITU IN LABORATORY
	COMPUTER SCIENCE NUMERICAL METHODS SIMULATIONS	

Standard formulae

The function y(x) is a periodic function, the period is 2l

(1.1)

.

The complex form of the Fourier series is

(1.2)

The Fourier integral theorem is defined by the relation

(1.3)

The Fourier transform and the inverse are

(1.4)

APPROXIMATION OF FUNCTIONS BY SERIES

Let us approximate y(x) by a linear combination of functions

(2.1)

by the least square method (the square error is minimum),

(2.2)

For a discrete set of data , j = 1,2, ..., r is given by a column matrix Y, the set of coefficients by a column matrix A, i=1,2, ..., n and a rectangular matrix Ф with elements corresponds to values at points , in i columns. In the discrete case the square error is

.	(2.3)

In the case n = r, Ф is a square matrix and if it is not singular, the error is minimum if the expression in the square brackets is zero. It follows A is the solution of the following linear equation

.	(2.4)

In the case r > n we have more equations than unknowns the differentiation with respect to leads to the following expression

.	(2.5)

It corresponds to the following matrix equation

.	(2.6)

The relative error is defined by.

Let us consider the approximation of measured damped free vibrations by the solution of the linear theory. Thus for the case of one degree of freedom the square error

(2.7)

should be minimum. In the integral y(t) is known a and b are displacements of the terms, f is the frequency in Hz and η is the damping coefficient in 1/s. It is convenient to introduce the following notation .

The necessary condition that the square error is minimum corresponds to the condition that all partial derivatives are zero

.	(2.8)

The problem is solved usually by iterations. The choice of estimated initial values is very important to obtain a convergent procedure.

A simple and traditional method of approach is to expand the function under the integral in Taylor series and to consider three terms only.

.	(2.9)

Upon substitution and integration the function E is a function of known (assumed) initial values and the unknown increments . It is a second order power series in the increments. To minimize the value of the square error the partial derivatives with respect to should be zero. This condition leads to four linear equations that, if the determinant is not singular, yield the values of the increments. The new value of the estimates is equal to the sum of the previous estimates and the calculated increments .

If the procedure converges, the step by step procedure leads to small values of the square error, that are satisfactory for the considered problem.

PROBABILITY THEORY AND
RANDOM VARIABLES

The probability space is denoted by and its elements ω are samples or experimental outcomes. Certain subsets (collection of outcomes) are called events Λ.

The set theory notation is used. If A and B are two sets, then is their union, is their intersection, is the complement of B with respect to A, the empty set is denoted by . If , that is the sets are disjoint then the events are mutually exclusive.

For a class of events we assign probabilities to events via a probability function Pr(.). That is, to each event we assign a number Pr(), called the probability of . The probability function satisfies the following axioms:

;

.

The class of events has to be defined. In defining the class of events we want the set operations (unions, intersections, complements) to yield sets that are also events. A class of sets having these properties is called a Borel field. A class F of ω sets is called a Borel field if

;

The triplet is called an experiment.

Example: , Borel field .

Example, rolling a die once . ,

Borel field ,

, is not a Borel field,

because .

RANDOM VARIABLES

A real finite-valued function defined on is called a (real) random variable if, for every real number x, the inequality defines an ω set whose probability is defined. A random variable is a Borel measurable function.

For a random variable the function

,	(3.1)

is defined for all real x and is called the distribution function of the random variable x. A random variable x is called discrete if there exists a mass function such that

.	(3.2)

A random variable is called continuous if there exists a density function such that

.	(3.3)

If the number of points at which is not differentiable is countable then

,	(3.4)

at all x at which the derivative exists.

The expectation, average, mean or first moment of a continuous random variable is defined by

.	(3.5)

The nth moment of x is defined by

.	(3.6)

The second moment is called the mean square value.

The nth central moment of x is define by

.	(3.7)

The second central moment is called the variance of x.

Example: rolling a die. We can define a probability space

, a random variable x(ω), and probability Pr (ω) as for example given in the table

	'1'	'2'	'3'	'4'	'5'	'6'
x(ω)	-30	-20	-10	+10	+20	+30
Pr(ω)	1/6	1/6	1/6	1/6	1/6	1/6

Let us introduce subsets (events, corresponding to odd or even numbers). The corresponding Borel field is and the corresponding probabilities are elements of the row matrix. This is also the case of tossing a coin with two events heads or tails.

Example: choice of a random phase ψ from the interval in continuum, thus. Let us define the random variable and the probability function is

The corresponding distribution function is:

,

,

.

The standard uniform density function results by differentiating the distribution function with respect to the random variable x = ψ.

A random variable x is Gaussian or normally distributed if its density function is given by

,	(3.8)

where (m - mean and - variance).

The normal distribution function is

.	(3.9)

It is convenient to specify a random variable by its characteristic function defined as the Fourier transform of the density function

.	(3.10)

The nth moment of the random variable is

.	(3.11)

It is easy to verify the following relation

.	(3.12)

Thus when the characteristic function is calculated it is easy to calculate the values of the set of nth moments by differentiation.

The characteristic function for a Gaussian density with expectation is

.	(3.13)

( The assumption that m = 0 is not a serious loss of generality because from all samples we can subtract the value of m to obtain a random variable with mean values equal zero.)

Introducing the change of variables in the integral the integration finally yields the following expression

.	(3.14)

Differentiation and substitution into the expressions for the nth central moments lead to the conclusion that for odd values of n the central moments are zero and it is easy to establish a relation between the expression for consecutive even numbers of n. The final result is that for a normal distribution the nth central moment is equal to

.	(3.15)

The even central moments grow without bound when n tends to infinity.

JOINTLY DISTRIBUTED RANDOM VARIABLES

Let us consider two random variables x and y. The two sets and are events with probabilities

(3.16)

where are distribution functions of the random variables x and y. The intersection of these two sets

(3.17)

is an event. The probability of this event is the joint distribution function of the jointly distributed random variables x and y

.	(3.18)

Example: temperature measured x(ω) at 6 at night and y(ω) at 12 (at the same day). It is possible to consider the random variables x and y separately or to look at the pair x(ω), y(ω) jointly. To estimate the joint probability one has to consider a two dimensional problem related to the x, y plane.

In general the continuous random variables defined on the same probability space are said to be jointly distributed. They may be characterized by their joint distribution function

,	(3.19)

where

,	(3.20)

or by their joint density function

.	(3.21)

For the differentiable case

.	(3.22)

The marginal distribution function is defined by

.	(3.23)

The marginal density function is

.	(3.24)

The expectation of , , is given by

.	(3.25)

The second moment of is

.	(3.26)

Of great importance in applications is the covariance of which is defined by

.	(3.27)

The generalization of the higher moments and central moments from the case of one random variable to the joint variables is straightforward.

Two jointly distributed random variables are independent if any of the following equivalent conditions is satisfied

(3.28)

We say that are mutually independent if

.	(3.29)

4 RANDOM SEQUENCES IN MATRIX NOTATION

Let us write the jointly distributed random variables in form of a column matrix, (random vector)

.	(4.1)

Its transpose is a row matrix. The expectation is a column matrix with elements equal to the expectations of the random numbers .

The square matrix with elements is called the covariance matrix of the random vector X and denoted by

,	(4.2)

its elements are

.	(4.3)

The covariance matrix is symmetric. For a symmetric square matrix the eigenvalues are real and the eigenvectors are orthogonal. A matrix A is said to be positive definite if for all vectors . It is easy to see that the matrix has to be positive definite. If the values of the elements of the covariance matrix are estimated from observation it is necessary to transform the matrix to a symmetric form and change the terms so that all the eigenvalues are positive.

In matrix notation the density function of a n jointly normally distributed random variables is

(4.4)

where denotes the determinant and the inverse of the covariance matrix.

Let us now look at some simple examples of Gaussian random sequences

Example 4.1. Let us consider the simple case the expectation is a zero column matrix and the covariance matrix is proportional to the unit matrix I, . Thus it is a diagonal matrix. When substituted into the expression for jointly normally distributed random variables it follows that the elements of the sequence are mutually independent.

It follows that

where U is a Gaussian white noise sequence all elements have the same variance .

Example 4.2. Let us consider the random sequence defined by the following difference equation and initial value

It means the sequence has independent increments with equal variances.

The covariance matrix has the following form

It is easy to verify that the elements of the covariance matrix are given by the following expression

.

Let us generalize the results of Example 4.2 to the of not unit interval but of length . This results in a change in the notation of the coefficient. This leads to the following form of the difference equation:

,	(4.5)

and the expression for the element of the covariance matrix becomes

.	(4.6)

This random difference equation may be used to study the case when . The result is that the sequence tends to a continuous function with no derivative in any point.

Example 4.3. Let us consider the random sequence defined by the following difference equation and initial value

where is an element of a white noise sequence and the condition that the elements have equal variances has to be satisfied. The condition yields the relation

It is easy to verify that the covariance matrix has the following form

It is easy to verify that the general expressions for the elements of the covariance matrix are

.

The covariance matrix has same values on all diagonals. The values depend upon the distances of the points .

Let us generalize the results of the Example 4.3 by introduction of the following change of notation in parameters in the expressions for the elements of the covariance matrix It follows

.	(4.7)

In the new notations the random difference equation becomes

(4.8)

This form is suitable to study the influence of the value of on the behaviour of the solution. The final result is: when the difference equation tends to an Itô random differential equation with a solution that is a continuous function with no derivative at any point.

5 STOCHASTIC PROCESSES,
MEAN SQUARE CALCULUS

Let be a stochastic process. For any finite set, the joint distribution function of the random variables is called a finite dimensional distribution of the process. The stochastic process can be characterized by specifying the joint density function

.

Classification of Stochastic Processes

Continuous state space	Random sequence	Stochastic process Random function
Discrete state space	Discrete parameter chain	Continuous parameter chain
	Discrete parameter set	Continuous parameter set

Example 5.1. At times 0,1,2,… we toss a fair coin. For each time we define a random variable

We assume that the start is at and the random variable is defined as

.

It is clear that is a discrete parameter chain.

The most important are the first-order and the second order densities .

The mean value function (expected value) is defined as

(5.1)

The (auto) correlation function is defined as

.	(5.2)

The (auto) covariance function is defined as

.	(5.3)

STATIONARY PROCESSES

A process is strictly stationary if it has the same probability laws

(5.4)

for all finite sets and for any .

The process is strictly stationary of order k if it holds for . For k=1 the first order density is independent of time

.	(5.5)

Thus the expected value is constant is constant.

The second order density can depend only on , thus the correlation and covariance functions are

(5.6)

The stochastic process is said to be weakly stationary (or stationary in the wide sense, or covariance stationary) if it has finite second moments, the mean value function is constant and the correlation function is a function of distance τ.

The stochastic process is said to have strictly stationary increments if the process is strictly stationary for every .

Convergence of Random Sequences

There are a number of ways in which a sequence may converge as .

The sequence is said to converge to x with probability 1 if

(5.7)

for almost all realizations ω. (Does not hold perhaps for event A with probability Pr(A)=0.

The sequence is said to converge to x in probability if, for every

(5.8)

The sequence is said to converge to x in mean square if

(5.9)

We write .

The Cauchy criterion

,	(5.10)

is a necessary and sufficient condition for mean square convergence.

Mean Square Calculus

The random function is said to be continuous in mean square at if

.	(5.11)

The correlation function of the random function is . Let us introduce a random function and calculate the cross correlation function

,

and then the correlation function

.

Finally upon substitution

.

It follows that for

,	(5.12)

and thus the random function is mean square continuous if the correlation function is continuous at .

Theorem. The random function is mean square continuous at if, and only if is continuous at (t,t).

Let us consider a derivative of a stochastic process.

The problem is simple when we know how to calculate the realizations of the stationary process with a mean value and a covariance function . If we can calculate the derivatives of all realizations, then the problem of the derivative is the problem of derivatives of a family of deterministic functions.

In view of the Cauchy criterion a sufficient condition for mean square differentiability of is the convergence to zero with h, s of

.

The final result is the theorem.

Theorem. The random function is mean square differentiable at if, and only if exists at (t,t).

For stationary processes it is easier to prove similar theorems.

A stationary process is mean square continuous at if, and only if is continuous at .

A stationary process is mean square differentiable at if, and only if exist at .

For example let us consider the stationary random process with the following covariance function

.

The covariance function is continuous at τ=0, thus the random process is continuous. The covariance function is not differentiable for . (The right side and left side derivatives are different.). The second derivative does not exist, the random function is not differentiable.

Let us consider a Riemann integral of a stochastic process.

We introduce a partition of the interval ,

We denote .

The random function is mean square Riemann integrable if the following limit, which then defines the integral, exists

.	(5.13)

Theorem. The random function is mean square Riemann integrable over [a,b] only if, and if , is Riemann integrable over [a,b]x[a,b] (The Riemann integral exist).

If a random function is mean square integrable over [a,b] for every , then

(5.14)

is a random function of t defined on [a,b].

Theorem. The random function is mean square integrable over [a,b] for every . Then is mean square continuous on [a,b]. If is ms continuous on [a,b], then is ms differentiable on and .

The Fundamental Theorem on ms Calculus.

Let be ms Riemann integrable on [a,b]. Then

(5.15)

with probability one.

The Brownian Motion Process

A random function has independent increments if, for all finite sets the random variables (vectors)

are independent.

The process has stationary independent increments if, in addition has the same distribution as for every and every h.

A random function is a Brownian motion (Wiener or Wiener-Lévy) process if

It follows from (ii) that is also normally distributed for every It remains to specify the distribution of increments for all The mean, in view of (iii) is zero. From the definition it follows:

.

The function does not decrease when t increases. Thus the equation

has a solution

.	(5.16)

It is the only solution.

Finally the mean value of the increment is zero and its variance is .

Let us calculate the correlation (the covariance is the same) function. For

.

Therefore in general for the Brownian motion process the correlation function is

.	(5.17)

The correlation function is continuous for every. Thus the process is mean square continuous on. The second derivative exists at no, thus the process is ms differentiable nowhere. It is Riemann integrable on every t interval.

To compute the values we have to chose the time increment and and to write the covariance function in terms of the increment

.	(5.18)

This leads to the following form of the difference equation:

,	(5.19)

where is a term of the Gaussian white noise sequence with mean value zero and unit variance. The expressions (5.18) and (5.19) reduce the problem in continuum to the solution of a difference equation as described in chapter 4 (compare with the relations (4.6)and (4.5)).

It should be noted that , and this expression is an invariant with respect to the choice of the increment , and not the relation for the derivative in finite differences.

6 STOCHASTIC DIFFERENTIAL EQUATIONS

Let us start our consideration with the case of linear stochastic differential equation. Let us consider the following first order, homogeneous differential equation

.	(6.1)

Looking at the realizations it is a deterministic differential equation

.	(6.2)

The randomness is in the initial condition. Let us assume that the initial condition is a normally distributed random variable with mean value m = 0 and variance equal to . The solution of the differential equation is

.	(6.3)

The process is Gaussian. The mean value and correlation functions are

.	(6.4)

The variance is

.	(6.5)

As a second example let us consider the following set of two differential equations

.	(6.6)

The randomness is in the initial conditions. Let us assume that the initial condition for is a normally distributed random variable with mean value and variance equal to . The initial conditions for are independent, normally distributed; the mean value is , and the variance is . The solution of the differential equation is

.	(6.7)

The mean values of the stochastic processes are

,

and the elements of the correlation matrix are

(6.9)

Let us consider the following first order, not homogeneous differential equation with the initial value :

.	(6.10)

On the right side is a linear function in time with random independent coefficients that are normally distributed with zero mean values and variances. The realizations have to be computed from the corresponding deterministic equation:

.	(6.11)

The solution for a realization is

.	(6.12)

The solution for the random case (for the family of realizations) is

(6.13)

The mean value and the correlation functions are

.	(6.14)

The stochastic process is continuous and differentiable. The second mixed partial derivative is a constant, thus all derivatives exist.

7 ITÔ DIFFERENTIAL EQUATION

A very important special case is the linear stochastic differential equation

(7.1)

where are deterministic functions and is the differential of the Brownian motion process. This is a first order differential equation with an additive noise. A more general case is when the function g(t) is random , but then we have a product of random functions. The integral form of the differential equation is

,	(7.2)

where the first integral is the Riemann stochastic integral and the second in general is an Itô stochastic integral and in the above case

is a Wiener stochastic integral.

Let us introduce the Itô stochastic integral by plausible arguments. Rigorous proves need a lot of advanced mathematics. The stress will be on applications and the use of numerical methods. We are going to discuss the Itô stochastic integral

.	(7.3)

Assume is a scalar Brownian motion process with variance parameter and the function is defined on . Let us introduce a partition

, and a step function

(7.4)

where is independent of , and . The Itô integral is defined by

.	(7.5)

In view of the independence condition the expected value is zero

.	(7.6)

Let us take another function and the corresponding step function and consider the expected value

(7.7)

where .The last expression results from the independence conditions.

Now let us consider a sequence of step functions converging to the random function in the sense that

.	(7.8)

Finally the following theorem is true:

Let the random functions satisfy the condition that they are independent of , for all , and the conditions

.

Then their Itô integrals are well defined as

(7.9)

with basic properties

,

.

Let us consider the following linear stochastic differential equation with constant coefficients

(7.10)

where are constants and is the differential of the Brownian motion process with variance parameter equal to one. This is a first order differential equation with an additive noise. The general solution of the homogeneous differential equation is

,	(7.11)

where is the random initial value for . Addition of a particular solution of the non homogeneous equation yields the general solution

,	(7.12)

where the integral on the right side is an Itô (Wiener) integral. Let us assume that the initial value is normally distributed with a mean value equal to zero. The mean value of the Itô integral is zero. Thus the mean value of the random function is zero. The correlation (covariance) function is

(7.13)

From the basic formula for the Itô integrals it follows

.	(7.14)

Upon integration it follows

,	(7.15)

where .

When tends to infinity in such a way that is constant, then in the limit

(7.16)

The process is asymptotically stationary in the wide sense and because it is Gaussian it is strictly asymptotically too. If the variance of the initial condition C(0) is equal to P the process is strictly stationary.

In general the process is not stationary. The variance of the process is

.	(7.17)

The variance changes from at the time to the asymptotic value P when t tends to infinity. The derivative with respect to time is

.	(7.18)

It is easy to show that the differential equation for the evolution of the variance in our example is

.	(7.19)

8 LINEAR MATHEMATICAL MODEL WITHOUT
DOMINANT FREQUENCY

Let us consider the following set of differential equations with constant coefficients η and ω and the variance parameter σ = 1.

(8.1)

The first equation is a stochastic Itô differential equation. Its general solution may be easily written. The solution is a continuous Riemann integrable random function. The second equation is a stochastic differential equation. The solution is a once differentiable Riemann integrable random function. Finally the function is n times differentiable.

It is convenient to write the relations in a matrix notation

(8.2)

where are (n+1) column matrices

and ψ is a lower triangular matrix

The fundamental solution is satisfying the corresponding homogeneous equation and initial conditions

(8.3)

is expressed by the following matrix

The general solution is

.	(8.4)

The mean value vector satisfies the following equation

.	(8.5)

Let us denote the variance matrix of as

.	(8.6)

From the differential equation in matrix notation it follows

, and thus

.

Finally the evolution of is described by the following differential equation

(8.7)

(The third term on the right side is due to the Itô integral; it gives a contribution in the first equation only.)

If the asymptotic covariance matrix exists it is a solution of the following algebraic equation.

(8.8)

The matrix has the following structure

(8.9)

The solution of the differential equation for the evolution of the variance is

.	(8.10)

The final expression for the covariance matrix is

(8.11)

For example in the stationary case the expressions are

The spectral density function of a stationary in the wide sense is the Fourier transform of the autocorrelation function

,	(8.12)

where ω is the angular frequency.

The inverse transformation is

.	(8.13)

For example the spectral densities for the non, once and twice differentiable processes are

Let us discuss how to calculate realizations of the above considered process in a recursive formulation. For the considered stationary process, we divide the time T into equal time intervals , and we want to express the column matrix in terms of the values of . The process is stationary and thus we may consider one interval from 0 to . It is easy to verify that in our case the general solution

(8.14)

for one step from may be written in the following form

.	(8.15)

The expected value of a stationary process is constant and thus without loss of generality it may be assumed equal to zero. For a stationary process the variances for must have the same values. It follows due to the independence conditions that

.	(8.16)

All the terms correspond to symmetric matrices thus they may be represented by where is a lower triangular matrix.

.	(8.17)

The matrix may be calculated from the chol(M ) procedure or directly by the following representation (given for a 3x3 matrix M)

The matrix on the right side is symmetric and thus has 6 different elements. The matrices on the left have 6 unknown elements. First we multiply the second and third matrices. Then we multiply the first row by the first column and obtain the value of . Multiplication by the second and third columns leads to the values of . Then we multiply the second row by the second column and obtain the value of . Multiplication by the third columns yields the value of . Finally the multiplication of the third row by the third column leads to the last unknown value. It may happen that is zero. It means the matrix M is singular but this matrix corresponds to a covariance matrix and must be positive definite. In such a case we have to reduce the number of elements in the column matrix of the white noise sequence.

Finally the values of the realizations in one step may be calculated from the following recursive equation

,	(8.18)

where is a column matrix with Gaussian independent random numbers in n+1 rows.

To obtain a stationary series the initial conditions must correspond to jointly normally distributed random numbers with mean values equal to zero and covariance matrix equal to the asymptotic variance matrix. To compute the initial conditions it is convenient to represent the asymptotic variance matrix by the product of a lower triangular matrix p by its transpose

.	(8.19)

For example for a twice differentiable process the matrix p is

,

where U is a column matrix with Gaussian independent random numbers in three rows.

9 LINEAR MATHEMATICAL MODEL WITH
DOMINANT FREQUENCY

Let us compute two independent realizations of the stationary process described in the preceding paragraph and form a complex random function

.	(9.1)

Let us introduce a complex random function by the following relation

,	(9.2)

whereis a dominant angular frequency. The real and imaginary parts define two real random functions

(9.3)

The mean values of are zero.

The covariance functions of the real functions are

(9.4)

The function is even and the function is odd. They do not depend upon the value of t and therefore the processes are stationary in the wide sense and because they are Gaussian they are strictly stationary. The cross covariance function for τ = 0 is zero, thus the random variables are statistically independent.

The spectral density is the Fourier transform of the covariance function. Thus

(9.5)

For a real function the spectral density is defined in the interval and is an even function. Thus when only positive values of ω are considered it follows

.	(9.6)

The complex random function may be written in an exponential form

.	(9.7)

The absolute value is equal to

,	(9.8)

and the phase shift may be calculated from the relations

.	(9.9)

It should be noted that if at a time t the function then the curve is tangent to the curve . It may be easily seen that the functions are envelopes for both functions and .

The amplitude of the discussed random function with a dominant frequency is changing in time and due to the random phase shift the local angular frequency is changing in time too. For example if in an interval the phase shift may be approximated by a linear function

the local angular frequency is . Thus the distance between the down or up crossings is a random sequence.

It may be easily verified that the random variable has a Rayleigh distribution and the random variable a uniform distribution on the interval of length 2π. These random variables are independent. Thus

(9.10)

and the joint probability function is equal to the product of these functions.

The differential equations for the functions are not written in a suitable form, they do not correspond to differential equations with constant coefficients. The first two equations that correspond to two independent processes without dominant frequencies may be written in matrix notations

,	(9.11)

where are independent Brownian motion processes.

For example for the first two equations the differentials are

and thus

.	(9.12)

The first matrix on the right side is an orthogonal and normal matrix (its inverse is equal to the transpose and the determinant is equal one). Such a matrix represents rotations and will be denoted by . Multiplication of the initial equation by the orthogonal matrix and upon substitution yields the following final equation

,	(9.13)

where the property was used that an orthogonal and normal transformation of two independent increments of Brownian motion preserves their properties. Similar relations

,	(9.14)

hold for the other differential equations in the set that correspond to standard differential equations.

Let us consider a case of a twice differentiable function in matrix notations.

(9.15)

This set of linear differential equations has constant coefficients and therefore it is easy to solve it by standard methods.

Let us look at the fundamental solution of the homogeneous equation (α = 0). The solution for the first matrix differential equation

by the standard method with initial conditions at t = 0 is

.	(9.16)

The second matrix differential equation

has a general solution that is the sum of a general solution ( similar as in the previous case) and a particular solution of the non homogeneous case (the right side is a known matrix). It follows

.	(9.17)

The same simple procedure leads to the general solution of the third matrix homogeneous differential equation

.	(9.18)

If we denote by the matrix

,

the general solution may be written in the form of a block matrix

,	(9.19)

where 0 is a 2x2 matrix with elements equal to zeros.

The asymptotic variance matrix has the following structure in block matrix notation

(9.20)

To simulate a stationary process the initial conditions for t = 0 should be computed with the help of a lower triangular matrix p that satisfies the relation . For a twice differentiable function, in block matrix notation the matrix p is

,	(9.21)

where and is a column matrix with Gaussian independent random numbers in three rows.

The stationary random series may be computed from the following recurrence equation

,	(9.22)

where is computed from the relation

.	(9.23)

It should be noted that when the block matrix notation is used in P and φ it follows that the block matrix multiplication leads to the following relation for the case of a twice differentiable function

,	(9.24)

where are elements of the matrix for the case of the corresponding random process without a dominant frequency.

10 KALMAN FILTERS - BASIC RELATIONS

The introduction to the Kalman filter theory will be limited to the linear theory. The stress is on the applications. Thus relations will be justified by plausible arguments based on physics of the problem. The stress is not on rigorous mathematical proofs. The dynamical system is described by a set of linear Itô differential equation

(10.1)

where is the column matrix of n random functions, are respectively n x n and n x m deterministic, continuous matrix functions and is a column matrix of m independent Brownian motion processes with variance parameters equal to ones. The initial conditions correspond to Gaussian processes. In the second row the dimensions of the matrices are given for a twice differentiable process with a dominant frequency.

Discrete, linear observations are taken at times

,	(10.2)

where is the column matrix of s random observations, is a deterministic matrix and is the column matrix of s Gaussian white noise sequences with zero mean values and prescribed variances. It is assumed that the initial conditions, the Brownian motion processes and the white noise observations sequences are independent.

The basic theorems are true for all linear systems. For constant coefficients the solutions are obtained by standard methods that are very easy to apply. In practical applications the computations are often restricted to discrete set of times with constant spacing. The most important case in applications is a dynamical system described by a set of difference equations with equal time intervals

,	(10.3)

where the column matrix with n rows contains the values of the mathematical model, the matrix corresponds to a transformation in one step, is an matrix and is an Gaussian white noise sequence with mean values equal to zeros and variances equal to one (the variance may be included in the values of the elements of q). The observation model is the same as in the case of the continuous discrete model

.	(10.4)

In the example of a twice differentiable process with a dominant frequency the matrix X and h are

(10.5)

For two events A and B the conditional probability function of event A given the event B is

.	(10.6)

Let us now introduce an outline of the procedure in the Kalman filter formulations. The discussion will be based on the example of a twice differentiable discrete random process and the corresponding observation model.

Let us assume that at the time we have a good estimate of the investigated realization of the random vector denoted by based on all observations up to r. Now we want to compute the expected value at the time on the data from the observation model for the time intervals . This value may be calculated from the mathematical model as the expectation of the conditional probability as defined in the relation (10.6). Thus

.	(10.7)

Let us now compute the prediction of the variance of the error

The final result is

(10.8)

Up till now we discussed the predictions for the time

based on the mathematical model and the observation up till .

To compute the estimate the observation at the time has to be taken into account. It was proved by Kalman that the influence of this observation leads to the following relation

,	(10.9)

where K is the Kalman gain given by the following relation

.	(10.10)

where R is a diagonal matrix of variances of the matrix white noise sequence V(r) in the relation (10.4).

From the physics of the problem it is clear that the prediction has to be supplemented by the expected value of the observation error. In our example there is only one observation and it corresponds to the measurement of . In general the matrix R is a set of numbers that has to be related to the variances of the measured random processes. For example the standard deviation is 5% of the standard deviation of .

The procedure must be completed by the estimation of the variance of the prediction error when the observation at the time is taken into account

.	(10.11)

This procedure leads to the computation of all the necessary initial values for the next step . Now it is possible to follow the outlined procedure to calculate the values in the next step: Taking consecutive steps the computation of the realization may be completed. It is only necessary to have the initial values at the initial time t = 0. If there is no information available it is reasonable to assume. It means the initial values are equal to the expected values and the variance of the initial error is equal to the asymptotic variance of the process.

Tego nie damy i nie damy przykładu pwsemk04.

Now let us discuss how to estimate the initial conditions that is for example for a twice differential process with a dominant frequency the six values .

A simple case corresponds to an example of free linear oscillations described by the differential equation

,

where m is the mass, ε the coefficient of damping, the natural angular frequency, k the spring constant and the general solution is

where A, D are constants which are fixed by the initial conditions and is the angular frequency of free linear damped vibrations We may introduce a solution by a complex signal

.

It is easy to verify that the real and the complex function satisfy the differential equation. The first and second derivatives and their initial values respectively are

Now we have to look at the corresponding solution of the mathematical model with the corresponding initial conditions. We consider a twice differential realization of a random function and it is convenient to use the expression in matrix notation. ( Tu jest powtórka powołać sie na nr)

.

In computations of derivatives it has to remembered that is a orthogonal and normal matrix with properties

.

The comparison with the solution of the differential equation for damped vibrations shows that when we take a good estimate requires that . Thus we take into account in the observation model the given (the solution of the differential equation) and additionally . Thus the observations are expressed by a matrix .

1. Fourier Series

1 The script file pwopro03 gives examples of Fourier series expansions of periodic functions. The first example is a discontinuous skew symmetric function, the second a continuous function with a discontinuous derivative.

2 The script file pwoprs03 calculates the Fourier series coefficients for damped free vibrations looks at the influence of the number of points. The values of coefficients change but the amplitudes are similar.

2. Approximation

1 The script file pwoprp03 calculates the Fourier series of a periodic soliton as a function of time for different shifts in space x and compares with solutions of approximations by trigonometric functions. The solutions when the orthogonal properties of the set of functions are not used are better.

2 The script file pwoprq03 calculates the Fourier series of a periodic soliton and discusses the problem of taking s points less compared with the perfect period. The approximation is not sensitive to such errors.

3 The script file pwoprr03 calculates the Fourier series of a periodic soliton and discusses the problem when additionally two and three periods are considered. Such a change does not change a solution for a perfect period.

4 The script file pwoprt03 approximates damped free vibrations by the expression for linear not damped oscillations and damped vibrations. The iterations converge.

5 Script file pwoprv03 is used to study the measured damped vibrations. The column matrix of measured values is denoted by Y, its elements by Y(r). We introduce the following column matrices E(r)=exp(-x3*r*dt), C(r)=cos(2*pi*x4*r*dt),S(r)=sin(2*pi*x4*r*dt), R(r)=r*dt; We introduce the following notations: EC=E.*C, ES=E.*S, REC=R.*EC, RES=R.*REC, RREC=R.*REC, RRES=R.*RES. The column matrix of errors is EY=Y-x1*EC-x2*ES, the mean square error is a scalar EJ=EY'*EY; The variance of Y is V=mean(Y'*Y) and thus the relative error is RE=sqrt(EJ/V).

3. Probability

4. Random sequences in matrix notation

1 The script file pwsema04 gives examples of simple random sequences: 1. White Gaussian Sequence, 2. Brownian Motion Sequence, 3. not differentiable stationary process.

5. Stochastic processes, mean square calculus

6.Stochastic differential equations

7. Itô differential equations

8. Linear mathematical model without dominant frequency

1 The script file pwsemb04 calculates examples of once and twice differentiable processes with no dominant frequency

2. The script file pwsemd04 depicts the correlation functions and the spectral densities of the not, once and twice differentiable processes

9. Linear mathematical model with dominant frequency

1 The script file pwsemc04 calculates examples of once and twice differentiable processes with dominant frequency, envelopes and derivatives.

2. The script file pwsemf04 depicts the correlation functions and the spectral densities of the non, once and twice differentiable processes.

3 The script file pwsemg04 computes examples of twice differentiable realizations with a dominant frequency.

10. Kalman Filters

1 The script file pwsemh04 computes examples of twice differentiable realizations with a dominant frequency computes matrices for Kalman filters and applies them.

2 The script file pwsemi04 is a program to discuss a random function and the estimated corresponding random sequence.

3 The script file pwsemj04 is a program to discuss a random function and the corresponding estimated random sequence and derivatives.

36

FOURIER SERIES

Your User Name

2004-01-02T17:44:00