# A New Method for Parameter Sensitivity Analysis of Lorenz Equations

Authors: Mehmet Ali Akinlar
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

A new method for parameter sensitivity analysis of Lorenz equations is presented. The sensitivity equations are derived based on the staggered methods. Experimental results indicate that it is possible to determine effects of parameters on model variables so that we can eliminate the less effective ones. Robustness can also be verified in some confidence intervals by simply looking at the corresponding phase portraits. This enables us to control the system. Although the stability properties of the Lorenz equations are studied extensively, to the best knowledge of the authors, the PSA of Lorenz equations has not been considered which is the main goal of this paper.


---

## Body

## 1. Introduction

Parameter sensitivity analysis (PSA) of large-scale differential algebraic systems is important in many engineering and scientific applications, including biology, chemistry, and economics. Problems such as population dynamics, network modeling, and chemical reactors coming from different branches of science have many parameters whose values may not be known accurately. Infinitesimal changes in most of these model input parameters change the future behavior of the systems partially or completely. Consequently, one can observe an uncontrolled and chaotic behavior of the system. In the present day, one has an opportunity to adjust these parameter values accordingly and make some list of parameters with respect to their effect on model. For instance, if a parameter is less effective than the other parameters, the designer of the model can eliminate that parameter. The analysis of this effectiveness is calledparameter sensitivity analysis. Consequently, algorithms which perform PSA in an efficient and rapid manner are invaluable to researchers in many fields.In this paper, a new method for parameter sensitivity analysis of Lorenz equations is presented. The sensitivity equations are derived based on the staggered methods. Experimental results indicate that it is possible to determine effects of parameters on model variables so that we can eliminate the less effective ones. Robustness can also be verified in some confidence intervals by simply looking at the corresponding phase portraits. This enables us to control the system. Although the stability properties of the Lorenz equations are studied extensively, to the best knowledge of the authors, the PSA of Lorenz equations has not been considered which is the main goal of this paper.The structure of this paper is as follows. In Section2, we overview the concept of parameter sensitivity analysis. In Sections 3 and 4, we study the chaotic and sensitivity analysis of Lorenz equations. We complete the paper by some simulation results.


## 2. Parameter Sensitivity Analysis

It is difficult to construct a model without any parameter. In fact, in reality, the problems coming from different branches of science such as engineering, biology, ecology, and meteorology have many parameters. With the help of faster computers of today, one has a chance to adjust them and make some list of parameters with respect to their effect on model. If a parameter is less effective than the others we, the designer of the model, can eliminate it. The analysis of this effectiveness is called “sensitivity.” When qualitative estimates of sensitivity are desired, a mathematical model of phenomena is desired or at least a relationship is required. Infinitesimal changes in all (or some) of the model input parameters change the future of mathematical design partially or completely (in some cases). The important thing here is the sensitivity of a single component compared to other input variables changed a little bit simultaneously. By saying single component we mean the parameters in the model whose values may not be accurately known. However, such a model brings the questions concerning stability, optimality, sensitivity, and so forth. In this work we concentrated on only the sensitivity analysis of a concrete example, namely, Lorenz equations.PSA generates essential information for parameter estimation, optimization, control, model simplification, and experimental design. In the literature, staggered direct method, simultaneous corrector method, adjoint method, and staggered corrector method are some of the well-known methods for parameter sensitivity analysis. We can give [1] as a general reference for most of these methods. Some popular software packages for the same task can be listed as ASAP, DASPK, and DASKADJOINT.In the theory of PSA, another important concept is the index structure which could be defined as the number of differentiations needed for transforming a DAE to an ODE. Intuitively, it is clear that all ODEs have index0. What defines the index is up to the constraints given in systems. For example, let us consider a simple predator-prey model
(1)y′=3x-y,x2+y=g(x),
where g is differentiable. Taking derivative of constraint equality, we get y′=3x-y, and 2xx′+y′=g′(x)⇒x′=(g′(x)-3x+y)/2x. Hence, new ODEs are given as
(2)y′=3x-y,x′=g′(x)-3x+y2x.
To obtain this, it takes one differentiation. Thus, the model has index 1.In order to capture the main idea of the PSA, let us consider the general form of the parameter-dependent DAEs given by(3)F(t,x,x˙,p)=0,x(0)=x0(p),
where x∈Rn and p∈Rp.It is not always the case but assume that we have index0 or 1 DAEs and convert (3) to explicit form of ODEs:
(4)x˙=f(t,x,p),x(0)=x0(p),
where (x,t,p)∈Rn×R×Rp.Sensitivity analysis requires the calculation of the term, namely,si, defined as the derivative of x with respect to pi; that is, si∶=∂x/∂pi. Since we are interested in partial derivatives, we can treat one parameter after another, while keeping the remaining ones fixed. Therefore, the derivative of (4) with respect to parameter pi is
(5)∂∂pi(∂x∂t)=∂∂pif(t,x,p)⟺∂∂t(∂x∂pi)=∂f∂pi+∂f∂x∂x∂pi.
Replacing si=∂x/∂pi into the right-hand side of (5), the ith sensitivity equation becomes
(6)s˙i=fpi+fxsi.
Since we have si=∂x/∂pi=[∂x1/∂pi,∂x2/∂pi,…,∂xn/∂pi]T=[si1,si2,…,sin]T, the previous process introduces n·p additional differential equations. Finally, defining a new variable
(7)z∶=[xs],
where
(8)s=[s1,s2,…,sp]=[s11,s12,…,s1n,s21,s22,…,s2n,…,sp1,sp2,…,spn],
we write
(9)z˙=[x˙s˙]=[f(t,x,p)∂fdpi+∂f∂pis]=[f(t,x,p)∂f∂p1+∂f∂p1s1∂f∂p2+∂f∂p2s2⋮∂f∂pp+∂f∂psp].
The initial condition takes the form
(10)∀i,si(0)=∂x0∂pi⟹z(0)=[x0∂x0∂pi]=[x0s11(0)s21(0)s31(0)⋮sn1(0)⋮s1p(0)s2p(0)⋮snp(0)].


## 3. Lorenz Equations

The Lorenz equations invented by E. N. Lorenz, a meteorologist and a pioneer of chaos theory, are typical examples of equations for system of differential algebraic equations that can be written as(11)[X′Y′Z′]=[σ(Y-X)rX-Y-XZ-bZ+XY],
where r>0 is Rayleigh number, σ is Prandtl number corresponding to temperature difference between two horizontal plates in convection problem, and b is a positive number. These equations arise in studies of convection and instability in planetary atmospheres, models of lasers and dynamos, and so forth. Although the stability and bifurcation properties of the Lorenz equations are studied in the literature [2, 3], to the best knowledge of the authors, the parameter sensitivity analysis of Lorenz equations has not been considered so far which is the main goal of this paper.The Lorenz equations are nonlinear due to the termsxy and xz. They are also symmetric equations, because the equations are invariant under (x,y)→(-x,-y). Thus, if (x(t),y(t),z(t)) is a solution of Lorenz equations, so is (-x(t),-y(t),z(t)). System of the Lorenz equations is dissipative. In other words, volumes in phase-space contract under the flow and σ and b are usually known as dissipation parameters. Next, we compute the fixed points of Lorenz equations. Letting each term of (11) be equal to 0, that is, letting x˙=y˙=z˙=0, we get the following identities:
(12)x=y,rx-y-xz=0,xy=bz,
which implies that (0,0,0) is a fixed point of the Lorenz equations. If x=y≠0, then from the appropriate ones of the previous equations we obtain z=x2/b,  rx-x-x(x2/b)=0 which gives us x2=b(r-1) by bearing in mind that x≠0. If r>1, then x*=y*=±b(r-1) and z*=r-1 that are the fixed-points of Lorenz equations. While r→1+, these equations generate a so-called pitchfork bifurcation. Now, we are making behavior analysis of these fixed points. Depending on the parameter r, we have three critical points; namely,
(13)x=y=z=0,a1:x=y=+[b(r-1)]1/2,a1:x=y=-[b(r-1)]1/2.In this study, behaviors dependent on initial conditions are not studied, and they are fixed asx0=0, y0=1, and z0=0. In the following figures different trajectories are given with respect to different r values.The behavior in Figure1 continues up to a value of r=24.08. After that, it becomes more complicated and chaotic; for example, for r=27 some periodic and aperiodic motions are observed as seen in Figure 2.Forr=3, projection of three-dimensional phase space onto z-y space (a) and corresponding solutions (b).
(a)
(b)Forr=27. The x in (a) and y in (b) are periodic, whereas z in (c) is aperiodic. (d) yz space for critical values a1 and a2.
(a)
(b)
(c)
(d)Further explanations of these and stability features of Lorenz equations might be seen, for instance, at [4]. In the next section, we study parameter sensitivity analysis of Lorenz equations.


## 4. Parameter Sensitivity Analysis of Lorenz Equations

Let us write the Lorenz equations having some initial conditions in the following way:(14)[x˙1x˙2x˙3]=[p1(x2-x1)p2x1-x2-x1x3-p3x3+x1x2]x1(0)=0,x2(0)=1,x3(0)=0.Our new variablez defined in (7) is given as follows:
(15)z=[x1x2x3s1s2s3],
where
(16){s1=dxdp1=[dx1dp1dx2dp1dx3dp1]T=[s11s12s13]Ts2=dxdp2=[dx1dp2dx2dp2dx3dp2]T=[s21s22s23]Ts1=dxdp3=[dx1dp3dx2dp3dx3dp3]T=[s31s32s33]T}.Finally, the sensitivity equations take the form(17)[z˙]=[x˙1x˙2x˙3s˙11s˙12s˙13s˙21s˙22s˙23s˙31s˙32s˙33]=[p1(x2-x1)p2x1-x2-x1x3-p3x3+x1x2x2-x1+p1(s12-s11)p2s11-s12-s11x3-x1s13-p3s13+s11x2+x1s12p1(s22-s21)x1+p2s21-s22-s21x3-x1s23-p3s23+s21x2+x1s22p1(s32-s31)p2s31-s32-s31x3-x1s33-x3-p3s33+s31x2+x1s32].Note that the initial conditions do not depend on parameters, so the new initial conditionssij=0  for all i,j=1,2,3. In other words, the initial conditions for the new variable are given as
(18)z0=[0,1,0,0,0,0,0,0,0,0,0,0]T.In the next section we are making some simulations in order to visualize the results.


## 5. Computational Results

In order to visualize the results, we made many different simulations as phase portraits, sensitivity analysis, and relation between the components. In this section, we present only some of the simulation results.After integrating the sensitivity equations, we getsi’s as a function of time so that one can analyze the change in the solution with respect to perturbations in the parameters. In these experiments, for time interval t=[0100], we solved the system by Matlab ODE solver, namely, ode45, based on an explicit Runge-Kutta Method. After integrating the sensitivity equations, we get si’s as a function of time so that one can analyze the sensitivity in phase portraits.Forr=30, the qualitative behaviors of the system might be seen in Figure 3.Figure 3
Forr=30, the behavior of the systems at x-z plane.Before entering the chaotic region which starts from the valuer=24.08, all nine sensitivity components demonstrate the same behavior. For example, taking r=20, we have well stable solutions, but the relation between the sensitivity variables is highly nonlinear, that is, a significantly important result for this well-known system. This is illustrated in Figure 4.(a)s13 and s31 components. (b) Phase portrait of s21 and s31.
(a)
(b)Remember thats21 and s31 represent sensitivity components of x1 with respect to p2 and p3. Figure 4(b) tells us that altering p2 and p3 in a simultaneous manner can affect the controllability of the system completely. In Figure 5, phase portrait for s31 and x1 is given.Figure 5
Phase portrait fors31 and x1.In the chaotic region, what happened to the sensitivity equations is that, first, they seemed to be very complicated and irregular when they are considered as a function of time as seen in Figure6.Solutions of sensitivity equations: (a)s22 and (b) s13.
(a)
(b)However, in the phase portraits ofs22 and s13, we obtained a completely linear relation between them as seen in Figure 7.(a)
(b)


## 6. Conclusion

When qualitative estimates of sensitivity are desired, a mathematical model of phenomena is desired or at least a relationship is required. However, such a model brings the questions concerning stability, optimality, sensitivity, and so forth. In this work, we concentrated only on the PSA of Lorenz equations. As we saw in the small application, it is possible to determine effects of parameters on model variables so that we can eliminate the less effective ones. Robustness can also be verified in some confidence intervals by just looking at the phase portraits. This enables us to control the system. This method is efficient if the number of variables is much more than that of the parameters. In a future work, we plan to study PSA for Van der Pool equations.



---

*Source: 537207-2013-09-16.xml*