Corso di Identificazione dei Modelli e Analisi dei Dati (AA 2009...

29
Corso di Identificazione dei Modelli e Analisi dei Dati (AA 2011-2012) Prof. Sergio Bittanti I esercitazione di laboratorio 21/11/2011 Marcello Farina

Transcript of Corso di Identificazione dei Modelli e Analisi dei Dati (AA 2009...

Corso di Identificazione dei Modelli e Analisi dei Dati (AA 2011-2012)

Prof. Sergio Bittanti

I esercitazione di laboratorio 21/11/2011

Marcello Farina

Sommario • Serie temporali e sistemi

– Generazione di sistemi e serie temporali – Analisi di sistemi – Analisi di serie temporali (“identificazione non parametrica”)

• Predittori

• Identificazione

– Identificazione parametrica – Selezione della struttura – Validazione dei modelli identificati

• Esercizi

Serie temporali e sistemi

• Generazione di sistemi a tempo discreto – LTI format (rappresentazione esterna)

>> system = tf(num,den,Ts)

– num: numeratore

– den: denominatore

– Ts: tempo di campionamento (-1 default)

>> system=tf([1 1],[1 -0.8],1)

Transfer function:

z + 1

-------

z - 0.8

Serie temporali e sistemi

• Generazione di segnali

– Risposta allo scalino >> T=[0:Ts:Tfinal];

>> [y,T]=step(sistema,Tfinal); % step response

>> step(sistema,Tfinal) % plot

– Risposta a un ingresso qualsiasi >> T=[0:Ts:Tfinal];

>> u=sin(0.1*T); % ingresso sinusoidale

>> [y,T]=lsim(sistema,u,T);

Serie temporali e sistemi – Rumore bianco

>> N=100; % numero di campioni

>> lambda=1; % deviazione standard

>> mu=1; % media

• I metodo >> u=mu+lambda*randn(N,1);

• II metodo >> u=wgn(M,N,P); % segnale a media nulla

M: numero di segnali

P: potenza del segnale (dBW) (10log10(lambda^2))

• Somma di segnali e WGN >> u=awgn(X,SNR,SIGPOWER);

X: segnale di partenza

SRN: signal-to-noise ratio (dB)

SIGPOWER: potenza del segnale X (dBW) (20log10(.))

- se SIGPOWER=‘measured’, AWGN misura la potenza el

segnale X.

Serie temporali e sistemi

– Risposta di un sistema (lineare) a segnali qualsiasi >> u=1+randn(3000,1);

>> [y,T]=lsim(sistema,u,[1:3000]);

– Il toolbox di identificazione si serve di istruzioni differenti per definire un modello (detto IDPOLY format) >> modello = idpoly(A,B,C,D,F,NoiseVariance,Ts)

A,B,C,D,F sono vettori (contengono coefficienti di

polinomi in ordine crescente di esponente di z-1)

Modello:

A(q) y(t) = [B(q)/F(q)] u(t-nk) + [C(q)/D(q)] e(t)

(Se ARMAX: F=[], D=[])

>> modello = idpoly(sistema)

Sistema è un sistema LTI (definito tramite ss, tf, etc.)

Serie temporali e sistemi

– Risposta di un sistema (lineare, in IDPOLY format) a segnali qualsiasi

>> y=sim(model,u);

– I dati possono essere forniti come IDDATA objects: >> y_id=sim(y,u,Ts);

u: input data

y: output data

Ts: sampling time

Serie temporali e sistemi

• Analisi del sistema (LTI e IDPOLY format):

– Zeri e poli: >> [Z,P,K]=zpkdata(sistema)

Z = [-1] % posizione degli zeri

P = [0.8000] % posizione dei poli

K = 1 % guadagno (non è pari a G(1)!)

– Diagramma di Bode: >> [A,phi] = bode(sistema,W)

W: valori della frequenza ‘w’ in cui calcolare ‘A(w)’ e

‘phi(w)’

A: amplificazione (A(w))

Phi: sfasamento (phi(w))

Serie temporali e sistemi

– Diagramma di Bode:

>> [A,phi] = bode(sistema,[0.1 1])

A(:,:,1) =9.1179

A(:,:,2) =1.9931

phi(:,:,1) =-24.2456

phi(:,:,2) =-78.5036

>> bode(sistema)

>> % disegna il diagramma

(modulo –dB- e fase)

Serie temporali e sistemi

• Analisi di serie temporali - m=mean(X) % Average or mean value.

– Y=detrend(X,’constant’) % Remove constant trend from data sets.

– Y=detrend(X,’linear’) % Remove linear trend from data sets.

– Funzione di covarianza (¿) >> gamma=covf(y,10) % tau=0,1,…,10-1

>> plot([0:9],gamma)

– Spettro ¡(!) – metodo del periodogramma >> periodogram(y)

>> [Pxx,w]=periodogram(y);

>> plot(w/pi,10*log10(Pxx))

° Attenzione: covf calcola la funzione di correlazione! BISOGNA PRIMA DEPOLARIZZARE IL SEGNALE

Serie temporali e sistemi

• IMAD_LAB simulink toolbox. Digitare:

>> imad_lab

Predittori

• Dato il sistema

il predittore a k passi si calcola come

dove Ek(z) e Rk(z) sono quoziente e resto della lunga divisione, e

(si calcolano “a mano”)

y(t) =B(z)

A(z)u(t¡ ¿) +

C(z)

A(z)e(t), e(t) =WN(0; ¸2)

y(t=t¡ k) =B(z)Ek(z)

C(z)u(t¡ ¿) +

Rk(z)

C(z)y(t)

=B(z)Ek(z)

C(z)u(t¡ ¿) +

~Rk(z)

C(z)y(t¡ k)

R(z) = ~R(z)z¡k

Predittori

C(z)

A(z)

B(z)

A(z)z¡¿

e(t)

u(t) y(t)z¡k

y(t¡ k)

B(z)Ek(z)

C(z)

~Rk(z)

C(z)

S Py(t=t¡ k)

Predittori

- Predittore

>> Yp = predict(model,data,k) % predittore a k passi.

model: IDPOLY format

data: IDDATA format

k: passi di predizione

– Errore di predizione

>> Epsilon = pe(model,data) % errore di predizione a k passi

– TEST DI ARDERSON PER L’ERRORE DI PREDIZIONE (regione di

confidenza del 99%)

>> resid(model,data)

Predittori

• Esempio:

>> modello=idpoly([1 0.5],[0 1 0.3],[1 0.6],[],[],1,1)

Discrete-time IDPOLY model: A(q)y(t) = B(q)u(t) + C(q)e(t)

A(q) = 1 + 0.5 q^-1

B(q) = q^-1 + 0.3 q^-2

C(q) = 1 + 0.6 q^-1

This model was not estimated from data.

Sampling interval: 1

>> dati=randn(100,2);

>> resid(modello,dati)

Identificazione >> ident % start the MATLAB identification toolbox GUI

Identificazione

• IMAD_LAB simulink

toolbox

Identificazione Altre istruzioni Matlab (2 Identification toolbox)

• Parametric model estimation – ar - AR-models of signals using various approaches.

– armax - Prediction error estimate of an ARMAX model.

– arx - LS-estimate of ARX-models.

– bj - Prediction error estimate of a Box-Jenkins model.

– init - Initialize (randomize) the parameters of a model.

– ivar - IV-estimates for the AR-part of a scalar time series.

– iv4 - Approximately optimal IV-estimates for ARX-models.

– n4sid - State-space model estimation using a sub-space method.

– nlarx - Prediction error estimate of a nonlinear ARX model.

– nlhw - Prediction error estimate of a Hammerstein-Wiener model.

– oe - Prediction error estimate of an output-error model.

– pem - Prediction error estimate of a general model.

Identificazione • Model structure selection

– aic - Compute Akaike's information criterion.

– fpe - Compute final prediction criterion.

– arxstruc - Loss functions for families of ARX-models.

– selstruc - Select model structures according to

various criteria.

– idss/setstruc - Set the structure matrices for idss

objects.

– struc - Generate typical structure matrices for

ARXSTRUC.

• Per maggiori dettagli si veda:

>> help ident

Esercizi

Esercizi

Esercizi

Esercizi

Esercizi

Esercizi

Esercizi

Esercizi

Esercizi

Esercizi