Tag Archives: ouragans

De la fréquence des ouragans aux Etats-Unis

Pour la fin du cours de statistique, on joue enfin avec des données. Après les taux de réussite au brevet des collèges, on a travaillé un peu sur la survenance des ouragans aux Etats-Unis. La base utilisée est la suivante

> base=read.table(
"http://freakonometrics.free.fr/ouragan.csv",
sep=";",header=TRUE,dec=",")
> base=base[1:207,1:7]

Le nombre annuel d’ouragan est la série suivante

>  TB <- table(base$Year)
>  years <- as.numeric(names(TB))
>  counts <- as.numeric(TB)
>  years0=(1900:2005)[which(!(1900:2005)%in%years)]
>  db <- data.frame(years=c(years,years0),
         counts=c(counts,rep(0,length(years0))))
>  X=db$counts
>  plot(db,type="h")

On peut comparer les fréquences observées avec ce que donnerait une loi de Poisson (paramètre obtenu par maximum de vraisemblance)

>  lambda=mean(X)
>  n=length(X)
>  cbind(X,dpois(0:6,lambda)*n)
  [,1]      [,2]
0   16 15.038430
1   34 29.367500
2   27 28.674870
3   13 18.665717
4    5  9.112744
5    6  3.559128
6    5  1.158396

On peut alors faire un test d’ajustement,

>  library(vcd)
>  fit=goodfit(X,"poisson")
>  plot(fit)

>  summary(fit)
 
	 Goodness-of-fit test for poisson distribution
 
                      X^2 df   P(> X^2)
Likelihood Ratio 14.17686  5 0.01452415

On voit que le test d’ajustement d’une loi Poisson n’est pas très concluant ici. Notons que ni la loi binomiale

>  fit=goodfit(X,"binomial")
>  plot(fit)
>  summary(fit)
 
	 Goodness-of-fit test for binomial distribution
 
                      X^2 df     P(> X^2)
Likelihood Ratio 46.21196  5 8.222959e-09

ni la loi binomiale négative ne semble marcher.

>  fit=goodfit(X,"nbinomial")
>  plot(fit)
>  summary(fit)
 
	 Goodness-of-fit test for nbinomial distribution
 
                      X^2 df   P(> X^2)
Likelihood Ratio 10.91507  4 0.02753523

On peut faire un test de stabilité temporelle du paramètre d’une loi de Poisson, à l’aide du test de Przyborowski & Wilenski (1940).  Un ancien package permettait de mettre en place ce test. Mais on peut heureusement utiliser une autre fonction, avec une version différente du test (on pourra aussi consulter Krishnamoorthy & Thomson (2004) pour aller plus loin). L’idée du premier test est assez simple: notons que si sur la première période (de durée ), on a un processus de Poisson de paramètre  et sur la seconde période (de durée ), alors

et

Un test de  s’écrit  Notons que la loi
conditionnelle de  sachant  suit une loi binomiale, donc le paramètre de probabilité est fonction de ce ratio. Plus précisément, conditionnellement à,  suit un loi binomiale  avec

A partir de là, il est facile de mettre en place une stratégie de test. Une autre stratégie est d’utiliser

>  seuil=1960
>  X1=db[db$years<seuil,"counts"]
>  X2=db[db$years>=seuil,"counts"]
>  N1=length(db$years<seuil)
>  N2=length(db$years>=seuil)
>  poisson.test(c(sum(X1),sum(X2)),c(N1,N2))
 
	Comparison of Poisson rates
 
data:  c(sum(X1), sum(X2)) time base: c(N1, N2)
count1 = 85, expected count1 = 103.5, p-value = 0.01216
alternative hypothesis: true rate ratio is not equal to 1
95 percent confidence interval:
 0.5218662 0.9265666
sample estimates:
rate ratio 
 0.6967213

 

Bref, on peut accepter que le nombre annuel d’ouragan n’est pas un processus de Poisson homogène dans le temps…

Autres Données pour le cours de Statistiques

Encore quelques données pour faire un TP de statistique. Ce sont des données de coûts d’ouragans, aux Etats-Unis,

> library(gdata)
> db=read.xls(
+     "http://sciencepolicy.colorado.edu/publications/special/public_data_may_2007.xls",
+     sheet=1)
> stupidcomma = function(x){
+  x=as.character(x)
+  for(i in 1:10){x=sub(",","",as.character(x))}
+  return(as.numeric(x))}
> base=db[,1:4]
> base$Base.Economic.Damage=
Vectorize(stupidcomma)(db$Base.Economic.Damage)
> base$Normalized.PL05=Vectorize(stupidcomma)(db$Normalized.PL05)
> base$Normalized.CL05=Vectorize(stupidcomma)(db$Normalized.CL05)

Lire la base est un peu pénible, à cause de virgules qui traînent comme séparateur de milliers. Mais peu importe, le petit code permet de nettoyer la base.

  • Fréquence d’ouragans

La première étape pourra être de travailler sur la fréquence. On pourra ajuster des lois de Poisson, ou binomiales, ou binomiales négatives, et faire des tests de stabilité du paramètre, avec le temps (on pourra  couper en 4 ou 5 blocs).

> TB <- table(base$Year)
> years <- as.numeric(names(TB))
> counts <- as.numeric(TB)
> years0=(1900:2005)[which(!(1900:2005)%in%years)]
> db <- data.frame(years=c(years,years0),
+       counts=c(counts,rep(0,length(years0))))
> db[88:93,]
   years counts
88  2003      3
89  2004      6
90  2005      6
91  1902      0
92  1905      0
93  1907      0
> plot(years,counts,type='h')

  • Coût des ouragans

Dans un second temps, on pourra travailler sur le coût des ouragans

> plot(base$Normalized.PL05/1e9,type="h")

On pourra tenter d’ajuster des lois, log normale, gamma, Pareto, etc. Faire des tests d’ajustement. Et tester la stabilité des paramètres, dans le temps (en coupant en 2 blocs, par exemple, ou davantage).

Arrêtez donc vos singeries !

Tout le monde connaît l’expérience où l’on compare les gains financiers obtenus par un traders expérimenté et un singe. Même si beaucoup de sites évoquent “la célèbre expérience“, je n ’en ai pas trouvé beaucoup – en France – qui donnent une référénce claire. Cette expérience  était née d’une théorie de Burton Malkiel dans son livre a random walk down wall street, qui affirmait “a blindfolded monkey throwing darts at a newspaper’s financial pages could select a portfolio that would do just as well as one carefully selected by experts“. En 1988, le Wall Street Journal a pris Burton Malkiel au mot, et a tenté l’expérience.Mais ils n’utilisèrent pas des singes, et ce furent des journalistes qui tiraient au hasard les actions sur lesquelles ils investiraient (alors que parallèlement, les mêmes sommes étaient confiées à des vrais inverstisseurs). Après 6 mois d’expérience, ils firent le bilan. “On October 7, 1998 the Journal presented the results of the 100th dartboard contest. So who won the most contests and by how much? The pros won 61 of the 100 contests versus the darts. That’s better than the 50% that would be expected in an efficient market. On the other hand, the pros losing 39% of the time to a bunch of darts certainly could be viewed as somewhat of an embarrassment for the pros. Additionally, the performance of the pros versus the Dow Jones Industrial Average was less impressive. The pros barely edged the DJIA by a margin of 51 to 49 contests. In other words, simply investing passively in the Dow, an investor would have beaten the picks of the pros in roughly half the contests (that is, without even considering transactions costs or taxes for taxable investors). The pro’s picks look more impressive when the actual returns of their stocks are compared with the dartboard and DJIA returns. The pros average gain was 10.8% versus 4.5% for the darts and 6.8% for the DJIA.” Il semble que l’Expansion en France ait aussi voulu tenter l’expérience.

La conclusion de cette expérience est qu’en moyenne, il est dur de battre le hasard ! Une autre conclusion peut aussi être qu’en mettant 500 singes dans une salle de marchés, au bout d’un certain temps on peut en trouver un qui sera aussi riche que Georges Soros, comme le disait Jean-Philippe Bouchaud. On peut parler de grandes déviations éventuellement…
Cette histoire peut prêter à sourire, effectivement, mais un think tank plutôt conservateur outre atlantique (le National Center for Public Policy Research) envisage de faire la même chose pour prédire les ouragans: au lieu d’écouter les climatologues, ils envisagent de faire appel à un singe (comme mentionné ici ou ).  Même si les prévisions des chercheurs ne sont pas toujours très fiables (en particulier sur les catastrophes comme les séismes ou les ouragants), elles reposent sur des éléments rationnels, avec souvent des intervalles de confiance…. Oui, je suis le premier à râler quand Météo France prévoit un beau week end et qu’au final il pleut, mais je suis agacé parces démagogues qui critiquent le travail des scientifiques (surtout qu’en l’occurence, je pense que le NCPPR confond  météo et climat, mais c’est un autre débat). Il y a un an, la communauté scientifique essayait de faire comprendre qu’elle n’était pas un parasite, voilà que maintenant il faut que l’on se batte pour faire comprendre qu’on fait des choses qu’un chimpanzé ne saurait pas faire ! On est décidément tombé bien bas…