Etude de la distribution des sauts d'amplitude dans une image numérique

Audran Le Baron

Description

Brique BASTA´Etude de cas :´Etude de la distribution des sauts d’amplitude dans une imagenum´eriqueAudran Le Baron17 d´ecembre 2002Fig. 1 – Image de travailLes donn´eesNoustravailleronstoutaulongdecette´etudedecassurl’imageaero.tiff,image en 256 niveaux de gris de taille 510×341. Nous chargeons cette imageet l’affichons par le code MatLab suivant. L’image est repr´esent´ee figure 1.a=imread(’aero.tiff’);imagesc(a);colormap gray;colorbar;axis equal;Nous travaillerons plus particuli`erement sur le vecteur des sauts d’ampli-tude horizontaux, not´e x, obtenu par le code suivant :a = double(a);[M, N] = size(a);h = a(:,2:end) - a(:,1:end-1); % Matrice des sauts horizontaux% (taille M x N-1)x = reshape(h, prod(size(h)), 1); % Matrice aplatie en vecteurLa figure 2 repr´esente les sauts d’amplitude horizontaux de cette image.1−50150 050 100 10050 1500 200−50 250300−100350−15040050 100 150 200 250 300 350 400 450 500Fig. 2 – Sauts d’amplitude horizontaux1 Mod´elisation par une gaussienneConsid´erons pour commencer le mod`ele param´etrique dans lequel la den-sit´edenotrevariableal´eatoire(valeurdusautd’amplitude)suituneloigaus-sienne : 21 (x−μ)g(x) = √ exp −22 2σ2πσ1.1 Comparaison d’histogrammesLe calcul de la moyenne et de la variance des observations se fait par lecode suivant :m = mean(x);v = var(x);Les r´esultats sont :μ = 0,01572σb = 337,3939Par le code suivant, on trace la distribution gaussienne dont les ...

Informations

• Publié par: Unje
• Nombre de lectures: 42
• Langue: Français

Extrait

´ Etude

de

la

Brique BASTA

´ Etude de cas :

distribution des sauts d’amplitude ´ri nume que

Audran Le Baron

17d´ecembre2002

dans

une

image

Fig. 1 – Image de travail

Les donnees ´ Noustravailleronstoutaulongdecettee´tudedecassurl’image aero.tiff , image en 256 niveaux de gris de taille 510 × 341. Nous chargeons cette image etl’affichonsparlecodeMatLabsuivant.L’imageestrepre´sent´eefigure1. a=imread(’aero.tiff’); imagesc(a); colormap gray; colorbar; axis equal; Noustravailleronsplusparticulie`rementsurlevecteurdessautsd’ampli- tudehorizontaux,note´ x , obtenu par le code suivant : a = double(a); [M, N] = size(a); h = a(:,2:end) - a(:,1:end-1); % Matrice des sauts horizontaux % (taille M x N-1) x = reshape(h, prod(size(h)), 1); % Matrice aplatie en vecteur La fi re 2 r ´ te les sauts d’amplitude horizontaux de cette image. gu epresen

1

Fig. 2 – Sauts d’amplitude horizontaux

1Mod´elisationparunegaussienne Consid´eronspourcommencerlemod`eleparam´etriquedanslequelladen- site´denotrevariableale´atoire(valeurdusautd’amplitude)suituneloigaus- sienne : g ( x ) = √ 2 1 πσ 2 exp  − ( x 2 − σ 2 µ ) 2  1.1 Comparaison d’histogrammes Le calcul de la moyenne et de la variance des observations se fait par le code suivant : m = mean(x); v = var(x); Le ´ ltats sont : s resu

µ = 0 , 0157 b σ 2 = 337 , 3939 Par le code suivant, on trace la distribution gaussienne dont les para- me`tressontlesmoyenneetvarianceempiriques,avecl’histogrammedesob- servations.Ler´esultatn’estpastr`essatisfaisantcommelemontrelafigure ?? .

2

Fig. 3 – Histogramme des observations et distribution gaussienne

x range = [-200:.1:200]; _ gauss_m_v = 1/(sqrt(2*pi*v)) * exp(-(x_range-m).^2 ./ (2*v)); hist(x, 100); hold on plot(x_range, length(x) * gauss_m_v, ’r’); hold off;

1.2 Comparaison de quantiles Soit F µ σ lafonctiondere´partitiond’unedistributiongaussiennedepa- , ram`etres( µ, σ ). Nous noterons F = F 1 , 1 . Le quantile d’ordre p de g µ,σ est le x tel que :

t − µ ) 2 p = F µ,σ ( x ) = Z − x ∞ √ 2 1 πσ 2 exp  − ( 2 σ 2  dt 1 x − µ t 2 dt = Z −∞ σ exp − 2 σ = σ 1 F  x σ − µ 

3

Fig. 4 – Quantiles d’une loi gaussienne en fonction des quantiles d’une autre loi gaussienne

Ainsi, si q µ,σ ( p ) est le quantile d’ordre p d’une gausienne ( µ, σ ) et q ( p ) celuid’unegaussiennecentre´er´eduite,ona: q µ,σ ( p ) = σq ( σp ) + µ C’est pourquoi, lorsque l’on trace les quantiles d’une variable gaussienne par rapport aux quantiles d’une seconde variable gaussienne, on obtient une droite ( Cf. figure 4). Ce n’est pas le cas lorsque l’on trace les quantiles des observations ( Cf. figure 5) en fonction des quantile d’une loi normale par le code suivant, comme le montre la figure 6. g_m_v = m + sqrt(v) * randn(length(x), 1); qqnorm(g1); % C’est bien une droite

qqnorm(x); % Ce n’est PAS une droite On note qu’il y a trop de poids dans la queue de la gaussienne et pas assez au ✭ centre ✮ .

4

Fig. male

6

–

Quantiles

Fig. 5 – Quantiles des observations

des observations en fonction des quantiles

5

d’une

loi

nor-

2Mod´elisationparuneloi ✭ laplaciennege´- ne´ralisee ✮ ´ Changeonsdemode`leparam´etriquepourprendreceluid’uneloideLa- placege´n´eralise´e,dedensit´e: αη = h ( x ) 2Γ(1 /α ) exp( − | ηx | α ) 2.1Estimationparlam´ethodedesmoments Unepremi`ereme´thoded’estimationestlam´ethodedesmoments.Nous avonsdeuxparam`etresscalaires( α et η ), donc deux moments qui sont ici la varianceetlekurtosis,donne´spar: σ 2 = η 2 Γ Γ   3 α α 1   Γ  κ = α 1  Γ  α 5  Γ  3 α  2 Le code suivant permet de tracer κ ( α ) et de calculer le kurtosis empirique κ . b alpha_vect = [25:100]/50; k = kurt(alpha_vect); plot(alpha_vect, k); k_emp = kurtosis(x); κ ( α )estrepre´sente´figure7etlecalculdukurtosisempiriquedonn b e κ = 17 , 8968.Doncd’apre`slacourbe, α se trouverait entre 0 , 5 et 0 , 6. D’o`ulare´solutionnume´riquesuivante: alpha_vect = [0.5:.0001:0.6]; [trash, indice] = min(abs(kurt(alpha_vect) - k_emp)) alpha = alpha_vect(indice) eta = 1 / sqrt(v * gamma(1/alpha) / gamma(3/alpha)) Lere´sultatdececalculdonne: α = 0 , 5672 η = 0 , 3469

6

en ; d 0]20 ap_ltesi [=egnar_:1.:002- mma(2*gapha)1/ala=plalecat(/ahe*an_r).gelp^a);haxe*)a-(pe(sbx*at t(2*pi*v))*exp(-agsu_s_m=v/1s(rq 01,x(tsih ;))v*2(./^2).-mgean_r(x x(*)gnht,eelargnt(x_ n plo oldo 0); h r_x(egna )’r’ tolp pllae,acnsdee_it ;) ohdlfo m_v_’,’gx)*gauss,length( Lesuivcodertiaftnanese´rpeleurrsteraegemmˆhpqieuC(.fgfiru8e) l’histogrammedesbosavrenoital,sstdiburiontiusgaee`aouv´netrsien e´´cnorp -c it ales ristdilaetteenedgecalpaLednoitubis´eeave´en´eral`mteer e slcseapar laaretm´imstsp´enemo x .stedohmsed

Fig. 7–Repr´esentationde κ ( α )

7

i=1:100 :1]1 f ; ro .0100:0.an_r=[ge halp; p ))ahpla/1,egnaa^.) daleC.snoitavres pr ullccaundeanemiteluqnadnsetcoiesobuesdiriqsemp ecedletttnauselieeac on n f deoiplLaeixuocpm ;f fiAdnmeracelesqarer,ont gean_r(xmaam*g).(sba(cnir_x*.ate (x)]:max arti ; rep1=2/itnoisng(*+1 letielsdoial _r . xegnaim[=)x(n1.0: ´ealabledesapproixamitnodsseuqna

Fig. 8 – Observations, gaussienne et Laplace (moments)

[temp, indice] = min(abs(repartition - p_range(i))); quant_emp(i) = x_range(indice); end plot(quantile(x, p_range), quant_emp); Onconstateunecertaineame´liorationparrapportaumode`legaussien (voir figure 9). Il y a cependant toujours trop de poids dans la queue de la distribution de Laplace.

2.2Estimationparlam´ethodedumaximumdevrai- semblance Uneautremethoded’estimationdesparam`etresestlam´ethodedumaxi-´ mumdevraisemblance.Ici,lesdeuxe´quationsdevraisemblancesont: n 1 /α + Ψ(1 /α ) /α 2 − n 1 X log( n | X i | ) | ηX i | α i =1 1 η =
α i = n X 1 | X i | α ! − α n Parite´rationssuccessivesetgrˆaceaucodesuivant,oncalcule˜ α et η ˜.

8

Fig. 9 – Quantiles de la loi de Laplace en fonction des quantiles empiriques (moments)

% Suppression des zeros [trash1 trash2 y] = find(x);

% Initialisation n = length(x);

alpha = 0.5672; eta = 0.3469;

pas = 0.01; ecart = 0.5;

logvrais = -Inf;

% Boucle alpha_vect = [alpha - ecart : pas : alpha + ecart]; for i=1:30 [trash indice] = min(abs(eq4(y, alpha_vect, eta, n))); alpha = alpha_vect(indice); eta = power(alpha / n * sum (power(abs(x), alpha)), -1/alpha); newlogvrais = sum(log(laplace(x, alpha, eta)));

9

´ Fig. 10–Evolutiondelalog-vraisemblancea`chaqueit´eration

i alpha eta newlogvrais (newlogvrais - logvrais)/abs(newlogvrais) logvrais = newlogvrais; logvraisvect(i) = logvrais; end ; plot(logvraisvect); L’e´volutiondelalog-vraisemblaceaufuret`amesuredesit´erationsest ´nte´figure10.Aubutde18it´erations,lalog-vraisemblancesestabilise represe o pour les valeurs :

˜ 0 , 4072 α = ˜ 1 , 7715 η = Onpeutalors,pourcesnouvellesvaleursdesparam`etres,retracerladis- tribution et les quantiles. Ceci est fait aux figures 11, 12. On constat ` e a nouveauuneame´lioration,maislemode`len’esttoujourspasparfait.

10