HYDROLOGIE STATISTIQUE ENSEEIHT Hydraulique et Mécanique des

Document Sample
HYDROLOGIE STATISTIQUE ENSEEIHT Hydraulique et Mécanique des Powered By Docstoc
					Hydrologie statistique 2007-2008 M. Ababou

                                     Idris Abdou

                                  Yohan Allouche

                            Renaud Champredonde




        HYDROLOGIE STATISTIQUE
         ENSEEIHT -- Hydraulique et
           Mécanique des Fluides
        3Hy SE , Mastère Hy (2007-08)


  Contrôle d' Hydrologie Statistique




                                                   1/27
Hydrologie statistique 2007-2008 M. Ababou

Table des matières :


Exo 1. Exercice de reconstitution de données par régression
linéaire (pluies mensuelles à Mens et Roissard, bassin du
Drac, Alpes).....................................................................................3



TD 1 : Hydrologie statistique – analyse univariee – gumbel &
poisson : Crues annuelle & crues extrêmes de la Garonne à
Toulouse........................................................................................10
        Extension du TD ..............................................................13



TD2 : Exercice 2 : A.C.P. (Analyse en Composantes
Principales)....................................................................................16



TD3 - Bureau d’Etudes « Hyd.Stat.», Jan.2007. Analyse
statistique de données de pluies & de débits............................22




                                                                                              2/27
Hydrologie statistique 2007-2008 M. Ababou

 Exo 1. Exercice de reconstitution de données par
 régression linéaire (pluies mensuelles à Mens et
         Roissard, bassin du Drac, Alpes)
Reconstitution des pluies mensuelles de Mens à partir de celles de
Roissard:
    ● Mars 1946: vérifier et compléter la réponse du polycopié, si possible en refaisant
       vous-même le calcul de reconstitution, et donner aussi l'intervalle de confiance à
       90% (à partir d'une table de la F de répartition de Gauss);

   ●    Mars 1940: procéder de même... obtenir (i) la pluie reconstituée, et (ii) l'intervalle
        de confiance à 90

Pluies mensuelles observées sur deux stations de jaugeages dans les Alpes.

Le tableau de pluies mensuelles fourni sur ces deux stations alpines comportent des
données manquantes en mars 1946 et 1940 à la station de Mens.




       Illustration 1: Histogramme de Mens et Roissard en fonction du temps



a) Si on avait pas la deuxième station on aurait corréler Y(t) de mens avec le temps

b) Ici comme on a une deuxième station (Roissard) on profite pour faire une estimation
aléatoire.
Mens étant conditionné par Roissard on calque en terme de pluviométrie l’évolution de la
station de Mens sur celle de Roissard. Ainsi on estime que quand la pluviométrie est
importante sur Roissard, elle l’est aussi sur Mens et vice versa. Ce qui nous donne de
tendance.




                                                                                         3/27
Hydrologie statistique 2007-2008 M. Ababou




             Illustration 2: Tendance des pluviométries de Mens et Roissard



   a) Analyse des données par extrapolation graphique
Histogramme de Mens et Roissard en fonction du temps :

                                              Histogramme données Roissard

                                        200
                   Pluviométrie en mm




                                        150

                                        100

                                        50

                                         0
                                        19 8
                                        19 9
                                        19 0
                                        19 1
                                        19 2
                                        19 3
                                        19 4
                                        19 5
                                        19 6
                                        19 7
                                        19 8
                                        19 9
                                        19 1
                                        19 2
                                        19 3
                                        19 4
                                        19 5
                                        19 7
                                           76
                                           2
                                           2
                                           3
                                           3
                                           3
                                           3
                                           3
                                           3
                                           3
                                           3
                                           3
                                           3
                                           4
                                           4
                                           4
                                           4
                                           4
                                           4
                                        19




                                                            Années




                                                                              4/27
Hydrologie statistique 2007-2008 M. Ababou

                                                                         Histogramme données Mens
                                                  Histogramme données Mens - Roissard


                                       Pluviométrie en mm
                                                            160
                                                            140
                                                            120
                                                            100
                           200                               80
      Pluviométrie en mm




                                                             60
                           180                               40
                           160                               20
                                                              0
                           140
                                                                 19 9
                                                                 19 0




                                                                 19 4
                                                                 19 5




                                                                 19 9
                                                                 19 1




                                                                 19 5
                                                                 19 7
                                                                 19 8




                                                                 19 1
                                                                 19 2
                                                                 19 3




                                                                 19 6
                                                                 19 7
                                                                 19 8




                                                                 19 2
                                                                 19 3
                                                                 19 4



                                                                   76
                           120
                                                                   2
                                                                   2
                                                                   3
                                                                   3
                                                                   3
                                                                   3
                                                                   3
                                                                   3
                                                                   3
                                                                   3
                                                                   3
                                                                   3
                                                                   4
                                                                   4
                                                                   4
                                                                   4
                                                                   4
                                                                   4
                                                             19




                                                                                                                               Série1
                           100                                                          Années
                            80                                                                                                 Série2
                            60
                            40
                            20
                             0
                           19 8
                           19 9
                           19 0
                           19 1
                           19 2
                           19 3
                           19 4
                           19 5
                           19 6
                           19 7
                           19 8
                           19 9
                           19 0
                           19 1
                           19 2
                           19 3
                           19 4
                           19 5
                           19 6
                           19 7
                              76
                              2
                              2
                              3
                              3
                              3
                              3
                              3
                              3
                              3
                              3
                              3
                              3
                              4
                              4
                              4
                              4
                              4
                              4
                              4
                              4
                           19




                                                                                    Années



                                                                     Regression linéaire de Mens / Roissard
                                 Racine (pluvio Roissard)




                                                            16                                      y = 0,7367x + 3,2264
                                                            14                                           R2 = 0,5861
                                                            12
                                                            10                                             Série1
                                                             8
                                                             6                                             Linéaire (Série1)
                                                             4
                                                             2
                                                             0
                                                                 0          5           10       15
                                                                         Racine (pluvio Mens)




   a) Si on avait pas la deuxième station on aurait corréler Y(t) de mens avec le temps
   b) Ici comme on a une deuxième station (Roissard) on profite pour faire une estimation
      aléatoire


   Mens étant conditionnée par Roissard on calque en terme de pluviométrie l’évolution de la
   station de Mens sur celle de Roissard. Ainsi on estime que quand la pluviométrie est importante
   sur Roissard, elle l’est aussi sur Mens et vice versa. Ce qui nous donne de tendance.




                                                                                                                                        5/27
Hydrologie statistique 2007-2008 M. Ababou

   b) Résolution analytique :
Le travail consiste à la reconstitution des données de pluies à Mens à partir de celle de Roissard à la
même période de l’année.
D’abord on annule les deux points manquants et on applique l’analyse sur les points connus.

Comme on connaît les données sur la station N°2 alors on fait la régression linéaire de X1/X2.
Soit    Y= X1= √ (P1) : variable à expliquer (Mens)
        X= X2= √ (P2) : variable explicative (Roissard)

On trace le tableau des pluviométries donné


                 Y=√(P1)
                                                                   Y= a X +b = 0.974*X – 0.800
                                                                   a = ρ* (σ1 / σ2)
                                                                   b = E(Y) – E(a X) =my – a*mx
                                                e

                                                              X= sqrt(P2)



                         P1(mm)= [0.974 sqrt(P2) – 0.800] ²

                       A.N     Pour P_Roissard = 60 mm on a P_Mens = 45,49 mm

                                     P_Roissard = 50 mm on a P_Mens = 37,05 mm




   Méthode

       On prends l’approximation Gaussienne la meilleure

   Moments de Mens X= √ P1 : m1 = 6.7            en √ mm      et    σ1 = 2.9   en √ mm      Moments de
   Roissard Y= √ P2 : m2 = 7.7 en √ mm et σ2 = 2.8 en √ mm                              a = ρ* (σ1 / σ2) =
   0.974 et b = E(Y) – E(a X) =my – a*mx = -0.8 √ mm

   La droite de régression introduit une erreur ( dite d’estimation)

   Par définition ^Y= a X + b est la moyenne conditionnée sachant que ^Y= E (Y/X) : théorème
   de Bayes et c’est le meilleur estimateur ( Min   σe ²)           si X, Y Gaussien.

                                                                                                      6/27
Hydrologie statistique 2007-2008 M. Ababou

                          Erreur       e = Y - ^Y est une V.A
                   σe = √ (1 - ρ) * σy
                         ρ_x;y = 0.94 ( Mens – Roissard )


ε rms = σe    / σy = √ (1 - ρ²)

                  A.N         ε rms = 0.34   Root means square error
                              Erreur relative à la normalisation

Intervalle de confiance hachurée à I90% directement lié à σe :
L’objectif est de caractériser une intervalle de confiance autour de la moyenne telle que la variable
aléatoire ait une probabilité P = 90% d’appartenir à cette intervalle à partir d’une table de la la
fonction de répartition F(u) d’une variable gaussienne centrée réduite [N (0,1)]

Table de la FdR de Gauss

F(0)         F(0.25)      F(0.52)       F(0.84)    F(1.28)      F(1.64)   F(2.32)     F(2.57)
0.50         0.60         0.70          0.80       0.90         0.95      0.99        0.995
L’intervalle à 90% de Probabilité (I90 %) est :
Proba ( U≤ 1.64) = 0.95
Proba ( U≥ -1.64) = 0.05 par symétrie de la loi
Donc
Proba (-1.68 ≤ U ≤ 1.64 ) ⇒ I90 % = [-1.68 ; +1.68 ] pour la valeur aléatoire centrée réduite
N(0,1)
Or     X = mx + σx*U
       U = 1.64        pour    I90 % = ^Y +U* σe

On obtient donc la variable aléatoire X gaussienne N[mx, σx^2]
L’objectif est de générer un des points ( étoiles bleues) qui appartienne à l’intervalle I90% qu’on va
corréler par σe puis reconstituer Yn :


                                    Yn_ reconstitué = a Xn +b + (σe * G)

G est une réplique d’ une V.A N[0,1]

N[0,1] : Générateur d’une V.A gaussienne réduite
Cette méthode est intéressante car σx et ρ interviennent alors que la première méthode (dite de
glissement) introduit beaucoup d’erreur (notamment dans le calcul de variance).




                                                                                                  7/27
Hydrologie statistique 2007-2008 M. Ababou




      Y=√ (P1)




                                             e




                                                 X= √ (P2)




                                                             8/27
Hydrologie statistique 2007-2008 M. Ababou

Figure 1 . Tableau pluies mensuelles de mars avril à Mens & Roissard (BV du Drac, Alpes)
                   Pluies Mensuelles en 2 stations d'un Bassin Versant du Drac (de 1928 à 1947, et en 1976)


                                    S1 Mens                                            S2 Roissard


          Années                Mars                       Avril                     Mars                     Avril


            1928                    61                        84                        44                     132


            1929                     7                        65                         3                      79


            1930                   109                        53                       135                     115


            1931                    90                        40                       116                      57


            1932                    59                        67                       101                      89


            1933                    33                        21                        83                      44


            1934                    74                       135                        88                     130


            1935                    41                        18                        91                     131


            1936                    56                       132                        64                     132


            1937                   143                        56                       188                      78


            1938                     3                        19                         3                       7


            1939                    53                        91                        86                      92


          1940                37.05                           x                       50                      112

            1941                    45                        83                        55                     117


            1942                    19                        23                        40                      42


            1943                     8                        25                        12                      35


            1944                    19                        30                        20                      30


            1945                    19                        17                        18                      18


          1946                45.49                           x                       60                       44

            1947                   103                        35                       134                      31


            1976                    57                        60                        62                      65




                                                                                                                      9/27
Hydrologie statistique 2007-2008 M. Ababou

HYDROLOGIE STATISTIQUE – TD1:
              ANALYSE UNIVARIEE – GUMBEL & POISSON :
 CRUES ANNUELLE & CRUES EXTREMES DE LA GARONNE A TOULOUSE

   1. Quelle est la variable hydrologique étudiée (expliquez le terme
      crue)?

La variable hydrologique étudiée est la hauteur d'eau. En utilisant la courbe de tarrage on
peut obtenir facilement les débits.

  2. Utilisez la FdR proposée pour obtenir la crue annuelle centennale
(expliquez). Question subsidiaire: est-ce une loi de Gumbel ?
(paramètres=?)

 Pour cela il convient au préalable de déterminer s’il s’agit d’une loi de Gumbel et s’il
 en est le cas, d’identifier les paramètres [ α, β ] qui caractérisent la loi de Gumbel.

          a) La représentation graphique laisse (à priori) à penser qu’il pourrait s’agir d’une loi
             de Gumbel. En effet nous avons là :
             une représentation à l’envers du graphe de la FdR: H=f [ FdR (H) ]
             On lit sur l’axe des ordonnées la hauteur et sur l’axe des abscisses la FdR(H).
             On constate également que ces deux axes ont subi une transformation qui reste à
             déterminer.
          b) Pour confirmer qu’il s’agit bien d’une loi de Gumbel nous proposons la méthode qui
             consiste à tester 3 points appartenant à la séquence des hauteurs tracés. Si ces trois
             points sont alignés par suite de la transformation alors nous pourrons affirmer qu’il
             s’agit d’une loi de Gumbel.

          Au préalable il convient de caractériser la transformation effectuée sur les axes.
          En faisant l’hypothèse que la sus transformation est une relation affine
          (type y= ax +b) lorsque F(H) est une loi de Gumbel


      1ère Etape : Détermination des paramètres [ α, β ]
                   Remarque [ α, β ] sont en mètre car la hauteur H est en mètre
                   F(H)= exp[-exp[- (H-α) /β]] on pose Y = (H-α) /β
      On détermine la pente (en prenant deux points éloignés sur la droite) et en déduit       β=
      0.9119
      Sachant β et un couple de point (x =F(H), y= H) ,
                      on en déduit α= H - β*[- ln [- ln F]] = 1.9822

      2ème étape: Détermination des paramètres (mH , σH) :
                moments et écart type

                                                                                                10/27
Hydrologie statistique 2007-2008 M. Ababou

                             σH = 1.28 * β= 1.1672
                             mH = α + 0.45* σH = 2.5075


      3ème étape : Test des points choisis
                  On prends un trio des points ( H1=2m, H2= 4, H3=6m) et détermine F(Hi) .
                 On constate que ces points se retrouvent bien sur le diagramme de Gumbel
                  donné dont l’équation est de type Hi=β*F(Hi) + α

      c) La crue annuelle centennale :
                  On utilise la loi de Gumbel H= - β*[- ln [- ln F]] + α = 6.2 m



   3°) La probabilité d’observer n crues supérieures à la crue centennale sur
   une période d’observation de 225 années
   Il s’agit là d’un évènement rare donc on utilisera la loi de Poisson sur un temps de retour
   Tr=100 ans. Cette loi représente la probabilité d’observer n dépassement de la crue de temps de
   retour Tr=100 précédemment calculé (H100). Voir 2°)
   La loi de poisson s’écrit :

                             Proba(k=n) = ((µ*TD)^n / n!) * exp (-µ*TD)
                             µ =1/ TR
                             TR = 100 ans et TD = 225 ans


       a) Au moins une :
                                             P0 = 0.105
                                             Proba( k≥1) = 1 - P0 = 89,4%

       b) Au moins deux:
                                             P1=0.237
                                             Proba( k≥2) = 1 - P0 - P1 = 65,7%
       c) Au moins trois :
                                             P2 = 0.266
                                             Proba( k≥3) =1 - P0 - P1 - P2 = 39.1%




                                                                                             11/27
Hydrologie statistique 2007-2008 M. Ababou

   Fichier.m Matlab du TD1 avec Résultat

   clear all;
   close all;

   %
   % Question 0 : Est-ce une loi de Gumbel ?
   %
   alpha = 4;
   beta = 1;
   h1 = 2;
   h2 = 4;
   h3 = 6;
   f1 = exp(-exp(-(h1-alpha)/beta));
   f2 = exp(-exp(-(h2-alpha)/beta));
   f3 = exp(-exp(-(h3-alpha)/beta));
   f_transfo1 = -log(-log(f1));
   f_transfo2 = -log(-log(f2));
   f_transfo3 = -log(-log(f3));
   f_transfo = [f_transfo1,f_transfo2,f_transfo3];

   %
   % Question 1 : Identifier alpha et beta
   %
   % Pente réelle en la calculant sur le graphique
   beta_reel = ((8.25-2)/(-log(-log(0.998))+log(-log(0.15))))

   % Expression de alpha
   alpha_reel = 3.35-beta_reel*(-log(-log(0.80)))

   % D’où la hauteur de la crue centennale
   hauteur_crue = ((-log(-log(0.99)))*beta_reel)+alpha_reel




                                                                12/27
Hydrologie statistique 2007-2008 M. Ababou

Extension TD 1 : Résolution du problème par la méthode Numérique :

Etape 1 : Numérisation des données de la Garonne

Etape 2 : Mise en place des deux vecteurs T_date et H_crue

Etape 3 : Calcul des moments d'ordre 1 des hauteurs :

Etape 4 : Classement des hauteurs par ordre croissant

Etape 5 : Fonction de répartition de la hauteur

         Expression de Hazen

Etape 6 : Calcul de la hauteur centennale

Etape 7 : Ajustement par une loi de Gumbel - relation paramètres/moments de la loi de
Gumbel

Etape 8 : Probabilité de non dépassement de la hauteur centennale (crue extrême)




                                                                                   13/27
Hydrologie statistique 2007-2008 M. Ababou

Fichier.m Matlab de Extension du TD1 avec Résultat

clear all;
close all;

%
% Question extension
%

% Chargement des données
load -ascii 'numerisation_donnees_garonne';
donnees = numerisation_donnees_garonne;
date = donnees(:,1);
hauteurs = donnees(:,2);
debits = donnees(:,3);
N =length(hauteurs);

% Calcul des moments d'ordre 1 des hauteurs :
moy_h = mean(hauteurs);
i = 1;
moment = zeros(N,1);
for i=1:N;
   moment(i) = mean(hauteurs(i)-moy_h);
end

% Classement des hauteurs par ordre croissant
[hauteurs_c i_c] = sort(hauteurs);
date_c = date(i_c);

% Fonction de répartition de la hauteur
% Expression de Hazen
j = 1;
fx = zeros(N,1);
for j = 1:N;
    fx(j)=(j-0.5)/N;
end
plot(hauteurs_c,fx)
grid on
xlabel('Hauteurs')
ylabel('Fx')
title('Fonction de répartition empirique de la Hauteur (Formule de Hazen)')

% Calcul de la hauteur centennale
ecart_type_h = std(hauteurs); % Ecart-type sur les hauteurs
alpha = moy_h-0.45*ecart_type_h;
beta = ecart_type_h/1.28;
hauteur_crue = ((-log(-log(0.99)))*beta)+alpha

                                                                              14/27
Hydrologie statistique 2007-2008 M. Ababou


% Ajustement par une loi de Gumbel - relation paramètres/moments de la loi
% de Gumbel
a = pi / (sqrt(6)*ecart_type_h);
h0 = moy_h - 0.045*ecart_type_h;

% Probalite de non dépassement de la hauteur centennale (crue extreme)
f_de_hauteur_crue = exp(-exp(-a*(hauteur_crue-h0)))


Ce script aboutit à la représentation suivante pour la fonction de répartition empirique de Hazen :

                                       F o n c t io n d e ré p a r t it io n e m p iriq u e d e la H a u t e u r ( F o rm u le d e H a z e n )
                        1


                     0 .9


                     0 .8


                     0 .7


                     0 .6


                     0 .5
                Fx




                     0 .4


                     0 .3


                     0 .2


                     0 .1


                        0
                            1      2              3                   4                  5                   6                   7               8   9
                                                                                   H a u te u rs



Nous obtenons ainsi les résultats suivants :

moy_h =
 3.0503

hauteur_crue =
  6.6054 mètres

a=
  1.1342

h0 =
  2.9994

f_de_hauteur_crue =             0.9834


                                                                                                                                                         15/27
Hydrologie statistique 2007-2008 M. Ababou

                        TD2-Exo.2 : A.C.P.
              (Analyse en Composantes Principales)

1.1 Question préalable : quelle est la signification de la variable
hydrologique analysée (débit « Q » -- ou débit spécifique « q ») ? A quel
type de normalisation des débits cela correspond-il ?

Il s'agit des écoulements mensuels en 6 stations hygrométriques (probablement des
exutoires de bassins versants. Dans ce cas ces mesures représentent des « quantités de
pluies » relevées en mm sur les exutoires du bassin versant, soit des mm par unité d'aire.

1.2 Moments simples. Calculer la moyenne, la variance et l’écart-type de
chaque variable (en utilisant directement les données, ou bien encore, les
sommes Σ données en annexe).
Scripte matlab :

load -ascii 'donnees.csv';

donnees = donnees;
annee = donnees(:,1);
observation = donnees(:,2);
naguilhes = donnees(:,3);
lanoux = donnees(:,4);
izourt = donnees(:,5);
gnioure = donnees(:,6);
caillaouas = donnees(:,7);
bleu = donnees(:,8);

% Somme de l'ensemble des valeurs pour chacune des stations
somme1=sum(naguilhes);
somme2=sum(lanoux);
somme3=sum(izourt);
somme4=sum(gnioure);
somme5=sum(caillaouas);
somme6=sum(bleu);

% Moyenne des valeurs pour chacune des stations
moyenne1=mean(naguilhes);
moyenne2=mean(lanoux);
moyenne3=mean(izourt);
moyenne4=mean(gnioure);
moyenne5=mean(caillaouas);
moyenne6=mean(bleu);

% Ecart-type des valeurs pour chacune des stations
ecarttype1=std(naguilhes);
ecarttype2=std(lanoux);
ecarttype3=std(izourt);
ecarttype4=std(gnioure);


                                                                                    16/27
Hydrologie statistique 2007-2008 M. Ababou

ecarttype5=std(caillaouas);
ecarttype6=std(bleu);


% On centre les valeurs
for i=1:23;
    centre1(i,1)=(naguilhes(i,1) - moyenne1)/ecarttype1;
    centre2(i,1)=(lanoux(i,1) - moyenne2)/ecarttype2;
    centre3(i,1)=(izourt(i,1) - moyenne3)/ecarttype3;
    centre4(i,1)=(gnioure(i,1) - moyenne4)/ecarttype4;
    centre5(i,1)=(caillaouas(i,1) - moyenne5)/ecarttype5;
    centre6(i,1)=(bleu(i,1) - moyenne6)/ecarttype6;
end

% Variance des valeurs pour chacune des stations
variance1=var(naguilhes);
variance2=var(lanoux);
variance3=var(izourt);
variance4=var(gnioure);
variance5=var(caillaouas);
variance6=var(bleu);

Les valeurs obtenues par ces différentes commandes sont regroupées dans le tableau ci
dessous :

   Station     Naguilhes     Lanoux          Izourt   Gnioure   Caillaouas    Bleu
 Moyenne           377,6      298,5          394,2     401.2      327.2      176.7
 Ecart Type         93.2      71.7            78.9      85.3       115        97.8
  Variance         8696.4    5143.5          6228.9    7281.1     13229      9569.1

1.3 Matrice de corrélation. Calculer la matrice de corrélation (i.e., la
matrice de covariance des variables réduites). Remarques ?
Script Matlab :

X=[232       180       450      450          391      163
228          155       355      337          271      110
416          344       391      376          306      125
479          370       503      490          387      234
323          250       358      334          293      162
379          260       288      269          432      351
423          325       476      505          380      144
154          141       215      197          137      37
523          400       567      590          516      337
440          340       337      364          318      137
478          370       412      441          518      314
431          329       365      386          313      241
359          294       313      358          274      160
295          271       318      305          208      104
464          360       381      415          597      406

                                                                                 17/27
Hydrologie statistique 2007-2008 M. Ababou

366        285          451         428       228       139
472        353          478         489       377       223
383        310          396         404       215       66
370        320          423         449       242       95
417        359          403         447       372       181
334        238          393         400       197       87
447        370          471         459       348       170
273        242          322         335       205       78]
X=

 232    180    450     450   391   163
 228    155    355     337   271   110
 416    344    391     376   306   125
 479    370    503     490   387   234
 323    250    358     334   293   162
 379    260    288     269   432   351
 423    325    476     505   380   144
 154    141    215     197   137    37
 523    400    567     590   516   337
 440    340    337     364   318   137
 478    370    412     441   518   314
 431    329    365     386   313   241
 359    294    313     358   274   160
 295    271    318     305   208   104
 464    360    381     415   597   406
 366    285    451     428   228   139
 472    353    478     489   377   223
 383    310    396     404   215    66
 370    320    423     449   242    95
 417    359    403     447   372   181
 334    238    393     400   197    87
 447    370    471     459   348   170
 273    242    322     335   205    78

>> Rqq=corr(X)

Rqq =

  1.0000      0.9627    0.6363     0.7030   0.6708   0.6455
  0.9627      1.0000    0.5998     0.6817   0.5775   0.5277
  0.6363      0.5998    1.0000     0.9651   0.4732   0.3119
  0.7030      0.6817    0.9651     1.0000   0.5347   0.3566
  0.6708      0.5775    0.4732     0.5347   1.0000   0.9194
  0.6455      0.5277    0.3119     0.3566   0.9194   1.0000

>> Cqq=cov(X)

Cqq =


                                                              18/27
Hydrologie statistique 2007-2008 M. Ababou


 1.0e+004 *

  0.8696    0.6439      0.4683    0.5594   0.7195     0.5888
  0.6439    0.5144      0.3395    0.4172   0.4764     0.3702
  0.4683    0.3395      0.6229    0.6499   0.4296     0.2408
  0.5594    0.4172      0.6499    0.7281   0.5247     0.2976
  0.7195    0.4764      0.4296    0.5247   1.3229     1.0344
  0.5888    0.3702      0.2408    0.2976   1.0344     0.9569



1.4 Diagonalisation de la matrice de corrélation. Afin d’alléger les
calculs, on donne en annexe la matrice diagonale D et la matrice de
passage P. En déduire les valeurs propres, ainsi que les vecteurs propres ou
« composantes principales ».
>> [P,D]=eig(Rqq)

%matrice de passage

P=

  0.6457 -0.3174 0.2324 -0.4728 0.0129 0.4524
 -0.5854 0.1670 -0.2377 -0.6211 -0.0655 0.4277
 -0.3404 -0.4918 0.3251 0.3832 -0.4889 0.3883
  0.3124 0.6476 -0.1879 0.2833 -0.4431 0.4137
  0.0450 -0.3129 -0.6462 0.3498 0.4465 0.4010
 -0.1578 0.3351 0.5753 0.2033 0.6007 0.3601

D=

%Matrice diagonale

  0.0209    0    0    0    0    0
     0 0.0246    0    0    0    0
     0    0 0.0700    0    0    0
     0    0    0 0.5537    0    0
     0    0    0    0 1.1228    0
     0    0    0    0    0 4.2079

On obtient donc le tableau des Cps suivants :

     CP             1              2           3               4     5        6
   Valeur       4.2079           1.1228      0.5537       0.0700   0.0246   0.0209



                                                                                  19/27
Hydrologie statistique 2007-2008 M. Ababou


   2. Analyses et applications

2. Montrer que, dans le cas présent, la CP1 représente les six variables
   avec un poids à peu prés égal pour toutes.

La valeur CP1 est de 4,2079. Pour la représenter en terme de pourcentage d'information
on effectue :

%CP1=4,2079/6*100=70,13%

CP1 contient donc environ 70% de l'information portée par les variances des 6 stations
hydrologiques. Cependant cet axe ne sera pas un axe discriminant, en effet si on observe
la dernière colonne de la matrice de passage, representant les Pij de CP1, on remarque
que les valeurs sont très similaires. Il est donc fortement probable que les stations ne
soient pas discriminées par cet axe.

La discrimination des stations se fera en apportant de l'information complémentaire (CP2,
CP3 etc..).

2. Calculer le % de variance expliquée par les K premières CP, en faisant
   varier K de 1 à 6. En déduire que l’on ne perd que quelques %
   d’information en éliminant les CP4, CP5 et CP6.
On vient de calculer le % de variance pour CP1 (K=1).
Les % pour les CP suivantes sont :

CP2=1.1228/6*100=18,71%
CP3=0.5537 /6*100=9,23%
CP4=0.0700/6*100=1,17%
CP5=0.0246/6*100=0,41%
CP6=0.0209/6*100=0,35%

On obtient donc :

Somme CP1,CP2,CP3= 70,13+18,71+9,23=98,07%

Les trois premiers axes comprennent donc plus de 98% de l'information disponible, les
axes 4,5 et 6 semblent donc inutiles pour discriminer les stations. Il cependant signaler
que ces axes peuvent être utiles dans le cas d'étude de population très proches.

De plus il faut signaler que les axes à utiliser seront dépendants des objectifs de l'étude.
Si on cherche a discriminer spatialement des stations hydrologiques la selection des axes
choisis semble adaptée. Par contre si on essaye de determiner des types de stations,
l'axe CP1 sera peu utile, car représentatif de toutes les stations.


                                                                                      20/27
Hydrologie statistique 2007-2008 M. Ababou

2. La figure 1 représente les 6 stations de jaugeage de débits (variables
   1,…,6) dans le plan des (CP2,CP3). Y a-t-il des regroupements
   possibles ? Que pouvez en déduire ?




                      Naguilhes et Lanoux




                                                                      Izourt et Gnioure




                                                 Caillaouas et Bleu




   Figure 1 : Représentation graphique des 6 stations dans le plan (CP2,CP3)

Sur le graphique on observe très nettement trois couples de stations. Pour identifier ces
stations (comme reporté sur la figure), il suffit de repérer les “coordonnées” de chaque
station dans la matrice de passage.

Ce fort appariemment deux à deux des stations hydrologiques amènent à penser que ces
couples sont situés sur le même exutoire de bassin versant, ces stations sont donc
probablement relativement proche géographiquement l'une de l'autre à l'interieur d'un
couple.




                                                                                          21/27
Hydrologie statistique 2007-2008 M. Ababou

              TD3 - Bureau d’Etudes « Hyd.Stat.», Jan.2007
       ANALYSE STATISTIQUE DE DONNEES DE PLUIES & DE DÉBITS

1 Calculer et interpréter des fonctions de covariances temporelles (auto et
croisées) sur les chroniques hydrologiques données

L’objectif de cette partie consiste à étudier les fonctions de covariances temporelles sur
les chroniques hydrologiques de la source pyrénéenne Aliou. L’analyse de ces interactions
mettra en évidence le caractère karstique de cette source. La première étape s’intéresse à
la visualisation des données pluviométriques en parallèle avec l’hydrogramme associé
correspondant à une sortie du bassin versant. Le script Matlab permettant la visualisation
des variables hydrologiques est fournie ci-dessous :

clear all;

close all;

load -ascii 'PtQt_ALIOU_de75a95_3cols.txt';

tpq= PtQt_ALIOU_de75a95_3cols;

t=tpq(:,1);

p=tpq(:,2);

q=tpq(:,3);

subplot(4,1,1)

plot(p)

grid on

axis ij

xlabel('Temps en jours')

ylabel('Pluie en mm')

title('Hyétogramme')

hold on;

subplot(4,1,2)

plot(t,q)

xlabel('Temps en jours')

ylabel('Débit en m^3/s')


                                                                                    22/27
Hydrologie statistique 2007-2008 M. Ababou

title('Hydrogramme')

hold off;




                  Figure 2 : hydrogramme et hyétogramme d'Aliou


Ce premier graphique met en évidence une certaine corrélation visuelle entre les deux
variables tracées. En effet la quantité d’eau tombée, tracée sur le graphique du haut,
semble avoir une conséquence sur le débit à la sortie du bassin versant. C’est ainsi que
nous observons pour certains maxima de précipitation, un maxima sur l’hydrogramme. La
deuxième remarque s’attache à décrire un décalage entre ces deux maxima. En effet le
débit réagit plus tard sur le graphique car nous observons un maximum de débit atteint
après un maximum pluviométrique.

Cependant cette relation n’est pas présente sur la totalité des maxima pluviométriques et
certains pics de quantité d’eau tombée semblent n’avoir aucune conséquence sur le débit
à l’exutoire.

L’analyse ci-après permet de quantifier cette relation.

Pour exprimer la présence ou non d’une relation entre les deux variables hydrologiques
tracées, nous allons exprimer l’autocorrélation tout d’abord puis ensuite la corrélation entre
les deux variables.



                                                                                        23/27
Hydrologie statistique 2007-2008 M. Ababou

Le script ci-dessous nous a permis de représenter le graphique de la fonction
d’autocorrélation pour la variable hydrologique «précipitation», puis pour la fonction
d’autocorrélation de la variable «débit». Le dernier graphique représente la fonction de
covariance entre les deux variables hydrologiques.




clear all;

close all;

load -ascii 'QALIOU93_SEQ.txt';

load -ascii 'PALIOU93_SEQ.txt';

q= QALIOU93_SEQ; % Q pour l annee 93

p= PALIOU93_SEQ; % P pour l annee 93

maxlag = 350;

cpp=xcov(p,maxlag,'biased'); %'unbiased'

mcpp = length(cpp)

tau = [-maxlag:1:maxlag];

subplot(3,1,1)

plot(tau,cpp)

grid on

ymin = min(cpp);

ymax = max(cpp);

axis([-10 10 ymin ymax])

xlabel('Tau delai')

ylabel('Cpp')

title('Fonction d autocovariance')

cqq=xcov(q,maxlag,'biased'); %'unbiased'

subplot(3,1,2)

plot(tau,cqq)

                                                                                   24/27
Hydrologie statistique 2007-2008 M. Ababou

grid on

ymin = min(cqq);

ymax = max(cqq);

axis([-100 100 ymin ymax])

xlabel('Tau delai')

ylabel('Cqq')

title('Fonction d autocovariance')

hold on

cpq=xcov(q,p,maxlag,'biased'); %'unbiased'

subplot(3,1,3)

plot(tau,cpq)

grid on

ymin = min(cpq);

ymax = max(cpq);

axis([-100 200 ymin ymax])

xlabel('Tau delai')

ylabel('Cpq')

title('Fonction de covariance')

Ce script conduit aux trois graphiques ci-dessous :




                                                      25/27
Hydrologie statistique 2007-2008 M. Ababou


                                                 F o n c t io n d a u t o c o va ria n c e

              10
     Cpp




               5
               0
               -1 0     -8           -6    -4            -2          0           2         4      6          8    10
                                                               T a u d e la i
                                                 F o n c t io n d a u t o c o va ria n c e

               2
       C qq




               1
               0
              -1 0 0   -8 0      -6 0     -4 0          -2 0           0         20          40   60         80   100
                                                                 T a u d e la i
                                                     F o n c t io n d e c o va ria n c e
               3
               2
       Cpq




               1
               0
              -1 0 0          -5 0               0                   50                    100         150        200
                                                                T a u d e la i


      Figure 3 : Fonction d’autocorrélation pour la variable hydrologique
 «précipitation», puis pour la fonction d’autocorrélation de la variable «débit».
        Fonction de covariance entre les deux variables hydrologiques.

Sur le premier graphique, nous observons que la fonction d’autocovariance présente un
maximum centré sur 0. La valeur maximale atteinte par la fonction d’autocovariance est de
5 aux extrémités. Les valeurs sont toutes positives ce permet d’affirmer qu’il y a
uniquement corrélation positive. Aucun retard n’est observé.

Le second graphique présente la fonction d’autocorrélation de la variable hydrologique
«débit». De la même façon que pour les précipitations, nous observons que la fonction est
strictement positive. C’est qui permet d’annoncer une corrélation positive. Le maximum
atteint en 0 n’indique aucun retard.

Les résultats de ces deux graphiques sont cohérents puisqu’ils donnent des résultats
attendus à savoir que chaque variable hydrologique n’est corrélée avec elle-même
uniquement de façon positive. Ceci signifie que le mouvement de la série hydrologique est
identique pour chaque variable analysée, c'est-à-dire dans notre cas, deux fois le débit
«qq» et deux fois les précipitations «pp».

Après s’être penché sur ces deux premiers graphiques, nous pouvons maintenant

                                                                                                                        26/27
Hydrologie statistique 2007-2008 M. Ababou

analyser le graphique de la covariance entre les précipitations et le débit. Sur celui-ci,
nous observons que le premier pic de la fonction de covariance est atteint pour un
décalage compris entre 10 et 12 jours. A la suite de ce-dernier suit une première forte
décorrélation jusqu’à une soixantaine de jours, c'est-à-dire l’équivalent de 3 mois. Puis la
fonction de covariance croit de nouveau pour un second maxima moins fort que le
premier. A la suite de cette deuxième corrélation la fonction décroit fortement.

Si l’on s’intéresse aux délais négatifs, nous remarquons que la fonction devient négative
pour des délais de 60 jours. Ceci signifie qu’il y a corrélation négative ou encore anti-
corrélation. Ces valeurs négatives correspondent bien aux délais correspondant au
minimum de la fonction de covariance entre les deux pics.

De ces observations nous concluons qu’il existe donc une forte corrélation entre les deux
variables hydrologiques précipitations et débit dont celle-ci est maximale pour un délai de
10 à 12 jours. Ce délai étant perceptible pour des valeurs positives, nous en déduisons la
relation de causalité suivante : les précipitations ont une avance de 10 à 12 jours sur les
débits. C’est donc la variable précipitation qui est une cause du débit.

Ensuite, nous avons montré qu’il existe un durée pour laquelle les variables sont les moins
corrélées sur une période hydrologique cohérente (i.e approximativement 100 jours). Cette
période atteint une valeur de 60 jours soit 2 mois approximativement. Le deuxième de
corrélation aux alentours de 90 jours (du troisième mois) indique un regain de corrélation
entre ces deux variables. L’ensemble des explications peuvent être mise en relation avec
la nature karstique des variables hydrologiques analysées.

Le premier pic que l'on observe correspond certainement à la réaction rapide du bassin
versant karstique. L'eau passe essentiellement dans les fractures du bassin versant. Le
second pic correspond à la réponse lente du bassin. Ce pic permet d'estimer la réponse
poreuse du bassin versant. Les deux échelles de temps sont nettement perceptibles sur
ces deux pics. La réponse rapide avec le premier pic et la réaction du sous sol non
fracturé plus lente dans le milieu poreux.




                                                                The end...

                                                                                      27/27