User talk:Mimmo12

Metodo delta-sigma Modello matematico MULTIPARAMETRICO per lo studio delle anomalie di tipo fisico per la previsione dei terremoti del Prof. Domenico Caccamo Collaboratori: Prof. Fabio Barbieri, Prof. Caterina Laganà, Prof. Vincenza Pirrone.

Premessa “Delta-Sigma” è un algoritmo multiparametrico che permette di trovare le anomalie presenti prima di una forte replica con magnitudo M≥ 5.5, utilizzando parametri geofisici e statistici Fu il sismologo giapponese Fusakichi Omori che, prendendo in considerazione la sequenza sismica che si verificò dopo il terremoto di Nobi del 1891, riuscì a definire per primo la decrescenza del numero di repliche dopo la scossa principale, con la formula empirica nota come legge di Omori (detta anche legge iperbolica) e qui di seguito riportata:

n(t)= 1/tp             (1) con: n(t)= numero di repliche; p = costante caratteristica di quella sequenza (con valori genericamente compresi tra 1 ed 1.4); t = tempo misurato in giorni a partire dalla mainshock. Poiché la serie temporale del numero di scosse per giorno di una sequenza sismica, dal punto di vista matematico, può essere intesa come somma di un contributo deterministico, dato dal decadimento della frequenza di repliche secondo la  funzione intensità (t)=K(t+a)-b     (2) e di uno stocastico, dovuto alle fluttuazioni casuali attorno alla media fornita dalla (2), è, allora, possibile modellare il fenomeno del decadimento come un processo poissoniano non stazionario (Page,1968; Matsu’ura,1986), in cui la funzione intensità coincide con la formula di Omori modificata da Utsu nel 1961: (t)=n(t)=K(t+c)-p, (3) dove n(t) è il numero di repliche per unità di tempo a partire dalla scossa principale, calcolato nell’intervallo di tempo Δt; t è il tempo contato a partire dalla scossa principale; K è una costante che dipende dall’intensità della sequenza di repliche; c è un valore positivo che serve a compensare la scala dei tempi (approssimativamente da 1 ora ad 1 giorno); p è un valore, normalmente prossimo ad 1, che descrive il decadimento della sequenza nel tempo.

Metodo delta-sigma Il metodo Delta/Sigma ideato prof. Domenico Caccamo è stato progettato per studiare le anomalie di tipo fisico per la previsione delle forti repliche per una tipologia ben precisa di terremoti e precisamente: Terremoti crostali, vale a dire con profondità non superiore a 70Km; Terremoti che presentano una mainshock non inferiore a 7, con sequenze temporali decrescenti ( mainshock-aftershocks). La funzione matematica utilizzata nel modello è data da: n(t)=K×(t)-p+K1           (4)  (Cacamo D. 2002), in cui n(t) è il numero di eventi per unità di tempo a partire dalla scossa principale;

K è una costante che dipende dal numero totale degli eventi presenti nella sequenza; p definisce il decadimento della sequenza ed è normalmente prossimo ad 1 ; K1 è una costante che tiene conto della sismicità di fondo, che varia da area ad area e che serve a limitare i falsi allarmi (Caccamo et,al,2004a,b,2005,2007a,b,c,;D’Amico et al 2010).

Se consideriamo l’intervallo di tempo Δt, la media degli eventi per quell’ intervallo di tempo sarà   n(t)∙Δt  con deviazione standard  σ = √( n(t)∙Δt );  per un campionamento Δt = 1(giorno) il numero medio atteso di scosse giornaliere risulta essere  n(t), con una deviazione standard pari a σ =√(n(t)). Dagli studi fatti è emerso che le fluttuazioni stocastiche attorno alla media n(t) rientrano entro un range di ±2.5σ (≈ 99%), pertanto, dette fluttuazioni hanno una probabilità di accadimento inferiore all’1%. (Caccamo et. al., 2005, 2007; Bussetti, 1983). Data la generica serie completa di dati reali: Noss(tj)      (5) con j=1……d, d durata della sequenza osservata, inclusi gli ultimi 10 giorni in cui non si hanno più scosse, poiché le anomalie si presentano alcuni giorni prima dell’eventuale replica di magnitudo superiore a 5.5, al fine di evitare gli smorzamenti dovuti alle anomalie, che potrebbero essere incluse nelle operazioni di fitting, occorre inserire nell’ estrapolazione della serie calcolata uno shift s pari a 6, per cui si ha: ncalc(th)      (6) con th tempo di campionamento della serie calcolata h=a, a+1,a+2+……d, a=s+ q=10, numero minimo di punti utilizzati per l’estrapolazione. s=6 ,shift, q=2µ=4, µ=2 numero dei parametri k e p utilizzati nella (4); d = n1+n2    ultimo giorno della serie temporale, n1 numero di giorni pari al numero di repliche avvenute nelle prime 24 ore a partire dalla scossa principale n2 numero di giorni dopo n1,che contengono d. Il generico elemento {n_calc (t_g)}, si ottiene dal sottoinsieme {n_oss (t_e)}, con e =1,…, d-s, g numero d’intervalli di tempo in accordo col tempo di campionamento (che nel nostro caso è pari ad un giorno). Le differenze tra i valori calcolati e quelli osservati: Δ(th )= |〖n(th )〗_oss-〖n(th)〗_calc | si possono ritenere un’anomalia quando: Δ(t_h )≥2.5√(n_oss (t_j )  ) con √(n_oss (t_j )  )= σ Da quanto sopra detto, le anomalie trovate con la seguente formula: Δ/σ ≥ 2.5    (7) possono essere considerate come dei veri e propri precursori di forti repliche (Caccamo 2004°,b,2005,2007°,b,c;Bussetti1983;D’amico et al 2007). Concludendo:- Delta/sigma valuta dal punto di vista statistico il rapporto tra le differenze prese in valore assoluto dei valori dell’andamento temporale osservato e di quello calcolato, dette “delta”,con le deviazioni standard, dette “sigma”-.

Calcolo della magnitudo di completezza di una data sequenza. Per la previsione delle forti repliche occorre valutare anche la soglia di completezza dell’insieme dei dati ottenuti. A tal fine si usa il diagramma magnitudo- frequenza Gutenberg e Richter dato dalla legge: Y=〖log〗_10⁡N = a - b(M)      (8)

dove M è la magnitudo ed N è la frequenza calcolata nei primi 10 giorni, cioè la somma del numero di eventi con magnitudo uguale ad M calcolata in un intervallo di tempo Δt=10                                       (Caccamo et al,2002,2004°, b,2005,2007°, b,c). Nella (8)i parametri a e b sono delle costanti: a è legato al numero di eventi e si riferisce alla quantità di terremoti avvenuti in quella regione; b, in genere prossimo ad 1, è legato al coefficiente angolare ed al grado di fatturazione delle rocce. Gli studi fatti sulla previsione delle forti repliche hanno indotto a porre, nella Gutenberg e Richter (8), una perdita di dati non superiore al 10%.(Caccamo et al,2002,2004,2005,2007°,b,c,;D’Amico et al 2011), pertanto la (8) diviene: y=Y-10%Y          (9)  (Caccamo 2004).

Nel piano semi-logaritmico, (vedi fig. 1), le due rette di regressione, la (8) e la (9) vengono rappresentate con due colori differenti, la prima(8) di colore nero,la seconda(9) di colore rosso. Ciò premesso, la y calcolata con la(9) è quindi l’ordinata che deve essere considerata per determinare il minimo valore di magnitudo di soglia  Ms, che coinciderà con la magnitudo di completezza Mc.

Fig.1 Didascalia fig.1D. Caccamo. Diagramma magnitudo-frequenza di Gutenberg e Richter per i primi 10 giorni del terremoto del Cile con mainshock di M=8.8. La freccia indica il minimo valore di magnitudo trovato Ms=Mc=4.2.

L’esperienza sulle banche dati fino al 2013 ci aveva indotto ad accettare come dati attendibili, per la magnitudo di completezza Mc solo quelli con MS≥4.0; oggi, tenuto conto che la qualità dalle banche dati è notevolmente migliorata, possiamo affermare che non ha più senso porre questo limite per magnitudo MS di valore superiore a 3, per cui : Mc= MS (Caccamo 2017)

Per individuare la sequenza sia nello spazio che nel tempo utilizziamo la seguente procedura: calcoliamo la dimensione dell’area coinvolta mediante la (Utsu 1969)         log_10⁡L= 0.5∙Ms – 1.8; acquisiamo i dati contenuti in un riquadro di lato 3L, con centro nell’epicentro della mainshock; consideriamo un periodo di tempo pari ad un anno a partire dall’accadimento della mainshock; calcoliamo il “baricentro” della sequenza, considerando solo i dati relativi ai primi 10 giorni, perché risultano essere quelli rappresentativi dal punto di vista della completezza, in quanto la maggior parte delle repliche avviene proprio in questo periodo;

B_lat=(∑_(i=1)^n▒〖Lat〗_i )/n; 〖 B〗_lon=(∑_(i=1)^n▒〖Lon〗_i )/n;  (Caccamo D.2004)

dove Lati è la latitudine Loni è la longitudine dell’epicentro dalla i-esima replica ed n numero di repliche con M≥4.0 Il termine “baricentro” non deve essere inteso in senso fisico, piuttosto come punto geometrico dato dalla media aritmetica delle latitudini e longitudini di tutti i terremoti della sequenza avvenuti nei primi 10 giorni. Detta procedura, solo da noi utilizzata, serve a migliorare sensibilmente l’individuazione dell’area della sequenza delle repliche. La banca dati usata, per individuare la soglia di completezza, è quella dell’USGS (U.S. Geological Survey), in quanto questa presenta dei formati compatibili con quelli utilizzati per l’elaborazione di Delta/Sigma con matlab.

Uso del Matlab

Per la elaborazione dei dati si utilizza una serie di programmi fatti dal Prof. Fabio Barbieri con il matlab, che ci permette di risolvere problemi che coinvolgono l’uso dei minimi quadrati non lineari ed è, quindi, adatto per il modello Caccamo D.( 2002), che dal punto di vista matematico è una legge di potenza. Il programma utilizza metodi matematici di tipo iterativo, basati su algoritmi di ottimizzazione a larga scala, noti come metodo di Newton (Coleman e Li1994; Dennis1977). ( Caccamo 2004 a,b,2005,2007 a,b e c; D’Amico et al 2010) Ogni iterazione coinvolge la soluzione approssimata di un grande sistema lineare tramite il “metodo dei gradienti coniugati precondizionati” (LevenbergK1994; Marquardt D. 1963;More J 1977;Dennis J.E:JR1977;Coleman T.F., LI Y 1994,1996) (D.Caccamo 2002,2004,2005,2007a,b,c;Bussetti 1982, D’Amico 2007).

Per ogni sequenza analizzata, le figure in uscita dal programma “Delta-Sigma”sono, relativamente ai grafici : - Diagramma di Gutenberg-Richter dei primi 10 giorni; - Andamento temporale n(t); - Mappa della distribuzione epicentrale con localizzazione della mainshock   (triangolo) e del baricentro (cerchio); -  Mappa ipocentrale 3D; -  Mappa d’allineamento epicentrale con direzione di faglia; -  Andamento temporale giornaliero delle forti repliche (M > 5.5); - Numero cumulativo di eventi; - Energia cumulativa liberata durante l’intera sequenza; - Fit di n(t); - Fit di n(t) dei primi 10 giorni;

- Parametri frattali D0 e D2; relativamente ai file in excel : dati input; n(t) con Mc calcolata con la Gutenberg e Richter; dati del fit nei primi 10 giorni; eventi con M≥5.5; valori della Gutenberg e Richter. Delta/Sigma Parametri frattali D0 e D2

Variabili del programma Matlab I dati da inserire per i parametri variabili sono: numero di giorni entro i quali calcolare la magnitudo di completezza, fattore moltiplicativo per il calcolo del lato del box in funzione di L,	percentuale di errore accettabile sulla stima di Mc, valore di magnitudo minima di completezza stimato con la G-R, magnitudo di soglia, passo di campionamento in giorni, numero di intervalli di campionamento senza scosse di fine sequenza, passo di campionamento desiderato, shift in termini del passo di campionamento inserito al punto 8, costante di fondo.

Come applicativo riportiamo in tabella1  i dati relativi alle previsioni del terremoto del Cile con M=8.8 del 27/02/2010 con le relative percentuali di riuscita del sulle forti repliche ( tabella 1, Caccamo D.,pag.39, Atti del 32° Convegno Nazionale. Trieste 19-21 novembre 2013. Tema 1:Geodinamica), in tabella2 quelli  del Giappone 11 marzo 2011,    Mmain = 9.0 tabella 2, in tabella3 quelli del Perù 23 giugno 2001, Mmain = 8.4

Tabella 1 Cile 27 febbraio 2010, M = 8.8 Magnitudo delle repliche indagate	Numero repliche dopo 10° giorno inizio sequenza	Numero repliche precedute,entro 10 giorni, da Δ/σ > 2.5	Numero repliche non precedute,entro 10 giorni, da Δ/σ > 2.5	Percentuale di riuscita del metodo Δ/σ M≥5.5	39	33	6	84.6% M ≥ 5.6	33	27	6	81.0% M ≥ 5.7	28	24	4	85.4% M ≥ 5.8	26	225	1	96.3% M ≥5.9	20	20	0	100% M ≥ 6.0	11	11	0	100% M ≥ 6.1	9	9	0	100% M ≥ 6.2	6	6	0	100% M ≥ 6.3	7	7	0	100% M ≥ 6.4	5	5	0	100% M ≥ 6.5	5	5	0	100% M ≥ 6.6	5	5	0	100% M ≥ 6.7	4	4	0	100% M ≥ 6.8	2	2	0	100% M ≥ 6.9	2	2	0	100% M ≥7.0	1	1	0	100%

Tabella 2 Giappone 11 marzo 2011, Mmain = 9.0 Magnitudo delle repliche indagate	Numero repliche dopo il 10° giorno dall’inizio della sequenza	Numero repliche precedute, entro e non oltre i 10 giorni, da Δ/σ > 2.5	Numero repliche non precedute, entro e non oltre i 10 giorni, da Δ/σ > 2.5	Percentuale di riuscita del metodo Delta/Sigma

M > 5.5	67	44	23	65.7% M > 5.6	54	38	16	70.4% M > 5.7	42	31	11	73.8% M > 5.8	32	24	8	75% M > 5.9	26	20	6	76.9% M > 6.0	23	17	6	74.0% M > 6.1	16	11	5	68.8% M > 6.2	10	6	4	60% M > 6.3	7	6	1	85.7% M > 6.4	5	4	1	80% M > 6.5	5	4	1	80% M > 6.6	4	3	1	75% M > 6.7	2	2	0	100% M > 6.8	2	2	0	100% M > 6.9	2	2	0	100% M > 7.0	1	1	0	100%

Tabella 3

Perù 23 giugno 2001, Mmain = 8.4 Magnitudo delle repliche indagate	Numero repliche dopo il 10° giorno dall’inizio della sequenza	Numero repliche precedute, entro e non oltre i 10 giorni, da Δ/σ > 2.5	Numero repliche non precedute, entro e non oltre i 10 giorni, da Δ/σ > 2.5	Percentuale di riuscita del metodo Delta/Sigma (Δ/σ). M ≥ 5.3	16	14	2	87,5% M≥5.4	12	12	0	100% M≥ 5.5	6	6	0	100% M ≥ 5.6	6	6	0	100% M ≥ 5.7	4	4	0	100% M > 5.8	3	2	0	100% M > 5.9	3	2	0	100% M > 6.0	3	2	0	100% M > 6.1	3	2	0	100% M > 6.2	3	2	0	100% M > 6.3	3	2	0	100% M > 6.4	2	2	0	100% M > 6.5	2	2	0	100% M > 6.6	1	1	0	100% M > 6.7	1	1	0	100% M > 6.8	1	1	0	100% M > 6.9	1	1	0	100% M >7.0	1	1	0	100%

Concludendo Le forti repliche per magnitudo M> 6.7  presentano una percentuale di previsione pari al 100%, per i tre terremoti analizzati, escludendo i primi 10 giorni. Considerando anche i primi dieci giorni, giorni in cui si ha il maggior numero di repliche, in percentuale più del 50% di tutte le repliche considerate, le previsioni per scosse di M>5.5 si aggirano intorno all’90%. Il metodo è stato elaborato per effettuare previsioni in tempo reale, nel senso che le previsioni si possono ottenere inserendo le informazioni giorno per giorno a partire dal 11° giorno dopo la mainshock, con shifth s=4 e ∆t variabile. Per far questo bisogna essere collegati con una rete di stazioni di rilevamento dati che, fornendoli in tempo reale, ne possa  permettere l’immediata elaborazione.

Ulteriori approfondimenti sulla teoria dei modelli (Tesi di laurea del Prof. Fabio Barbieri)

Modellizzazioni.

Otsuka (1985, 1987) ha proposto una formula composta,, la quale esprime il prodotto tra la formula di Utsu 1961 e la funzione esponenziale legge di Omori. Per t grande, l’effetto della funzione esponenziale predomina. Il decadimento esponenziale dell’attività nei periodi finali di alcune sequenze di repliche è stato evidenziato da Utsu (1957), Mogi (1962), Watanabe e Kuroiso (1970) e Otsuka (1985). L’equazione  si rivela utile in alcuni casi limitati: nel modello considerato da Otsuka (1985), il valore di  è fortemente dipendente dal livello di magnitudo, cioè la diminuzione dell’attività è molto differente per le piccole e grandi repliche. Alcune funzioni mostrano una diminuzione proporzionale a  in un vasto dominio di t, ed una più rapida diminuzione per t molto grandi. Ad esempio, la funzione di Weibull, >0, >0, possiede tali proprietà se <<1, e può essere utilizzata per rappresentare l’attività delle repliche. La distribuzione gamma e la distribuzione log-normale hanno anche caratteristiche analoghe in un dominio limitato di valori dei parametri. Souriau et al. (1982) hanno utilizzato la distribuzione di Weibull per rappresentare l’andamento temporale di una sequenza di repliche avvenute nei Pirenei. Kisslinger (1993) ha utilizzato la funzione esponenziale per esprimere il numero di repliche non accadute sino ad un tempo t:,  ,  , dove   è il numero totale di repliche,  ,   è il numero cumulativo di repliche al tempo t. Un paragone dei valori AIC della formula di Omori modificata e dell’equazione  ,  ,  , indica che la prima si adatta meglio nella maggior parte dei casi se TS=0. Se TS è positivo, l’equazione,  ,  , diviene più rappresentativa del fenomeno nella metà dei casi, approssimativamente. Un modello “a trigger” è stato proposto da Vere-Jones e Davis (1966) nel corso di uno studio della sismicità della Nuova Zelanda. In questo modello la sismicità consiste di due tipi di eventi, cosiddetti primari e secondari. Gli eventi primari (scosse principali) si distribuiscono secondo un processo di Poisson con un tasso di accadimento costante . Ogni evento primario (che ha luogo all’istante ti) dà l’avvio ad una serie di eventi secondari (repliche) la cui intensità è espressa da:, per  , e  , per. Nella precedente relazione,  è una funzione normalizzata, in altri termini. Il numero totale A di eventi secondari innescati da un evento primario è una variabile aleatoria con media a e varianza v. Per applicare questo modello, la forma funzionale di  deve essere specificata, ma non è necessario identificare i singoli eventi come primari o secondari. Vere-Jones e Davis (1966) hanno provato ad utilizzare come  le funzioni di tipo esponenziale e la Omori modificata. Essi hanno mostrato che quest’ultima (p è stato posto pari a 1.25) è più appropriata per la sismicità della Nuova Zelanda, e ne hanno verificato la bontà comparando le curve spettrali teoriche e sperimentali. Il modello “Poisson generalizzato” proposto da Shlien e Toksoz (1970) costituisce un modello a trigger semplificato in cui. Tale modello ha riscosso successo nella rappresentazione di alcune proprietà della sismicità, ed è stato utilizzato da diversi ricercatori (ad es. Utsu, 1972; Sase, 1974; Hawkes e Adamopoulos, 1973). I primi due autori hanno utilizzato la funzione “Omori modificata”, mentre gli ultimi hanno introdotto la somma di due funzioni esponenziali le quali rappresentano le componenti a breve e lungo termine, rispettivamente. Nel modello a trigger, una replica non dà mai inizio ad una propria sequenza di repliche. I valori di magnitudo dei singoli eventi non sono presi in considerazione, cosicché l’entità della sequenza di repliche A non è correlata alla magnitudo della scossa principale. Hawkes e Adamopoulos (1973) hanno anche paragonato i modelli sopra menzionati ad un modello noto come “self-exciting point process model”, nel quale ogni evento può produrre eventi “figli”. Il tasso di accadimento all’istante t è espresso da, ove  e   rappresentano rispettivamente la sismicità a “rate” costante ed il tasso di attività “stimolato” da un evento originatosi al tempo ti. Gli autori hanno utilizzato la somma di due funzioni esponenziali per rappresentare la funzione. Il modello ETAS (Ogata, 1986, 1988, 1989, 1992, 1994) è un modello del tipo “self-exciting” (ad auto-eccitazione). Vengono considerati i terremoti di magnitudo Mz e più forti originatisi in una regione durante un periodo di tempo che va da un istante iniziale Ts all’istante finale Te. Il numero totale dei terremoti è indicato con N. L’attività sismica  include l’attività di fondo (poissoniana) con tasso di accadimento . Ogni terremoto è seguito da una sequenza di repliche, sebbene soltanto scosse di magnitudo Mz e più grandi siano incluse nell’insieme dei dati. L’attività delle repliche è rappresentata dalla formula di Omori modificata. Il tasso di accadimento al tempo t relativo alle repliche che seguono l’i-esimo terremoto verificatosi al tempo ti ed avente magnitudo Mi è dato dalla relazione, per  , ove K, , c e p sono costanti comuni a tutte le sequenze di repliche. Il numero cumulativo di terremoti all’istante t, successivo a Ts, è dato da    I cinque parametri , K, c,  e p rappresentano le caratteristiche dell’attività sismica della regione. Le stime di massima verosimiglianza dei cinque parametri possono essere ottenute massimizzando la funzione. La maggior parte dei valori di p e c ottenuti per vari insiemi di dati ricade negli intervalli 0.91.4, e 0.0030.3 giorni, rispettivamente (Ogata, 1986, 1988, 1992). Il parametro , variabile in genere tra 0.2 e 3.0, misura l’efficienza di una scossa nel generare una propria sequenza di repliche. L’attività di sciame presenta di solito valori di  alquanto bassi. La quantità  è grande nel caso di sequenze "scossa principale-repliche".