Catastrophes naturelles et population mondiale

Pour reprendre une question posée au téléphone par François, et par Eric dans un commentaire (ici), il est possible de comparer l’évolution du nombre de catastrophes et de l’évolution de la population (mondiale). Pour les données de population on peut en trouver sur wikipedia (ici), sur le site des nations unies () après 1950, ou encore sur le site du census.gov (). Ensuite on peut moyenner les doublons.

POP=read.table("http://freakonometrics.blog.free.fr/public/data/pop-mondial.txt",
sep=";",header=TRUE)
population=tapply(POP$POP,POP$YEAR,mean)
years=as.numeric(names(population))
reg=lm(population~bs(years,3))
lines(years,predict(reg),col="red")
plot(years,population)
baseline=predict(reg,newdata=data.frame(years=1900:2010))

Revenons à nos catastrophes “naturelles“. Par exemple les tremblements de terre. La première piste (évoquée ) était de lisser cette série pour l’utiliser en baseline (multiplicatif). On retrouve ici une série qui oscille autour de 1.

earthquake=read.table("http://freakonometrics.blog.free.fr/public/data/DataSet.csv",
sep=",",header=TRUE)
flood=read.table("http://freakonometrics.blog.free.fr/public/data/DataSet_2_.csv",
sep=",",header=TRUE)
storm=read.table("http://freakonometrics.blog.free.fr/public/data/DataSet_3_.csv",
sep=",",header=TRUE)
library(splines)
n=nrow(earthquake)
earthquake=earthquake[-((n-3):n),]
DATE=1900:2010
EARTHQUAKE=FLOOD=STORM=rep(0,length(DATE))
EARTHQUAKE[as.numeric(as.character(earthquake$Year))-1899]=
as.numeric(as.character(earthquake$Number))
plot(DATE,EARTHQUAKE,type="b")
reg=lm(EARTHQUAKE~bs(DATE,6))
lines(DATE,predict(reg),col="red")
plot(DATE,EARTHQUAKE/predict(reg),type="b")
abline(h=1,lty=2,col="blue")

Si on utilise cette fois la population mondiale comme baseline, on obtient

POP=read.table("http://freakonometrics.blog.free.fr/public/data/pop-mondial.txt",
sep=";",header=TRUE)
population=tapply(POP$POP,POP$YEAR,mean)
years=as.numeric(names(population))
reg=lm(population~bs(years,3))
lines(years,predict(reg),col="red")
plot(years,population)
baseline=predict(reg,newdata=
     data.frame(years=1900:2010))/1000000
plot(DATE,EARTHQUAKE/baseline,type="b")

La série de points présente une jolie rupture. Afin de trouver la date optimale de rupture, utilisons un test basé sur le critère d’Akaike.

Y=EARTHQUAKE/baseline
X=DATE
DUMMY=1910:2000
VAIC=rep(NA,length(DUMMY))
for(x in DUMMY){
IX=(X>x)
reg=lm(Y~1+IX)
VAIC[x-1909]=AIC(reg)
}
plot(DUMMY,VAIC)
k=which.min(VAIC)
DUMMY[k]

Manifestement, il y a une rupture en 1975 (que je ne sais trop comment interpréter).

Si on regarde une autre série, comme celle des inondations.

n=nrow(flood)
flood=flood[-((n-3):n),]
DATE=1900:2010
FLOOD[as.numeric(as.character(flood$Year))-1899]=
as.numeric(as.character(flood$Number))
plot(DATE,FLOOD,type="b")
reg=lm(FLOOD~bs(DATE,6))
lines(DATE,predict(reg),col="red")
population=tapply(POP$POP,POP$YEAR,
mean)
years=as.numeric(names(population))
reg=lm(population~bs(years,3))
baseline=predict(reg,newdata=data.frame(years=1900:2010))
plot(DATE,FLOOD/baseline,type="b")

On retrouve ici clairement le fait que la croissance de la population ne peut expliquer la croissance du nombre d’inondations.