Analisi dati sperimentali e principi di modellazione tramite identificazione

Carlo Alberto Avizzano – carlo@sssup.itCorso di Simulink L2 - A.A. 2008/2009

Introduzione Sensori e sistemi si dicono ideali se producono

comportamenti e leggono grandezze in maniera aderente a un modello teorico senza alcun errore sui parametri e di misura.

Durante i corsi di controllo si immagina solitamente di poter disporre di sensori e sistemi ideali (semplificazione).

Tuttavia eccetto che nelle simulazioni, questi sistemi non sono realistici e quando possibile il comportamento dei sensori e dei sistemi va modellato tramite una campagna di acquisizione dati.

Obiettivi

Si necessita pertanto di poter caratterizzare il comportamento di un sensore quando I segnali sorgenti non sono piu’ modellabili come valori certi.

La misura porta I sistemi a interagire con una realta` in cui tutti I segnali sono soggetti ad un disturbo detto rumore che ne altera il valore teorico

0 1 2 3 4 5 6 7 8 9 10-20

-15

-10

-5

0

5

10

X=0:0.1:10Y=sin(X)+X-0.2*X.*X+ X.*sin(X)Z=rand(1,101)*2;plot (X,Y,X,Y+Z-1,'.',X,Y+1,X,Y-1)grid

Obiettivi In molti casi sperimentali il

processo di design del controllore dipendera` dalla corretta identificazione del sistema in studio e da una sua modellazione con strutture e parametri adeguati

Si tratta quindi di elaborare dati dalla conoscenza a priori e quella sperimentale che possiamo testare sul sistema da controllare.

In tutti i casi rumore, vincoli sul segnale o altre fonti esterne potranno disturbare l’acquisizione corretta di dati dal modello.

Modellare il comportamento lineare/non lineare di un algebricita` (interpolazione/regressione)

Modellare il comportamento lineare(/non lineare) di un sistema dinamico (identificazione)

Modellare il comportamento sequenziale di un sistema complesso (machine learning)

Il training set In tutti i casi di cui sopra si partira` sempre da

un set di addestramento sulla base del quale potremo manipolare le informazioni a nostra conoscenza per renderle adeguate a quelle del modello.

Curva caratteristica La caratteristica di

funzionamento (o relazione di taratura statica o funzione di taratura), è la relazione che lega fra loro i valori della grandezza in uscita dallo strumento a quelli in ingresso, stabilendo una corrispondenza tra i valori di lettura dello strumento e i risultati della misurazione.

Nel funzionamento, la catena di misura è quindi caratterizzata dalla sua funzione di trasferimento, ovvero dal legame che unisce il segnale d’ingresso x con il segnale d’uscita y:

y= f(x)

Yi + uiYi

I i

curva di taratura

I

Yi - ui

Y

Approssimazioni di Linearità Vi sono diversi tipi di approssimazioni: linearità riferita allo zero (zero based

linearity): la retta di riferimento passa per l’estremo inferiore della curva di taratura, corrispondente all’estremo inferiore del campo di misura, ed è tracciata in modo da rendere minimo il più elevato (in valore assoluto) degli scostamenti;

linearità riferita agli estremi (end point linearity): la retta di riferimento congiunge i due estremi della curva di taratura corrispondenti ai due estremi del campo di misura;

linearità secondo i minimi quadrati (least squares linearity): la retta di riferimento è quella che corrisponde al valor minimo della somma dei quadrati degli scostamenti;

linearità indipendente (independent linearity): la retta di riferimento è quella che rende minimo il più elevato (in valore assoluto) degli scostamenti.

Misura Supponendo di

suddividere l’intervallo di misura in N punti, il diagramma di taratura di solito riporta al posto della nuvola di punti acquisiti una linea continua detta ‘curva di taratura’.

Retta di taratura Nel caso in cui un

sensore abbia un comportamento quasi lineare, risulta possibile definire una curva di taratura lineare

Y=Au + B

Il metodo di ricerca si fa tramite una tecnica definita dei minimi quadrati

Metodo dei minimi quadrati

I valori dei parametri si scelgono quindi in maniera tale da minimizzare la somma degli scarti quadratici.

Larisoluzione di questo problema si ottiene imponendo le due derivate parziali uguali a zero

Soluzione specifica

Cosa succederebbe in un caso del 2 ordine

Vi sono tre parametri e tre derivate parziali, che si risolvono con l’inversionedi una matrice 3x3.

0 1 2 3 4 5 6 7 8 9 100

10

20

30

40

50

60

Soluzione Matlab

X=0:0.1:10; Y=3.1+5.2*X; Ym=3.1+5.2*X+

(rand(1,101)-0.5)*4.5; [AB,C]=polyfit(X,Ym,1) Yr= AB(2)+AB(1)*X; plot (X,Y,X,Ym, '.',X,Yr); AB = 5.1200 3.4997

Approccio statistico

Il comando Polyfit Polyfit esegue una

interpolazione ai minimi quadrati di un polinomio di grado qualunque (1== lineare, 2== quadratico,…)

Polyfit restituisce sia I parametri dei coefficienti polinomiali sia la matrice di covarianza dei campioni rispetto alla retta data

Misura di variabili La misura di una variabile e’ una lettura soggetta

ad differenti tipologie di errori e che per sua natura viene gestita quindi tramite metodologie statistiche.

In particolare si introduce il concetto di funzione di densita’ distribuzione della lettura P(x) che descrive la probabilita` che a una lettura di un sensore corrisponda una determinato valore del misurando.

L’integrale della densita’ di probabilita’ e’ detta funzione cumulata T(x) e’ corrisponde alla probabilita’ che data la lettura x (evento) il misurando sia minore o uguale all’evento Pr(X<=x)

Distribuzioni Frequenti Errore

costantemente distribuito indipendente

Errore a rampa Errore Gaussiano o

normale (Teorema del limite centrale)

Esempio Quantizzazione e Integrali

Operatori

Riassunto di regole

Relazioni tra minimi quadrati e probabilita’ Perche’ si usa il

metodo dei minimi quadrati?

Che senso ha nell’ottica della misura che abbiamo descritto

Proviamo a calcolare la probabilita’ di N misura congiunte supposte indipendenti e cerchiamo quale sia la soluzione ottima (piu’ probabile)

Identificazione di algebricita’ non lineari Modelli parzialmente strutturati (esempio

Robot) Linearita` nei parametri?

No: Punto di stima Si?

Metodi di identificazione

lsqcurvefit

Data-fitting nonlineare ttramite i minimi quadrati

Xdata: vettore di ingressi (xdata)

Ydata: un vettore di osservazioni (ydata)

Trovare i parametri “theta” che meglio “adattano” la funzione F(x,xdata) alle osservazioni.

Identificazione (matlab)

[P,residual] = lsqcurvefit(‘F(xdata,Px…)’,P0,xdata,ydata,[lb,ub,options,P1,P2,…])

residual: matrice dei residui

P0: valore iniziale di P

Lb, ub: limiti dell’algoritmo

options: opzioni di minimizzazione (vedi help optimset)

P: vettore dei parametri

Identificazione (matlab)

‘F’: funzione modello memorizzata in un M-file (F.m)

function ydata = F(P,xdata,P1,P2,…)

‘F’ puo’ anche richiamare un solver ‘ode’ per risolvere sistemi differenziali

Identificazione di modelli dinamici

Identificazione non parametrica Modelli strutturali per l’identificazione Identificazione tramite modelli regressivi

Prerequisiti di identificazione

Un modello derivato da dati non puo’ fornire piu’ informazioni di quante contenute nei dati stessi!

Per l’identificazione di un sistema dinamico e’ necessario che i dati presentati siano sufficientemente descrittivi del modello da rappresentare.

In particolare, essendo la funzione di uscita correlata al segnale in ingresso, i dati di ingresso devono essere tali da sollecitare tutti i modi propri del sistema.

Ingressi comuni sono: random binary, pseudo random, treni di gradini/impulsi, sinusoidi a frequenza variabile, ingressi dipendenti dalla struttura del modello

Operazioni di Make-Up dei dati

Operazione Comando

Rimozione dei disturbi a bassa frequenza

Rimozione dei disturbi ad alta frequenza

Rimozione degli outliers

Fill-up di dati mancanti

Sottrazione della media

Detrend del miglior fit lineare

Filtraggio passa alto Filtraggio passa

basso

Detrend a spezzate

Ispezione dei comportamenti frequenziali (fft) ed analisi del tempo di campionamento

IDData Manipulation

iddata - IDDATA(Y,U,Ts) iddata - IDDATA(Y,U,Ts, ‘FREQS’, Freqs)

Basic Operations fft/ifft - Transform time to frequency and

viceversa. Collections – getexp, subsref, merge Manipulation – Resample, nkshift, misdata Test/Analysis – plot , isnlarx, spa (Spectrum),

IDData Manipulation

detrend - Remove trends from data sets. delayest - Estimate the time delay (dead

time) from data. feedback - Investigate feedback effects in

data sets. idfilt - Filter data through Butterworth

filters. idinput - Generates input signals for

identification. isreal - Check if a data set contains real

data.

Classi di modelli

LTIBlack-box(i.e. ARX)

Modelli nonParametrici

Identificazione non parametrica Covf(Z,D) – Stima della funzione di covarianza tra Z=[y u]

con ritardo D IMPULSE (IDDATA) – Stima della risposta all’impulso step (IDDATA) - Stima della risposta al gradino

[ir,R,cl]=Cra(IDDATA,[D,na,plot]) – Stima della risposta all impulso: ir

Modelli LTI tempodiscreti

Alcune definizioni Model = IDPOLY(A,B,C,D,F,Noise,Ts) Ts=0 Continuous Time A,B,C,D,F Polynoms B leading zeros == delays == nk MIMO => B,F = Matrices

Also possible Model = IDPOLY(SYS) Polydata reverse the computation [A,B,C,D,F] = Polydata(Model) Model puo’ essere simulato sim(Model,Udata),

impulse(Model)

Modelli ARX

Il comando ARX Model = arx(data,order) Order =[na nb nk] Na = numero poli Nb = numero zeri Nk ritardo in ingresso Se B(z) == 0 nb = 0 ar command

Errore di predizione (ARX)

Classi Armax Model = ARMAX(IDDATA,[na nb nc nk]) [na nb nc nk] are the orders and delays of the

ARMAX model, according to the reference equation

A(q) y(t) = B(q) u(t-nk) + C(q) e(t)