ARMAX: One Step Prediction IIT Bombay Consider 2' nd order - - PowerPoint PPT Presentation

Automation Lab IIT Bombay

System Identification 99 4/10/2012 System Identification 99

ARMAX: One Step Prediction

) ( 1 1 ) ( 1 ) ( ) 2 ( ) 1 ( ) ( ) 3 ( ) 2 ( ) 2 ( ) 1 ( ) ( 1 d with model ARMAX

• rder

nd 2' Consider

2 2 1 1 2 2 1 1 2 2 1 1 3 2 2 1 2 1 2 1 2 1

k e q a q a q c q c k u q a q a q b q b k y k e c k e c k e k u b k u b k y a k y a k y

− − − − − − − −

+ + + + + + + + = − + − + + − + − + − − − − = = Difficulties: Sequences {y(k)} and {u(k)} are known but {e(k)} is unknown Non-Linear in parameter model – optimum can’t be computed analytically Solution Strategy Problem solved numerically using nonlinear optimization procedures

Automation Lab IIT Bombay

System Identification 100 4/10/2012 System Identification 100

Invertability of Noise Model

∞ < − = =

∑ ∑

∞ = ∞ = − 1

~ ) ( ~ ) ( ) ( ) (

i i i i

h i k v h k v q H k e i.e. stable is (q) H 1

• Crucial Property of Noise Model

Noise model and its inverse have to be stable i.e. all its poles and zeros should be inside the unit circle (follows from spectral factorization theorem) Key problem in identification: Find stable and inversely stable such H(q) and a white noise sequence {e(k)} 1 h i.e. polynomial monic' ' always is H(q) : Note

0 =

∞ < − = =

∑ ∑

∞ = ∞ =

i.e. stable is H(q) ) ( ) ( ) ( ) (

i i i i

h i k e h k e q H k v

Automation Lab IIT Bombay

System Identification 101 4/10/2012 System Identification 101

Example: A Moving Average Process ∑ ∑

∞ = ∞ = − −

− − = < − = + =

1 1

• )

( ) ( e(k) v(k)

• f

ts measuremen from recovered be can e(k) and 1 c if ) ( 1 1 ) ( H Then,

i i i i i

i k v c q c cq q

Inversion of Noise Model plays a crucial role in the procedure for model identification

c

• q

at zero and q at pole a has cq 1 H(q) i.e. sequence noise white a is {e(k)} where 1)

• ce(k

e(k) v(k) process MA

• rder

first a Consider

1

• =

= + = + = + = q c q

Automation Lab IIT Bombay

System Identification 102 4/10/2012 System Identification 102

Example: an ARMA Process

... ) 3 ( 0325 . ) 2 ( 067 . ) 1 ( 13 . ) ( .... 0325 . 067 . 13 . 1 5 . 1 8 . 1 ) (

3 2 1 1 1

+ − − − + − − = + − + − = + − =

− − − − −

k v k v k v k v q q q q q q e(k) v(k)

• f

ts measuremen from recovered be can e(k) and H Then,

1

• 0.5
• q

at zero and 0.8 q at pole a has 0.8q 1 0.5q 1 H(q) i.e. sequence noise white a is {e(k)} where 1)

• 0.5e(k

e(k) v(k) process ARMA

• rder

first a Consider

1

• 1
• =

= − + = − + = + + − = 8 . 5 . ) 1 ( 8 . q q k v neglected be can rate) decay

• n

(depending n k some after terms Thus, increases k as h stable inversely and stable is model noise When : Note

k

= → 0 ~

Automation Lab IIT Bombay

System Identification 103 4/10/2012 System Identification 103

One Step Prediction

[ ]

∑

∞ = −

− − = − = −

1 1

) ( ~ ) ( 1 ) ( ) 1 | ( ˆ

i i

i k v h k v q H k k v

1)

• (k

time upto ts measuremen

• n

based v(k) predict to want we and 1)

• (k

t upto v(t)

• bserved

have we Suppose ≤

∑ ∑

∞ = ∞ =

− + = − =

1

) ( ) ( ) ( ) (

i i i i

i k e h k e i k e h k v

∑

∞ =

− = − −

1

mean zero has e(k) as ) ( ) 1 | ( ˆ 1)

• (k

upto n informatio

• n

based v(k)

• f

n expectatio l Conditiona : ) 1 | ( ˆ

i i

i k e h k k v k k v

[ ]

) ( ) ( 1 ) ( ) ( 1 ) ( ) ( ) ( ) 1 | ( ˆ k v q H q H k e q H k e k v k k v − = − = − = −

Automation Lab IIT Bombay

System Identification 104 4/10/2012 System Identification 104

One Step Output Prediction

[ ]

[ ]y(k)

1 ) ( ) ( ) ( ) 1 | ( ˆ ) (

• r

y(k) ) ( 1 ) ( ) ( ) ( ) 1 | ( ˆ have we g Rearrangin

1 1

− + = − − + = −

− −

q H k u q G k k y q H q H k u q G q H k k y

v(k) G(q)u(k) y(k) have We and 1)

• (k

t upto u(t) and y(t)

• bserved

have we Suppose + = ≤

) 1 | ( ) ( ) ( ) 1 | ( − + = − k k v k u q G k k y 1)

• (k

time upto n informatio

• n

based y(k) predict to want we and

[ ]

) ( ) ( 1 ) ( ) (

1

k v q H k u q G

−

− + =

G(q)u(k)

• y(k)

v(k) However, =

[ ][

]

G(q)u(k)

• y(k)

) ( 1 ) ( ) ( ) 1 | ( ˆ

1 q

H k u q G k k y

−

− + = −

Automation Lab IIT Bombay

System Identification 105 4/10/2012 System Identification 105

ARX: One Step Predictor

) 1 | ( ˆ ) ( ) ( as estimated be can instant k'th at Residual known are RHS in terms All : Advantage − − = k k y k y k e

) ( 1 1 ) ( 1 ) ( 1 d with model ARX

• rder

nd 2' Consider

2 2 1 1 2 2 1 1 3 2 2 1

k e q a q a k u q a q a q b q b k y ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + + + ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + + + = =

− − − − − −

) ( 1 ) ( 1 ) 1 | ( ˆ is model this for predictor ahead step One

2 2 1 1 3 2 2 1

k y q a q a k u q b q b k k y ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ − − + ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + = −

− − − −

) 3 ( ) 2 ( ) 2 ( ) 1 ( ) 1 | ( ˆ equation difference to equivalent is which

2 1 2 1

− + − + − − − − = − k u b k u b k y a k y a k k y

Automation Lab IIT Bombay

System Identification 106 4/10/2012 System Identification 106

ARMAX: One Step Predictor

) 1 | ( ˆ ) ( ) ( as estimated be can instant k'th at Residual − − = k k y k y k ε

) ( 1 1 ) ( 1 ) ( 1 d with model ARMAX

• rder

nd 2' Consider

2 2 1 1 2 2 1 1 2 2 1 1 3 2 2 1

k e q a q a q c q c k u q a q a q b q b k y ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + + + + + ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + + + = =

− − − − − − − −

) ( 1 ) ( ) ( ) ( 1 ) 1 | ( ˆ is model this for predictor ahead step One

2 2 1 1 2 2 2 1 1 1 2 2 1 1 3 2 2 1

k y q c q c q a c q a c k u q c q c q b q b k k y ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + + − + − + ⎥ ⎦ ⎤ ⎢ ⎣ ⎡ + + + = −

− − − − − − − −

( ) ( )

) 2 ( ) 1 ( ) 3 ( ) 2 ( ) 3 | 2 ( ˆ ) 2 | 1 ( ˆ ) 1 | ( ˆ equation difference to equivalent is which

2 2 1 1 2 1 2 1

− − + − − + − + − + − − − − − − = − k y a c k y a c k u b k u b k k y c k k y c k k y

Automation Lab IIT Bombay

System Identification 107 4/10/2012 System Identification 107

ARMAX: One Step Predictor

) 2 ( ) 1 | 2 ( ˆ ; ) 1 ( ) | 1 ( ˆ ; ), ( ) 1 | ( ˆ : Assume parameters model for guess and 0,1,2,...} k : {u(k) : Given y y y y y y = = = − =

( ) ( )

) 2 | 3 ( ˆ ) 3 ( ) 3 ( ) 1 ( ) 2 ( ) ( ) 1 ( ) | 1 ( ˆ ) 1 | 2 ( ˆ ) 2 | 3 ( ˆ

2 2 1 1 2 1 2 1

y y y a c y a c u b u b y c y c y − = − + − + + + − − = ε

( ) ( )

) 3 | 4 ( ˆ ) 4 ( ) 4 ( ) 2 ( ) 3 ( ) 1 ( ) 2 ( ) 1 | 2 ( ˆ ) 2 | 3 ( ˆ ) 3 | 4 ( ˆ

2 2 1 1 2 1 2 1

y y y a c y a c u b u b y c y c y − = − + − + + + − − = ε

( ) ( )

) 4 | 5 ( ˆ ) 5 ( ) 5 ( ) 3 ( ) 4 ( ) 2 ( ) 3 ( ) 2 | 3 ( ˆ ) 3 | 4 ( ˆ ) 4 | 5 ( ˆ

2 2 1 1 2 1 2 1

y y y a c y a c u b u b y c y c y − = − + − + + + − − = ε

Automation Lab IIT Bombay

System Identification 108 4/10/2012 System Identification 108

ARMAX: One Step Predictor

{u(k)}. and {y(k)} sequences using (k)} { sequence generate can we ), , , , , , (a parameters model given and, (2) (1) (0) guesses initial with prediction start can We

2 1 2 1 2 1

ε ε ε ε c c b b a = = =

[ ] [ ]

) 3 | 2 ( ˆ ) 2 ( ) 2 ( ) 2 | 1 ( ˆ ) 1 ( ) 1 ( instances previous at residuals using ely, Alternativ − − − − = − − − − − = − k k y k y k k k y k y k ε ε

) 1 | ( ˆ ) ( ) ( ) 2 ( ) 1 ( ) 3 ( ) 2 ( ) 2 ( ) 1 ( ) 1 | ( ˆ as predictor step

• ne

rearrange can we

2 1 2 1 2 1

− − = − + − + − + − + − − − − = − k k y k y k k c k c k u b k u b k y a k y a k k y ε ε ε

Automation Lab IIT Bombay

System Identification 109 4/10/2012 System Identification 109

ARMAX: One Step Predictor

) 2 ( ) 1 ( ) ( : Assume parameters model for guess and 0,1,2,...} k : {u(k) : Given = = = = ε ε ε

) 2 | 3 ( ˆ ) 3 ( ) 3 ( ) 1 ( ) 2 ( ) ( ) 1 ( ) 1 ( ) 2 ( ) 2 | 3 ( ˆ

2 1 2 1 2 1

y y c c u b u b y a y a y − = + + + + − − = ε ε ε ) 3 | 4 ( ˆ ) 4 ( ) 4 ( ) 2 ( ) 3 ( ) 1 ( ) 2 ( ) 2 ( ) 3 ( ) 3 | 4 ( ˆ

2 1 2 1 2 1

y y c c u b u b y a y a y − = + + + + − − = ε ε ε ) 4 | 5 ( ˆ ) 5 ( ) 5 ( ) 3 ( ) 4 ( ) 2 ( ) 3 ( ) 3 ( ) 4 ( ) 4 | 5 ( ˆ

2 1 2 1 2 1

y y c c u b u b y a y a y − = + + + + − − = ε ε ε

Automation Lab IIT Bombay

System Identification 110 4/10/2012 System Identification 110

Parameter Estimation

∑

=

=

N k

N c c b b a Min c c b b a

3 2 2 1 2 1 2 1 2 1 2 1 2 1

) (k, 1 ) , , , , , (a ) ˆ , ˆ , ˆ , ˆ , ˆ , a ˆ ( function

• bjective

minimizes that Find

• n

Optimizati by n Estmimatio Parameter θ ε θ

) 1 | ( ˆ ) ( ) ( ) 2 ( ) 1 ( ) 3 ( ) 2 ( ) 2 ( ) 1 ( ) 1 | ( ˆ To Subject

2 1 2 1 2 1

− − = − + − + − + − + − − − − = − k k y k y k k c k c k u b k u b k y a k y a k k y ε ε ε

Nonlinear optimization problem Important to give good initial guess

Automation Lab IIT Bombay

System Identification 111 4/10/2012 System Identification 111

Prediction Error Method ∑

=

=

N k k N

N Z V

2

) (k, 1 ) , ( function

• bjective

minimizes that Find Method Error Prediction by n Estmimatio Parameter θ ε θ θ

( ) { }

, .... 2 , 1 : u(k) y(K), Z set data Given

N

N k = =

) ( ) , ( ) ( ) , ( ) ( : Model k e q H k u q G k y θ θ + =

[ ]

) ( ) , ( 1 ) ( ) , ( ) , ( ) 1 | ( ˆ predictor step

• ne

Optimal

1 1

k y q H k u q G q H k k y θ θ θ

− −

− + = − ) , 1 | ( ˆ ) ( ) (k, as defined is error prediction step One θ θ ε − − = k k y k y

Automation Lab IIT Bombay

System Identification 112 4/10/2012 System Identification 112

Two Tank System: ARMAX Model

[ ]

[ ]

) , , , , , ( to respect with minimized is ) 1 | ( ˆ ) ( e(k) function

• bjective

that such ) , , , , , , ( Estimate

2 1 2 1 2 1 3 2 3 2 2 1 2 1 2 1

c c b b a a k k y k y c c b b a a

N k N k

∑ ∑

= =

− − = = Ψ Optimization formulation

3

• 3

10 2.6813 ˆ : Variance Estimated 10 4.3601 E{e(k)} : Mean Estimated × = × =

2

λ

Residual {e(k)} Statistics Identified Model Parameters

q 0.2501 q 0.8367

• 1

C(q) q 0.01154 q 0.001748 B(q) q 0.68 q 1.651

• 1

A(q)

2

• 1
• 3
• 2
• 2
• 1

+ = + = + =

Automation Lab IIT Bombay

System Identification 113 4/10/2012 System Identification 113

Two Tank System: ARMAX Model

50 100 150 200 250

• 1

1 y(k) ARMAX(2,2,2,2): Measured and Simulated Ouputs 50 100 150 200 250

• 0.05

0.05 e(k) Time

2’nd Order ARMAX Model

Automation Lab IIT Bombay

System Identification 114

ARMAX Model: Correlation Functions

5 10 15 20 25

• 0.5

0.5 1 1.5 ARMAX Model AMX(2,2,2,2): Auto-correlation function of residuals e(k) lag

• 30
• 20
• 10

10 20 30

• 0.4
• 0.2

0.2 0.4 Cross corr. function between input u(k) and residuals e(k) lag

Automation Lab IIT Bombay

System Identification 115 4/10/2012 System Identification 115

Histogram of {e(k)}

• 0.06
• 0.04
• 0.02

0.02 0.04 0.06 10 20 30 40 50 60 70

Histogram of innovatoion sequence e(k) e(k) Number of Samples

Automation Lab IIT Bombay

System Identification 116 4/10/2012 System Identification 116

Comparison of Model Predictions

50 100 150 200 250

• 1
• 0.8
• 0.6
• 0.4
• 0.2

0.2 0.4 0.6

Time y(k) Measured and simulated model output

ARMAX(2,2,2,2) ARX(6.6.2) Plant OE(2,2,2)

Best Fit (%) ARMAX : 76.45 ARX : 76.37 OE : 77.38

Automation Lab IIT Bombay

System Identification 117 4/10/2012 System Identification 117

Comparison of Models

10 20 30 40 50 0.2 0.4 Sampling Instant Step Response y(k)

5 10 15 20 25 30 35 40 45 0.01 0.02 0.03

Sampling Instant Impulse Response y(k)

OE(2,2,2) ARX(6,6,2) ARMAX(2,2,2)

Automation Lab IIT Bombay

System Identification 118 4/10/2012 System Identification 118

Comparison of Models: Nyquist Plots

• 0.1

0.1 0.2 0.3 0.4

• 0.6
• 0.4
• 0.2

0.2 0.4 0.6 OE(2,2,2) ARX(6,6,2) ARMAX(2,2,2)

Automation Lab IIT Bombay

System Identification 119 4/10/2012 System Identification 119

Steps in Model Development

Selection of model structure Planning of experiments for estimation of

unknown model parameters

Design of input perturbation sequences Open loop / closed loop experimentation

Estimation of model parameters from

experimental data using optimization techniques

Model validation

Prediction capabilities Steady state behavior

Automation Lab IIT Bombay

System Identification 120 4/10/2012 System Identification 120

Model Structure Selection

Issues in Model Selection

• Process application (batch / continuous)
• Time scale of operation
• Type of application

(scheduling/optimization/MPC/Fault Diagnosis)

• Availability of physical knowledge / historical

data

• Development time and efforts

Model granularity decides how well we can make control / planning moves or diagnose / analyze process behavior