Legge lognormale a tre parametri

Legge lognormale a tre parametri Paolo Maccallini Settembre 2023 Generalizzo la definizione di legge lognormale, introducendo un terzo parametro, comunemente detto threshold e indicato b in questo documento, che permette di traslare lungo l’asse delle ascisse il sostegno della legge. Il terzo parametro consente spesso una maggiore aderenza ai dati sperimentali, come nel caso della potenza cardiaca normalizzata (cardiac power index). Risolvo il problema della determinazione dei tre parametri con il metodo dei momenti. 1 Legge lognormale DEF. 1 (legge di lognormale). Una v.a. aleatoria 𝑋 segue la legge lognormale di parametri 𝜇, 𝜎 quando la sua densità si scrive: 1 𝑓𝑋 (𝑡) = {√2𝜋𝜎𝑡 0, Eq. 1 𝑒 1 ln 𝑡−𝜇 2 ) 𝜎 −2 ( , 𝑝𝑒𝑟 𝑡 > 0 𝑝𝑒𝑟 𝑡 ≤ 0 dove 𝜇 è un parametro reale, e 𝜎 è un parametro reale non negativo. Si sconsideri ora che 1 √2𝜋𝜎 ∫ 𝑡 −∞ 𝑒 1 ln 𝑥−𝜇 2 ) 𝜎 −2( 𝑥 𝑑𝑥 = 1 √2𝜋𝜎 𝑡 ∫ 𝑒 0 1 ln 𝑥−𝜇 2 ) 𝜎 −2 ( 𝑑(ln 𝑥 ) = Dunque 𝑓𝑋 è una legge e la sua funzione di ripartizione è Eq. 2 𝐹𝑋 (𝑡) = { ln 𝑡−𝜇 Φ𝜇,𝜎 (ln 𝑡) = Φ ( 0, 𝜎 1 √2𝜋𝜎 ∫ ln 𝑡 −∞ 𝑒 1 𝑧−𝜇 2 ) 𝜎 −2( 𝑑𝑧 = Φ𝜇,𝜎 (ln 𝑡) ) , 𝑝𝑒𝑟 𝑡 > 0 𝑝𝑒𝑟 𝑡 ≤ 0 PROP. 1 (legge lognormale e legge normale). Se 𝑋~𝑁(𝜇, 𝜎 2 ), allora la v.a. 𝑌 = 𝑒 𝑋 è una legge lognormale di parametri 𝜇, 𝜎. DIM. Per 𝑦 > 0 si ha 𝐹𝑌 (𝑦) = 𝑃(𝑌 < 𝑦) = 𝑃(ln 𝑌 < ln 𝑦) = 𝑃(𝑋 < ln 𝑦) = Φ𝜇,𝜎 (ln 𝑦) che è la funzione di ripartizione della legge lognormale ∎ ESEMPIO 1 (mediana). La mediana della legge lognormale è data da 𝐹𝑋 (𝑚 ) = 1 ln 𝑚 − 𝜇 1 ln 𝑚 − 𝜇 ⇔ Φ( )= ⇔ = 0 ⇔ 𝑚 = 𝑒𝜇 2 𝜎 2 𝜎 La mediana di una legge lognormale dunque dipende solo da 𝜇. Questo lo si desume anche dalla Figura 1 dove sono riportate leggi e funzioni di ripartizione per tre valori di 𝜎, con 𝜇 costante. Si vede che le tre funzioni di ripartizione si incontrano tutte nel punto di coordinate (𝑒, 0.5), proprio perché condividono la stessa mediana. CODICE 1 (legge lognormale). Il seguente codice traccia il grafico della densità e della funzione di ripartizione della legge lognormale di parametri 𝜇, 𝜎. La funzione di ripartizione si ottiene per 2 integrazione numerica, con metodo di cavalieri. La densità attinge il suo massimo in 𝑡𝑀 = 𝑒 𝜇−𝜎 , infatti la derivate di 𝑓𝑋 si scrive 𝑑𝑓𝑋 (𝑡) 1 𝑒 = [− 𝑑𝑡 √2𝜋𝜎 1 ln 𝑡−𝜇 2 ) 𝜎 −2( 𝑡2 ⇔ − 𝑒 1 ln 𝑡−𝜇 2 ) 𝜎 −2( 2𝑡 1 ln 𝑡−𝜇 2 − ( ) 𝑒 2 𝜎 ln 𝑡 − 𝜇 2 ln 𝑡 − 𝜇 ( )] = − (1 + )=0⇔ 𝜎 𝜎2 𝑡𝜎 √2𝜋𝜎𝑡 2 ln 𝑡 − 𝜇 2 = −1 ⇔ ln 𝑡 = 𝜇 − 𝜎 2 ⇔ 𝑡 = 𝑒 𝜇−𝜎 𝜎2 L’andamento della legge lognormale è riportato in Figura 1 per 𝜇 = 1 e per tre diversi valori di 𝜎. Figura 1. Densità (a sinistra) e funzione di ripartizione per la legge lognormale con 𝜇 = 1, in corrispondenza di tre valori di sigma. % file name = legge_lognormale % date of creation = 11/05/2019 clear all delta = 0.01; index = 601 sigma (1)= sqrt(2.); sigma (2)= 1.; sigma (3)= 1./sqrt(2.); media (1)= 1.; media (2)= 1.; media (3)= 1.; for j=1:3 moltiplicatore(j) = 1/(sigma(j)*sqrt(2.*pi)); n(j,1) = delta; for i=2:1:index; n(j,i) = delta + n(j,i-1); end for i=1:1:index f(j,i) = ( e^( (-0.5)*( ( log(n(j,i)) -media(j) )/sigma(j) )^2. ) )/n(j,i); end integrale(1) = 0.; integrale(3) = integrale(1) + delta*( f(j,1) + ( 4*f(j,2) ) + f(j,3) )/3; integrale(2) = integrale(3)*0.5; 1 for k=2:1:index-2 integrale(k+2) = integrale(k) + delta*( f(j,k) + ( 4*f(j,k+1) ) + f(j,k+2) )/3; end for i=1:1:index F(j,i) = moltiplicatore(j)*integrale(i); end end figure(1) plot(n(1,:),moltiplicatore(1)*f(1,:),'-k','Linewidth', 1) hold on plot(n(2,:),moltiplicatore(2)*f(2,:),'-.k','Linewidth', 3) hold on plot(n(3,:),moltiplicatore(3)*f(3,:),'.k','Linewidth', 1) legend('{\mu}=1 {\sigma}=1.4142','{\mu}=1 {\sigma}=1','{\mu}=1 {\sigma}=0.7071', "location", 'northwest') xlabel('x') ylabel('f(x)') axis ([0. index*delta 0. 0.4]) grid on figure(2) plot(n(1,:),F(1,:),'-k','Linewidth', 1) hold on plot(n(2,:),F(2,:),'-.k','Linewidth', 3) hold on plot(n(3,:),F(3,:),'.k','Linewidth', 1) legend('{\mu}=1 {\sigma}=1.4142','{\mu}=1 {\sigma}=1','{\mu}=1 {\sigma}=0.7071', "location", 'northwest') xlabel('x') ylabel('F(x)') grid on axis ([0. index*delta 0. 1.]) PROP. 2 (legge lognormale). La speranza matematica e la varianza della legge lognormale di parametri 𝜇, 𝜎 valgono rispettivamente 𝜎2 𝐸 [𝑋] = 𝑒 𝜇+ 2 , 𝑉𝑎𝑟 [𝑋] = 9𝜎2 2 2 𝐸 [𝑋 2 ] = 𝑒 2𝜇+2𝜎 , 2 2 𝑒 2𝜇+𝜎 (𝑒 𝜎 DIM. Per il momento di ordine 𝑛 si ha − 1), 𝐸 [𝑋 𝑛 ] = Si consideri ora la sostituzione 1 √2𝜋𝜎 𝐸 [𝑋 3 ] = 𝑒 3𝜇+ 𝛾 [𝑋 ] = 𝑒 ∫ +∞ 0 𝑡 𝜎2 𝜇+ 2 , 2 𝑒 2𝜇+𝜎 (𝑒 3𝜎 − 1) − 3 𝜎3 1 ln 𝑡−𝜇 2 ) 𝜎 𝑛−1 −2 ( 𝑒 2 𝑛2 𝜎2 2 𝐸 [𝑋 𝑛 ] = 𝑒 𝑛𝜇+ 𝑑𝑡 ln 𝑡 − 𝜇 = 𝑥 ⇒ ln 𝑡 = 𝜎𝑥 + 𝜇 ⇒ 𝑡 = 𝑒 𝜎𝑥+𝜇 ⇒ 𝑑𝑡 = 𝜎𝑒 𝜎𝑥+𝜇 𝑑𝑥 𝜎 Quindi – cambiano opportunamente l’estremo di integrazione inferiore – abbiamo 𝐸 [𝑋 𝑛 ] = 1 √2𝜋 ∫ +∞ −∞ 1 2 (𝑒 𝜎𝑥+𝜇 )𝑛−1 𝑒 −2𝑥 𝑒 𝜎𝑥+𝜇 𝑑𝑥 = 2 1 √2𝜋 ∫ +∞ −∞ 𝑒 1 −(2𝑥2 −𝑛𝜎𝑥−𝑛𝜇) 𝑑𝑥 1 Consideriamo la ulteriore sostituzione 𝑥 2 = 𝑧 2 ⟺ 𝑥 = 𝑧√2 e abbiamo 2 𝐸 [𝑋 𝑛] +∞ Consideriamo ora l’integrale ∫−∞ 𝑒 −(𝑥 = 𝑒 𝑛𝜇 √𝜋 2 +𝜂𝑥) 𝐸 [𝑋 𝑛 ] = ∫ +∞ −∞ 2 −𝑛 𝑒 −(𝑧 𝜂2 √2𝜎𝑧) 𝑑𝑧 𝑑𝑥 = √𝜋𝑒 4 e abbiamo 𝑒 𝑛𝜇 √𝜋 √𝜋𝑒 𝑛2 2𝜎2 4 𝑛2 𝜎2 2 = 𝑒 𝑛𝜇+ Quindi in particolare abbiamo i primi tre momenti indicati nella tesi. Per la varianza abbiamo poi: 2 2 2 2 2 𝑉𝑎𝑟 [𝑋] = 𝐸 [𝑋 2 ] − 𝐸 2 [𝑋] = 𝑒 2𝜇+2𝜎 − 𝑒 2𝜇+𝜎 = 𝑒 2𝜇+𝜎 (𝑒 𝜎 − 𝑒 𝜎 ) Per l’indice di asimmetria si ha 𝜎2 𝜎2 𝜎2 2 𝜎2 2𝜇+4𝜎 − 3 − 𝑒 2𝜇+𝜎 𝐸 [𝑋 3 ] − 3𝐸 [𝑋] − 𝐸 3 [𝑋] 𝑒 3𝜇+9 2 − 3𝑒 𝜇+ 2 − 𝑒 3𝜇+3 2 𝜇+ 2 𝑒 𝛾 [𝑋 ] = = = 𝑒 𝜎3 𝜎3 𝜎3 2 Con un’altra semplice manipolazione si ottiene la tesi ∎ OSSERVAZIONE 1. Può essere utile ricavare i parametri della legge lognormale in funzione di media e varianza: 2 2 Quindi 𝑒 𝜎 = 𝑉𝑎𝑟[𝑋] 𝐸2 [𝑋] 𝑉𝑎𝑟[𝑋] + 1 ⇒ 𝜎 2 = ln ( 𝐸2[𝑋] 2𝜇 = ln 𝐸 2 [𝑋] − ln ( Si conclude dunque che 𝜎2 2 2 2 2𝜇 + 𝜎 2 = ln 𝐸 2 [𝑋] ⇒ 𝑉𝑎𝑟 [𝑋] = 𝑒 2𝜇+𝜎 (𝑒 𝜎 − 1) = 𝐸 2 [𝑋](𝑒 𝜎 − 1) + 1) da cui si ricava 𝑉𝑎𝑟 [𝑋] 𝐸 4 [𝑋 ] + 1) = ln 𝐸 2 [𝑋 ] 𝑉𝑎𝑟 [𝑋] + 𝐸 2 [𝑋] 𝑉𝑎𝑟 [𝑋] + 1) , = ln ( 2 𝐸 [𝑋 ] 𝜇 = ln Legge lognormale a tre parametri 𝐸 2 [𝑋 ] √𝑉𝑎𝑟 [𝑋] + 𝐸 2 [𝑋] DEF. 2 (legge di lognormale a tre parametri). Sia 𝑋 una legge lognormale di parametri 𝜇, 𝜎 (DEF. 1). Si consideri la legge di 𝑌 = 𝑋 + 𝑏, detta legge lognormale a tre parametri. La sua legge si scrive 𝑓𝑌 (𝑦) = 1 √2𝜋𝜎(𝑦 − 𝑏 ) 𝑒 2 1 ln(𝑦−𝑏) −𝜇 ) 𝜎 2 − ( per 𝑦 > 𝑏, zero altrimenti. Il valore di 𝑏 può essere positivo, nullo, negativo e rappresenta la traslazione sull’asse delle ascisse del sostegno di 𝑓𝑋 [1]. Valgono le relazioni seguenti: 3 𝑌 = 𝑋 + 𝑏 = 𝑒 𝑍 + 𝑏, 𝑍~𝑁(𝜇, 𝜎 2 ), 𝑋~𝐿𝑁(𝜇, 𝜎 2 ), PROP. 3 (media e varianza). Si dimostra che 𝜎2 𝑌~𝐿𝑁(𝜇, 𝜎 2 , 𝑏 ) 2 𝐸 [𝑌 ] = 𝐸 [𝑋] + 𝑏 = 𝑒 𝜇+ 2 + 𝑏, 2 𝑉𝑎𝑟 [𝑌 ] = 𝑉𝑎𝑟 [𝑋] = 𝑒 2𝜇+𝜎 (𝑒 𝜎 − 1) DIM. Segue immediatamente da PROP. 2 ∎ PROP. 4 (momenti). Si dimostra che 𝐸 [𝑌 𝑛] = ∑ DIM. Abbiamo 𝐸 [𝑌 𝑛 ] = 𝐸 [(𝑋 + 𝑏 )𝑛 ] = 𝐸 [∑ 𝑛2 𝜎2 2 E dunque, essendo 𝐸 [𝑋 𝑘 ] = 𝑒 𝑛𝜇+ 𝑛 𝑘2 𝜎2 𝑛 ( ) 𝑏 𝑛−𝑘 𝑒 𝑘𝜇+ 2 𝑘=0 𝑘 𝑛 𝑛 𝑛 ( ) 𝑏 𝑛−𝑘 𝐸 [𝑋 𝑘 ] ( ) 𝑋 𝑘 𝑏 𝑛−𝑘 ] = ∑ 𝑘 𝑘=0 𝑘 𝑘=0 𝑛 abbiamo la tesi ∎ PROP. 5 (indice di asimmetria). Si dimostra che 𝜎2 𝜎2 𝜎2 2 𝑒 3𝜇+9 2 + 3𝑏𝑒 2𝜇+2𝜎 − 3𝑒 𝜇+ 2 − 3𝑏 − 𝑒 3𝜇+3 2 − 3𝑏𝑒 2𝜇+𝜎 𝛾 [𝑌 ] = 𝜎3 2 DIM. si può considerare Eq. 6 e scrivere Essendo poi 𝛾 [𝑌 ] = 𝐸 [𝑋 3 ] + 3𝑏𝐸 [𝑋 2 ] − 3𝐸 [𝑋] − 3𝑏 − 𝐸 3 [𝑋] − 3𝑏𝐸 2 [𝑋] 𝜎3 2 2 2 2 2 𝜎 2 𝜎 3 𝜎 3 3 3 3 𝐸 [𝑌 3] = ( ) 𝑏 3 + ( ) 𝑏 2 𝑒 𝜇+ 2 + ( ) 𝑏𝑒 2𝜇+ 2 + ( ) 𝑒 3𝜇+ 2 = 1 0 2 3 segue la tesi ∎ = 𝑏3 + 𝜎2 3𝑏 2 𝑒 𝜇+ 2 + 3𝑏𝑒 22 𝜎2 2 2𝜇+ +𝑒 32 𝜎2 2 3𝜇+ PROP. 6 (mediana). Si dimostra che la mediana della legge 𝐿𝑁(𝜇, 𝜎 2 , 𝑏 ) è 𝑚𝑦 = 𝑏 + 𝑒 𝜇 . DIM. Partendo dalla definizione di mediana abbiamo 2 2 𝑚𝑦 𝑚𝑦 1 ln(𝑦−𝑏) −𝜇 1 ln(𝑦−𝑏)−𝜇 1 1 1 ) ( ) −2( − 𝜎 𝜎 𝑑𝑦 = ∫ 𝑑𝑦 =∫ 𝑒 𝑒 2 2 √2𝜋𝜎(𝑦 − 𝑏 ) −∞ √2𝜋𝜎(𝑦 − 𝑏 ) 𝑏 4 Operando la sostituzione 𝑥 = 𝑦 − 𝑏 abbiamo 𝑚𝑦 −𝑏 1 ln 𝑥−𝜇 2 ln 𝑥 − 𝜇 𝑚𝑦−𝑏 1 ln(𝑚𝑦 − 𝑏) − 𝜇 1 −2( 𝜎 ) 𝑑𝑥 = Φ ( =∫ 𝑒 )| ) = Φ( 𝜎 𝜎 2 √2𝜋𝜎𝑥 0 0 Si osservi che la funzione Φ(𝑥 ) attinge il valore nullo per 𝑥 = 0 e dunque 𝑚𝑦 = 𝑏 + 𝑒 𝜇 , che generalizza quanto trovato per la legge lognormale ∎ 2 2 2 OSSERVAZIONE 2 (media e mediana). Poiché 𝑒 𝜎 > 1, segue 𝑒 𝜇 𝑒 𝜎 > 𝑒 𝜇 , quindi 𝑏 + 𝑒 𝜇 𝑒 𝜎 > 𝑏 + 𝑒 𝜇 . Ovvero si ha sempre 𝐸 [𝑌 ] > 𝑚𝑌 . PROP. 7 (metodo dei momenti). Ci si chiede quali siano i parametri di 𝑌~𝐿𝑁(𝜇, 𝜎 2 , 𝑏 ) in funzione di 𝐸 [𝑌 ], 𝑉𝑎𝑟 [𝑌 ], 𝑚𝑦. Ciò equivale a chiedersi quale sia la legge lognormale a tre parametri che meglio si adatta a una distribuzione empirica assegnata. Si dimostra che 𝑉𝑎𝑟 [𝑌 ] 𝑏 ≅ 𝑚𝑦 − , 2(𝐸 [𝑌 ] − 𝑚𝑦 ) 𝜇 = ln (𝐸 [ 𝑌 ] − 𝑏 )2 √𝑉𝑎𝑟[𝑌 ] + (𝐸 [𝑌 ] − 𝑏 )2 𝜎 2 = ln ( , 𝑉𝑎𝑟 [𝑌 ] + 1) (𝐸 [ 𝑌 ] − 𝑏 )2 DIM. Si osserva che i parametri 𝜇, 𝜎 2 di 𝑌 sono gli stessi di 𝑋 = 𝑌 − 𝑏, pertanto OSSERVAZIONE 1 porge la tesi, non appena si considera che 𝐸 [𝑋] = 𝐸 [𝑌 ] − 𝑏. Possiamo poi scrivere 𝑚𝑦 = 𝑏 + 𝑒 𝜇 = 𝑏 + (𝐸 [ 𝑌 ] − 𝑏 )2 √𝑉𝑎𝑟 [𝑌 ] + (𝐸 [𝑌 ] − 𝑏 )2 2 ⇒ (𝑚𝑦 − 𝑏) = (𝐸 [ 𝑌 ] − 𝑏 )4 𝑉𝑎𝑟 [𝑌 ] + (𝐸 [𝑌 ] − 𝑏 )2 Quindi la dislocazione 𝑏 è la soluzione della equazione di quarto grado seguente: Eq. 3 2 (𝑚𝑦 − 𝑏) [𝑉𝑎𝑟 [𝑌 ] + (𝐸 [𝑌 ] − 𝑏 )2 ] − (𝐸 [𝑌 ] − 𝑏 )4 = 0 Questa equazione ammette due soluzioni complesse e una reale che tuttavia presenta una espressione molto complicata (risolto con Wolfram Mathematica). Seguendo le argomentazioni in [1], 2 consideriamo lo sviluppo di McLaurin di 𝑒 𝜎 , ovvero 2 𝑒𝜎 = 1 + 𝜎2 + 2 𝜎4 𝜎6 + +⋯ 3! 2 Se ci limitiamo alla approssimazione 𝑒 𝜎 = 1 + 𝜎 2 , per i momenti abbiamo 𝐸 [𝑌 ] ≅ 𝑒 𝜇 (1 + 𝜎2 ) + 𝑏, 2 𝑉𝑎𝑟 [𝑌 ] ≅ 𝑒 2𝜇 (1 + 𝜎 2 )𝜎 2 ≅ 𝑒 2𝜇 𝜎 2 Sostituendo 𝑒 𝜇 = 𝑚𝑦 − 𝑏 in queste equazioni, si ha 𝐸 [𝑌 ] ≅ (𝑚𝑦 − 𝑏 ) (1 + 𝜎2 ) + 𝑏, 2 Ricavando poi 𝜎 2 dalla seconda e sostituendo nella prima, si ha 5 2 𝑉𝑎𝑟 [𝑌 ] ≅ (𝑚𝑦 − 𝑏 ) 𝜎 2 𝐸 [𝑌 ] ≅ (𝑚𝑦 − 𝑏 ) (1 + 𝑉𝑎𝑟 [𝑌 ] 2(𝑚𝑦 − 𝑏 ) 2 ) + 𝑏 = 𝑚𝑦 + 𝑉𝑎𝑟 [𝑌 ] 2(𝑚𝑦 − 𝑏 ) ⇒ 2(𝑚𝑦 − 𝑏 )𝐸[𝑌] ≅ 𝑚𝑦 2(𝑚𝑦 − 𝑏 ) + 𝑉𝑎𝑟 [𝑌 ] ⇒ 𝑚𝑦 − 𝑏 ≅ E si trova così la formula nella tesi ∎ ⇒ 𝑉𝑎𝑟 [𝑌 ] 2(𝐸 [𝑌 ] − 𝑚𝑦 ) Figura 2. Legge lognormale a tre parametri con media 346.56, deviazione standard 152.19, e mediana 318.5 (mmHgL/min/𝑚2 ). Generata da CODICE 2. Si ha sia la densità ottenuta calcolando 𝑏 direttamente da Eq. 3, in modo numerico (linea continua), sia quella che si ottiene utilizzando la formula approssimata in PROP. 7, che in figura è indicata 𝑏0 (linea tratteggiata). ESEMPIO 2 (cardiac power index). L’indice di potenza cardiaca, dato dal prodotto tra la pressione arteriosa media e l’indice cardiaco (output cardiaco normalizzato rispetto la superficie corporea) segue una legge lognormale a tre parametri negli esseri umani [2]. I dati empirici presentano media 346.56, deviazione standard 152.19, e mediana 318.5 (mmHgL/min/𝑚 2 ). Si chiede di identificare i parametri 𝜇, 𝜎, 𝑏 della legge lognormale a tre parametri che possiede la medesima media, varianza, e mediana dei dati empirici. Utilizzando CODICE 2 si ottiene Figura 2, dove si ha sia la densità ottenuta calcolando 𝑏 direttamente da Eq. 3 (in modo numerico), sia quella che si ottiene utilizzando la formula approssimata in PROP. 7 (che in figura è indicata 𝑏0 ). Nel primo caso il programma restituisce 𝑏0 = 37.37, 𝜇 = 5.66, 𝜎 = 0.39; nel secondo caso sia ha 𝑏 = 115.53, 𝜇 = 5.31, 𝜎 = 0.51. Si vede pertanto che l’approssimazione non deve essere usata in modo acritico ed è comunque preferibile procedere con una soluzione numerica della equazione di quarto grado in 𝑏. CODICE 2 (metodo dei momenti). Il seguente codice in Octave ricava la densità lognormale a tre parametri a partire da media, varianza, e mediana assegnate. Genera Figura 2. Il metodo traccia il grafico sia ricavando 𝑏 direttamente da Eq. 3 (in modo numerico), sia quella che si ottiene utilizzando la formula approssimata in PROP. 7 (che in figura è indicata 𝑏0 ). 6 % file name = three_parameter_lognormal_var % date of creation = 27/09/2023 % % It plots the 3-parameter lognormal law from assigned mean, variance and median. % This script allows for negative values of b, b0. % clear all close all pkg load statistics % EY=343.35; VY=1.15e+04; my=342.27; steps=200; y(1)=my-2*sqrt(VY); y(steps)=my+5*sqrt(VY); delta_y=(y(steps)-y(1))/(steps-1); y=[y(1):delta_y:y(steps)]; % function [b,b0] = solve_b (EY,VY,my) b0 = my - 0.5*VY/(EY-my); b = fsolve(@(b)((my-b)^2)*(VY + (EY-b)^2) - (EY-b)^4,b0) endfunction % function [feY,b] = density_Y (y,EY,VY,my) [b,b0] = solve_b (EY,VY,my); if ((y-b)<=0) feY=0; else EYb=EY-b; MU=log((EYb^2)/sqrt(VY+EYb^2)); SIGMA=sqrt(log(1+VY/EYb^2)); A=-0.5*((log(y-b)-MU)/SIGMA)^2; C=1/(SIGMA*(y-b)*sqrt(2*pi)); feY=C*exp(A); endif endfunction % function [feY,b0] = density_Y0 (y,EY,VY,my) [b,b0] = solve_b (EY,VY,my); if ((y-b0)<=0) feY=0; else EYb=EY-b0; MU=log((EYb^2)/sqrt(VY+EYb^2)); SIGMA=sqrt(log(1+VY/EYb^2)); A=-0.5*((log(y-b0)-MU)/SIGMA)^2; C=1/(SIGMA*(y-b0)*sqrt(2*pi)); feY=C*exp(A); endif endfunction % for i=1:length(y) [density(i),b]=density_Y(y(i),EY,VY,my); endfor plot(y,density,'-k','Linewidth',1) 7 hold on for i=1:length(y) [density2(i),b0]=density_Y0(y(i),EY,VY,my); endfor plot(y,density2,'--k','Linewidth',1) title(strcat("b=",num2str(b),", b_{0}=",num2str(b0)),"fontsize",15) legend ('b','b_{0}','location','northeast',"fontsize",15) grid on 3 Appendice DEF. 3 (indice di asimmetria). Detto anche indice di skewness, è il momento centrato di ordine tre: Eq. 4 𝛾 [𝑋 ] = 𝐸[(𝑋−𝐸[𝑋])3] 𝜎3 = +∞ ∫−∞ (𝑥−𝐸[𝑋])3𝑓𝑋 (𝑋)𝑑𝑥 𝜎3 dove 𝜎 2 = 𝑉𝑎𝑟 [𝑋]. Come vedremo con il calcolo di questo indice per le densità note, un indice positivo indica la presenza di una coda a destra del grafico della densità, un indice negative indica la presenza di una coda a sinistra. Per il calcolo di 𝛾 si osservi che 𝛾𝜎 3 = 𝐸[𝑋 3 − 3𝑋 2 𝐸 [𝑋] + 3𝑋𝐸 2 [𝑋] − 𝐸 3 [𝑋]] = 𝐸 [𝑋 3 ] − 3𝐸 [𝑋]𝐸 [𝑋 2 ] + 3𝐸 [𝑋]𝐸 2 [𝑋] − 𝐸 3 [𝑋] = Eq. 5 𝛾 [𝑋 ] = = 𝐸 [𝑋 3 ] − 3𝐸 [𝑋](𝐸 [𝑋 2 ] − 𝐸 2 [𝑋]) − 𝐸 3 [𝑋] ⇒ 𝐸[𝑋 3 ]−3𝐸[𝑋]−𝐸3[𝑋] 𝜎3 Si osservi che definita la v.a. 𝑌 = 𝑎𝑋 + 𝑏 si ha 𝐸 [𝑌 3 ] = 𝐸 [𝑎3 𝑋 3 + 3𝑎2 𝑋 2 𝑏 + 3𝑎𝑋𝑏 2 + 𝑏 3 ] = 𝑎3 𝐸 [𝑋 3 ] + 3𝑎2 𝑏𝐸 [𝑋 2 ] + 3𝑎𝑏 2 𝐸 [𝑋] + 𝑏 3 −3𝐸 [𝑌 ] = −3𝑎𝐸 [𝑋] − 3𝑏 −𝐸 3 [𝑌 ] = −𝑎3 𝐸 3 [𝑋] − 3𝑎2 𝑏𝐸 2 [𝑋] − 3𝑎𝑏 2 𝐸 [𝑋] − 𝑏 3 Sommando membro a membro Eq. 6 4 𝛾[𝑌 ]𝜎 3 = 𝑎3 𝐸 [𝑋 3 ] + 3𝑎2 𝑏𝐸 [𝑋 2 ] − 3𝑎𝐸 [𝑋] − 3𝑏 − 𝑎3 𝐸 3 [𝑋] − 3𝑎2 𝑏𝐸 2 [𝑋] ⇒ 𝛾[𝑎𝑋 + 𝑏 ] = Riferimenti 𝑎 3𝐸[𝑋 3 ]+3𝑎2𝑏𝐸[𝑋 2]−3𝑎𝐸[𝑋]−3𝑏−𝑎 3𝐸 3 [𝑋]−3𝑎2𝑏𝐸2[𝑋] 𝜎3 1. Sangal BP. et Biswas AK. “The 3 parameter Lognormal Distribution and its Application in Hydrology” 1970 2. Maccallini Paolo, “Probability Density of Left-Ventricular Mechanical Power: derivation from metadata” 2023 8