This week, I discovered a picture on http://waitbutwhy.com/, which represent a (so-called) typical human life, in weeks,

I found that interesting. But the first problem is that I don’t understand the limit, below: 90 years, that’s not the average life length. That’s not what you should expect to live when you get born. The second problem is that it cannot be as static as it might seem, when you look at the picture. I mean, life expectancy at age 0 is not the same as life expectancy at age 30, or 50. So I did try to make an animated graph, using prospective life tables. Here a code to generate life tables, at different period, for a French population (I distinguish, here male and female)

library(demography)
france.fcast <- forecast(france.LC1,h=100)
L2 <- lifetable(france.fcast)
ex2=L2$ex L1=lifetable(fr.mort,series="female") ex1=L1$ex
exF=cbind(ex1,ex2)
france.fcast <- forecast(france.LC1,h=100)
L2 <- lifetable(france.fcast)
ex2=L2$ex L1=lifetable(fr.mort,series="male") ex1=L1$ex
exM=cbind(ex1,ex2)
Y=colnames(exF)

Based on those lifetables, we can extract remaining life expectancy, at various ages (say, for instance 50, 51, 52, etc), for someone born on some given year (say 1950). Based on those expected remaining lifetimes, we can plot

picture=function(yearborn=1950,age=50){
k=which(Y==yearborn)
M=diag(exM[,k+0:100])
F=diag(exF[,k+0:100])
par(mfrow=c(1,2))
va=0:(52*100-1)
plot(va%%52,va%/%52,cex=.6,pch=15,col=c("light yellow","light blue","white")[1+
(va>=age*52)*1+(va>(age+M[age+1])*52)*1],ylim=c(100,0),axes=FALSE,xlab="Week",
ylab="Age",main=paste("Man, born on ",yearborn,
", age ",age,sep=""))
axis(1)
axis(2)
plot(va%%52,va%/%52,cex=.6,pch=15,col=c("light yellow","pink","white")[1+
(va>=age*52)*1+(va>(age+F[age+1])*52)*1],ylim=c(100,0),axes=FALSE,xlab="Week",
ylab="Age",main=paste("Woman, born on ",yearborn,
", age ",age,sep=""))
axis(1)
axis(2)}

For instance, if we want the graph above, for someone age 30, born in 1980, we use

picture(1980,30)

Now, if we run a code to get an animated gif, we can get, for someone born in 1950,

and for someone born in 2000

Now, if I could get historical datasets, with the average time spent in schools, ages of retirement, etc, I guess I could add it on the graph. But that’s another story…

# Régression linéaire, quelques codes

Un rapide billet pour mettre en ligne les codes utilisés la semaine passée, complétant les codes des transparents. On travaille toujours sur la même base, ou on cherche à prévoir une distance de freinage d’un véhicule, tenant compte de la vitesse du véhicule.

> plot(cars)
> reg=lm(dist~speed,data=cars)
> summary(reg)

Call:
lm(formula = dist ~ speed, data = cars)

Residuals:
Min      1Q  Median      3Q     Max
-29.069  -9.525  -2.272   9.215  43.201

Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -17.5791     6.7584  -2.601   0.0123 *
speed         3.9324     0.4155   9.464 1.49e-12 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 15.38 on 48 degrees of freedom
Multiple R-squared:  0.6511,	Adjusted R-squared:  0.6438
F-statistic: 89.57 on 1 and 48 DF,  p-value: 1.49e-12

Pour faire plusieurs prévisions, à la main, on peut utiliser le code suivant (la boucle permet de faire des prévisions pour plusieurs valeurs)

> for(x in seq(3,30,by=.25)){
+ b0=coef(reg)[1]
+ b1=coef(reg)[2]
+ Yx=b0+b1*x
+ V=vcov(reg)
+ Vx=V[1,1]+2*V[1,2]*x+V[2,2]*x^2
+ IC1=Yx+c(-1,+1)*1.96*sqrt(Vx)
+ s=summary(reg)$sigma + IC2=Yx+c(-1,+1)*1.96*s + points(x,Yx,pch=19,col="red") + points(c(x,x),IC1,pch=3,col="blue") + points(c(x,x),IC2,pch=3,col="purple")} On avait ensuite fait une régression linéaire sur une sous-base, avec 20 observations tirées au hasard > I=sample(1:50,size=20) > reg=lm(dist~speed,data=cars[I,]) Le but était de visualiser l’impact du nombre d’observation sur la qualité de la régression > summary(reg) Call: lm(formula = dist ~ speed, data = cars[I, ]) Residuals: Min 1Q Median 3Q Max -23.529 -7.998 -5.394 11.634 39.348 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -20.7408 9.4639 -2.192 0.0418 * speed 4.2247 0.6129 6.893 1.91e-06 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 16.62 on 18 degrees of freedom Multiple R-squared: 0.7252, Adjusted R-squared: 0.71 F-statistic: 47.51 on 1 and 18 DF, p-value: 1.91e-06 > for(x in seq(3,30,by=.25)){ + b0=coef(reg)[1] + b1=coef(reg)[2] + Yx=b0+b1*x + V=vcov(reg) + Vx=V[1,1]+2*V[1,2]*x+V[2,2]*x^2 + IC=Yx+c(-1,+1)*1.96*sqrt(Vx) + points(x,Yx,pch=19,col="purple") + points(c(x,x),IC,pch=3,col="green")} Notons qu’il est possible d’utiliser des fonctions de R pour faire des prévisions, avec des intervalles de confiance > predict(reg, + newdata=data.frame(speed=c(15,25)),interval= "confidence") fit lwr upr 1 42.62976 34.75450 50.50502 2 84.87677 68.92746 100.82607 > predict(reg, + newdata=data.frame(speed=15),interval= "prediction") fit lwr upr 1 42.62976 6.836077 78.42344 Quand on a plus d’une variable explicative, c’est plus compliqué de “visualiser” la régression > chicago=read.table("http://freakonometrics.free.fr/chicago.txt", + header=TRUE,sep=";") > Y=chicago$Fire
>  X1=chicago$X_1 > X2=chicago$X_2
>  X3=chicago$X_3 > base=data.frame(Y,X1,X2,X3) > plot(X2,X3) > reg=lm(Y~X2+X3,data=base) > y=function(x2,x3) predict(reg,newdata=data.frame(X2=x2,X3=x3)) > VX2=seq(0,80) > VX3=seq(5,25) > VY=outer(VX2,VX3,y) > image(VX2,VX3,VY) > contour(VX2,VX3,VY,add=TRUE) qui correspond à un plan de régression > persp(VX2,VX3,VY,theta=30,ticktype=detailed) On reviendra plus en détails sur ce point, mais il est possible de faire des régressions non linéaires assez facilement, à partir de ce modèle linéaire. On avait commencé par un modèle linéaire sur le logarithme de la distance > plot(cars$speed,log(cars$dist)) > reg1=lm(log(dist)~speed,data=cars) > abline(reg1,col="red") (on le verra, ce n’est pas fini, car on n’a pas ici de prévision sur la distance, juste sur son logarithme… mais promis, on en reparlera) ou sur la racine carrée > plot(cars$speed,sqrt(cars$dist)) > reg1=lm(sqrt(dist)~speed,data=cars) > abline(reg1,col="red") Au lieu de transformer la variable d’intérêt, on peut aussi transformer la variable explicative. On peut pendre des puissances, ou des fonctions simples, mais aussi mettre des ruptures. On avait commencé par une variable indicatrice, > plot(cars$speed,cars$dist) > s=10 > abline(v=s,col="green") > regs=lm(dist~speed+I(speed>s),data=cars) > summary(regs) Call: lm(formula = dist ~ speed + I(speed > s), data = cars) Residuals: Min 1Q Median 3Q Max -29.472 -9.559 -2.088 7.456 44.412 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -17.2964 6.7709 -2.555 0.0139 * speed 4.3140 0.5762 7.487 1.5e-09 *** I(speed > s)TRUE -7.5116 7.8511 -0.957 0.3436 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 15.39 on 47 degrees of freedom Multiple R-squared: 0.6577, Adjusted R-squared: 0.6432 F-statistic: 45.16 on 2 and 47 DF, p-value: 1.141e-11 Mais on peut aussi mettre des fonctions afin d’avoir un modèle linéaire par morceaux, tout en étant continu > plot(cars) > s=15 > abline(v=s,col="green") > positive=function(x) ifelse(x>0,x,0) > regs=lm(dist~speed+positive(speed-s),data=cars) > summary(regs) Call: lm(formula = dist ~ speed + positive(speed - s), data = cars) Residuals: Min 1Q Median 3Q Max -29.502 -9.513 -2.413 5.195 45.391 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -7.6519 10.6254 -0.720 0.47500 speed 3.0186 0.8627 3.499 0.00103 ** positive(speed - s) 1.7562 1.4551 1.207 0.23350 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 15.31 on 47 degrees of freedom Multiple R-squared: 0.6616, Adjusted R-squared: 0.6472 F-statistic: 45.94 on 2 and 47 DF, p-value: 8.761e-12 On a ici une rupture, mais on pourrait imaginer en avoir plusieurs > nouvellebase=data.frame(speed=5:25) > y=predict(regs,newdata=nouvellebase) > lines(5:25,y,col="red") > > plot(cars$speed,cars$dist) > s1=10 > s2=20 > abline(v=c(s1,s2),col="green") > positive=function(x) ifelse(x>0,x,0) > regs=lm(dist~speed+positive(speed-s1)+positive(speed-s2),data=cars) > summary(regs) Call: lm(formula = dist ~ speed + positive(speed - s1) + positive(speed - s2), data = cars) Residuals: Min 1Q Median 3Q Max -24.374 -9.475 -2.625 6.639 43.914 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) -7.6305 16.2941 -0.468 0.6418 speed 3.0630 1.8238 1.679 0.0998 . positive(speed - s1) 0.2087 2.2453 0.093 0.9263 positive(speed - s2) 4.2812 2.2843 1.874 0.0673 . --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 15 on 46 degrees of freedom Multiple R-squared: 0.6821, Adjusted R-squared: 0.6613 F-statistic: 32.89 on 3 and 46 DF, p-value: 1.643e-11 Comme vu en cours, le test de significativité des deux derniers coefficients ne veut pas dire que la pente est nulle, mais qu’elle est significativement différente de cette obtenue sur la zone de gauche (avant les deux seuils). # Toronto, Statistical Conference I will be in Toronto this week, for the SSC annual conference (Canadian Statistical Society), to present some recent work, with Emmanuel Flachaire. Because of some administrative duty last week, I am a bit late, so I won’t be able to upload the slides before the talk. Sorry about that… More to come when I’ll be back in Montréal ! # De la difficulté de faire des prévisions (quand on a peu de données) Depuis plusieurs mois, on observe un engouement (probablement légitime) pour le big data. Si beaucoup peut être fait pour utiliser les volumes énormes de données à la disposition des assureurs, il convient de garder en mémoire que dans beaucoup de cas, les données sont rares et que la technologie ne devrait pas pouvoir changer grand chose. Le manque de données (fiables) crée une variabilité importante. • Loi des grands nombres, approximations et statistique asymptotique Tout comme les assureurs, les statisticiens aiment les gros volumes de données. Pour les assureurs, les gros portefeuilles sont considérés comme moins incertains, et pour les mêmes raisons, les statisticiens disposent, avec des grosses bases de données, d’estimateurs moins volatiles (on parle d’ailleurs du risque associé à un estimateur). Supposons que l’on suive un portefeuille de 1047 assurés, pendant un an, et que 159 ait réclamé un sinistre auprès de la compagnie. La probabilité qu’un assuré ait un sinistre est $p=\frac{159}{1047}\sim 15.2\%$ Si on admet que l’ont dispose d’assez d’observations, on peut utiliser une approximation Gaussienne du taux de réclamation, et obtenir un intervalle de confiance à 95%, de telle sorte qu’il y a 95 chances sur 100 que le taux de réclamation appartienne à l’intervalle1 $\left[p\pm\sqrt{\frac{p(1-p)}{n}}\right]\sim [13\%;17\%]$ On utilise ici l’approximation2 de la loi binomiale (le nombre de personnes qui réclame suit une loi binomiale) par une loi Gaussienne, comme sur la Figure ci-dessous, Avec cette approximation, on peut valoriser des contrats de réassurance. On peut ainsi approcher les probabilités d’avoir des dérives de la sinistralité dans notre portefeuille. Dans cet exemple, on notera que la probabilité d’avoir un taux de réclamation excédant 20% est de l’ordre de 0.0014%. • Prédire des “cygnes noirs” Supposons qu’un assureur décide de vendre à des collectivités locales des contrats d’assurance en cas de marée noire. Dans un modèle actuariel standard, la prime à payer devrait être $p\cdot c\cdot i$, où $p$ est la probabilité qu’une marée noire survienne par kilomètre de côte, $c$ la longueur de la côte de la région assurée, et $i$ le montant de l’indemnité versée. Après 5 ans d’expérience, l’assureur souhaite faire le bilan de son produit. Sauf qu’aucune marée noire n’a été observée. Que peut-il dire sur son tarif ? Cet exemple (à peine) fictif est a rapprocher d’une question posée dans les années 1950 à l’actuaire L.H. Longley-Cook: est-il possible de prédire (ou d’estimer) la probabilité d’observer, une année donnée, une collision en plein vol, entre deux avions ? Il n’y avait jamais eu de (grave) collision d’avions commerciaux lorsque la question avait été posée. Et sans aucune expérience passée, les statisticiens ne savaient trop quoi répondre. Pourtant, Longley-Cook avait prédit ” anything from 0 to 4 […] collisions over the next ten years“, comme le raconte McGrayne (2012). Deux ans plus tard, 128 personnes perdirent la vie au dessus de Grand Canyon lors d’une telle collision, et quatre an plus tard, ce sont 133 personnes qui périrent au dessus de New York. • La réponse Bayésienne L’idée de Longley-Cook était d’utiliser des techniques bayésienne en actuariat. Comme le notait Liu et al. (1996)Statistical methods with a Bayesian flavor […] have long been used in the insurance industry”. Dans Charpentier (2007), nous avions vu l’apport des techniques de crédibilité et les liens avec la statistique bayésienne. Et comme le raconte McGrayne (2012), Arthur Bailey a joué un rôle essentiel, précisément en proposant des techniques dans le cas d’évènements (très) rares. [Arthur] Bailey spent his first year in New York [in 1918] trying to prove to himself that ‘all of the fancy actuarial [Bayesian] procedures of the casualty business were mathematically unsound.’ After a year of intense mental struggle, however, he realized to his consternation that actuarial sledgehammering worked. He even preferred it to the elegance of frequentism. He positively liked formulae that described ‘actual data . . . I realized that the hard-shelled underwriters were recognizing certain facts of life neglected by the statistical theorists.’ He wanted to give more weight to a large volume of data than to the frequentists small sample; doing so felt surprisingly ’logical and reasonable.’ He concluded that only a ‘suicidal’ actuary would use Fishers method of maximum likelihood, which assigned a zero probability to nonevents. Since many businesses file no insurance claims at all, Fishers method would produce premiums too low to cover future losses. Pour rappel, le théorème de Bayes permet d’écrire $\mathbb{P}(A\vert B)=\frac{\mathbb{P}(A\cap B)}{\mathbb{P}(B)}=\frac{\mathbb{P}(A)}{\mathbb{P}(B)}\mathbb{P}(B\vert A)\propto \mathbb{P}(A) \mathbb{P}(B\vert A)$ Si on suppose disposer d’un échantillon d’observations $\{x_1,\cdots,x_n\}$ de loi $F_\theta$, où, classiquement $\theta$ est vu comme un paramètre inconnu que l’on cherche à estimer. Ici, on va supposer qu’il s’agit d’une variable aléatoire dont on va pouvoir obtenir une distribution (dite à posteriori) compte tenu des observations dont on dispose. On écrit alors $\mathbb{P}(\Theta=\theta\vert \underline{\boldsymbol{x}}_n)\propto \underbrace{\mathbb{P}(\Theta=\theta)}_{\text{loi a priori}}\cdot \underbrace{\mathbb{P}(\underline{\boldsymbol{x}}_n\vert \Theta=\theta)}_{\text{vraisemblance}}$ où le premier terme est une loi a priori que l’on se donne sur le paramètre , et le second est la vraisemblance de l’échantillon. Afin d’illustrer un peu, reprenons notre exemple initial, où on dispose de 5 assurés, et qu’aucun n’ai réclamé de sinistre. Que peut-on dire sur la probabilité qu’un assuré ait un sinistre. Avec les approches classiques, la probabilité estimée serait nulle. Ce qui est gênant pour calculer ensuite une prime convenable. La réponse Bayésienne consiste à suppose que la probabilité de réclamer un sinistre est une variable aléatoire, dont on peut calculer la loi a posteriori, compte tenu du fait que l’on observé – pour l’instant – cinq observations nulle. En utilisant la relation précédente, la densité a posteriori est ici $\pi(\theta\vert \underline{\boldsymbol{x}}_n)\propto \underbrace{\pi(\theta)}_{\text{loi a priori}}\cdot \underbrace{\mathbb{P}(\underline{\boldsymbol{x}}_n\vert\theta)}_{\text{vraisemblance}}$ où le second terme est la probabilité d’avoir cinq fois 0 avec une loi de Bernoulli, ce qui vaut $(1-\theta)^5$, et où le premier terme est une loi que l’on se donne, a priori pour $\Theta$. Une loi qui pourrait sembler naturelle serait la loi uniforme, pouvant être perçu comme non informative. Dans ce cas, la densité de $\Theta$, a posteriori, est proportionnelle à $(1-\theta)^5$. Comme le montre la Figure ci-dessous, on a ici beaucoup de chance que $\Theta$ prenne une valeur proche de 0 (ce qui a du sens, car on n’a jamais observé de réclamation, jusqu’à présent). A priori, il a avait 95 chances sur 100 que la probabilité soit inférieure à 95% (c’est l’idée de la loi uniforme comme loi a priori). Après cinq expériences sans sinistre, a posteriori, il y a 95 chances sur 100 que la probabilité de réclamer soit inférieure à 40%. Si nous avions observé dix fois la valeur nulle, il y aurait, a posteriori, 95 chances sur 100 que la probabilité de ne pas réclamer soit inférieure à 25%. La Figure ci-dessous montre la déformation de la loi a posteriori, au fur et à mesure de la mise à jour, avec la loi a priori à gauche, puis a posteriori sachant qu’une observation est arrivée (nulle), puis une seconde, puis une troisième, etc. Sur la figure ci-dessous, il s’agit à nouveau de la déformation de la loi a posteriori, avec une loi a priori différente (en l’occurrence, la loi dite non-informative de Jeffreys) • Les outils statistiques pour le small data De même que le big data pose des problèmes aux techniques statistiques classiques, et permet de découvrir les techniques d’apprentissage, le small data pose des problèmes à laquelle la statistique classique n’apporte aucune réponse. Les actuaires redécouvrent les techniques bayésiennes, mais aussi d’autres outils qui permettent malgré tout d’apporter des réponses quantitatives pertinentes. Pour quantifier l’incertitude lors du calcul des provisions pour sinistres à payer, les actuaires utilisent des techniques de rééchantillonage (bootstrap) depuis une vingtaine d’années (rappelé dans Charpentier, Devineau & Nessi (2010)) pour tenir compte du peu de données dans les triangles de liquidation. Si le big data ouvre des perspectives fascinantes à explorer pour les actuaires (en particulier en souscription), il convient de garder en mémoire que dans beaucoup de cas, les actuaires doivent quantifier des risques avec très peu de données. L’avenir n’est aucunement au big data ! L’avenir est au contraire dans la connaissance de techniques permettant de traiter à la fois des gros volumes de données, mais aussi de faire face à un insuffisance de données. 1 Dans le cas des événements (très) rares, l’approximation par une loi de Poisson sera nettement meilleure, mais ici, la probabilité est de l’ordre de 15%, ce qui est trop élevé. 2 On retrouve ici les intervalles de confiance mentionnés désormais dans les sondages. # Allez les Bleus ! In almost three weeks, the (FIFA) World Cup will start, in Brazil. I have to admit that I am not a big fan of soccer, so I will not talk to much about it. Actually, I wanted to talk about colors, and variations on some colors. For instance, there are a lot of blues. In order to visualize standard blues, let us consider the following figure, inspired by the well known chart of R colors, BLUES=colors()[grep("blue",colors())] RGBblues=col2rgb(BLUES) library(grDevices) HSVblues=rgb2hsv( RGBblues[1,], RGBblues[2,], RGBblues[3,]) HueOrderBlue=order( HSVblues[1,], HSVblues[2,], HSVblues[3,] ) SetTextContrastColor=function(color) ifelse( mean(col2rgb(color)) > 127, "black", "white") TextContrastColor=unlist( lapply(BLUES, SetTextContrastColor) ) c=11 l=6 plot(0, type="n", ylab="", xlab="",axes=FALSE, ylim=c(0,11), xlim=c(0,6)) for (j in 1:11){ for (i in 1:6){ k=(j-1)*6 + i rect(i-1,j-1,i,j, border=NA, col=BLUES[ HueOrderBlue[k] ]) text(i-.5,j-.5,paste(BLUES[k]), cex=0.75, col=TextContrastColor[ HueOrderBlue[k] ])}} All the color names that contain “blue” in it are here. Having the choice between several possible colors is interesting, but it can also be interesting to get a palette of blue colors, What we can get is the following library(RColorBrewer) blues=colorRampPalette(brewer.pal(9,"Blues"))(100) In order to illustrate the use of palette colors, consider some data, on soccer players (officially registered). The dataset – lic-2012-v1.csv – can be downloaded from http://data.gouv.fr/fr/dataset/… (I will also use a dataset we have on location of all towns, in France, with latitudes and longitudes) base1=read.csv( "http://freakonometrics.free.fr/popfr19752010.csv", header=TRUE) base1$cp=base1$dep*1000+base1$com
base2=base2[base2$fed_2012==111,] names(base2)[1]="cp" base2$cp=as.numeric(as.character(base2$cp)) The problem with France (I should probably say one of the many problems) is that regions and departements are not well coded, in the standard functions. To explain where départements are, let us use the dept.rda file, and then, we can get a matching between R names, and standard (administrative) ones, base21=base2[,c("cp","l_2012","pop_2010")] base21$dpt=trunc(base21$cp/1000) library(maps) load("dept.rda") base21$nomdpt=dept$dept[match(as.numeric(base21$dpt),dept$CP)] L=aggregate(base21$l_2012,by=list(Category=base21$nomdpt),FUN=sum) P=aggregate(base21$pop_2010,by=list(Category=base21$nomdpt),FUN=sum) base=data.frame(D=P$Category,Y=L$x/P$x,C=trunc(L$x/P$x/.0006))
france=map(database="france")
matche=match.map(france,base$D,exact=TRUE) map(database="france", fill=TRUE,col=blues[base$C[matche]],resolution=0)

Here are the rates of soccer players (with respect to the total population) It is also possible to look at rate not by département, but by town,

base10=base1[,c("cp","long","lat","pop_2010")]
base20=base2[,c("cp","l_2012")]
base=merge(base10,base20)
Y=base$l_2012/base$pop_2010
QY=as.numeric(cut(Y,c(0,quantile(Y,(1:99)/100),10),labels=1:100))
library(maps)
map("france",xlim=c(-1,1),ylim=c(46,48))
points(base$long,base$lat,cex=.4,pch=19,col=blues[QY])

The darker the dot, the more player, We can also zoom in, to get a better understanding, in the northern part of France, for instance, or in the Southern part,

We can obtain a map which is not (too) far away from the one mentioned a few months ago on http://slate.fr/france/78502/.

# Modèles de Prévision, suite

La seconde partie de transparents pour le cours ACT6420 est en ligne. Il s’agit d’un document qui permettra à tout le monde garder des traces, mais je n’entends pas forcément les suivre… Par exemple, ce soir, je propose de jouer un peu avec les données d’incendies à Chicago

> chicago=read.table("http://freakonometrics.free.fr/chicago.txt",
> Y=chicago$Fire > X1=chicago$X_1
> X2=chicago$X_2 > X3=chicago$X_3
> base=data.frame(Y,X1,X2,X3)

ou une base sur les états américains (cf un ancien billet)

> US=read.table("http://freakonometrics.free.fr/US.txt", + header=TRUE,sep=";")

Sinon, la semaine passée, on avait évoqué la matrice chapeau (dont on reparlera davantage lorsque l’on évoquera la phase de diagnostique). L’idée était de noter que

$\widehat{\boldsymbol{Y}} = H \boldsymbol{Y}$

avec

$H = \boldsymbol{X} \left(\boldsymbol{X}^\top \boldsymbol{X} \right)^{-1} \boldsymbol{X}^\top$

On peut visualiser cette matrice avec le code suivant

> X=cbind(1,cars\$speed)
> M=X %*% solve( t(X) %*% X) %*% t(X)
> Q=matrix(as.numeric(cut(M,c(-10,quantile(M,(1:5)/6),10),labels=1:6)),50,50)
> plot(1:50,1:50,col="white",axes=FALSE,xlab="",ylab="",ylim=c(50,0))
> library(RColorBrewer)
> couleur=rev(brewer.pal(6,"RdBu"))
> for(i in 1:50){for(j in 1:50){ points(i,j,col=couleur[Q[i,j]],pch=15,cex=1)}}
> axis(1)
> axis(2)

avec des cases rouges pour les valeurs importantes (et positive) et des case bleues pour des valeurs importantes (mais négatives)

# Πάντες ἄνθρωποι τοῦ εἰδέναι ὀρέγονται φύσει

J’essaye continuellement de me convaincre qu’Aristote avait raison lorsqu’il affirmait que tous les hommes désiraient naturellement savoir. C’est cette motivation qui fait que j’aime toujours enseigner, et partager des choses que j’ai pu apprendre, sur mon blog. Ces dernières semaines, j’ai été amené à plusieurs reprises à repenser mon activité de blogger, après qu’un ami qui suivait la blogosphère depuis plusieurs années me disait qu’il allait se lancer (sur son propre blog) et m’interrogeait sur ma pratique, qu’un autre qui bloguait sous un pseudo allait tenter de bloguer aussi sous son vrai nom, ou après quelques discussions sur Twitter sur une tentative de définition du terme “blogs scientifiques“. Ayant promis de revenir un peu sur ma conception des blogs académiques et scientifiques, j’ai lu tout ce qui me passait sous la main, afin de réfléchir encore un peu sur ma pratique, et sur ma position vis à vis de la vulgarisation. Comme toujours, c’est un peu brouillon, mais c’est le principe de mon blog.

Lorsque j’ai migré mon blog sur la plateforme d’hypotheses.org, j’aimais l’idée d’avoir un carnet de recherche ouvert. Car mon blog n’est rien d’autre que cela. Et mon blog est à l’image de mes cahiers, dans lequel je note des idées que je peux avoir dans le métro en allant chercher ma fille chez sa gardienne, des preuves que je dois faire en cours, des idées de visualisation glanées dans des exposés, ou dans des articles, etc. Il n’a aucune autre prétention. Ou au contraire, ce blog a une prétention incroyable: celui de redonner la possibilité, à tous ceux qui le souhaitent, de s’approprier des outils et des techniques d’ordinaire laissés aux mains des “experts”. Comme disait Gustave Flaubert, dans Bouvard et Pécuchet, “la science est un monopole aux mains des riches. Elle exclut le peuple…“. J’espère que ce blog aidera a réparer cela.

• Définir le terme “vulgarisation

Comme le notait Marcel Thouin, “selon la définition usuelle, qui est peut-être illusoire, la vulgarisation scientifique consiste à rendre la connaissance scientifique accessible“. Il ne s’agit pourtant pas d’open-access, mais plus précisément, selon  Marcel Thouin (rappelant la définition du dictionnaire) la vulgarisation consiste à “adapter un ensemble de connaissances scientifiques et techniques de manière à les rendre accessibles à un lecteur non spécialiste” (on retrouvera des propos similaire chez Henri Canonge). Mais le terme “lecteur non spécialiste” est vague. Pablo Jensen parle lui du “citoyen” mais sans forcément poser de définition claire des termes (pour sa défense, l’article mentionné est un article de vulgarisation). Patrice Cohen-Seat (cité par Roger Clausse et repris dans la thèse de Baudouin Jurdant) disait que vulgariser, c’est “transposer à l’usage de tous ce qui a été conçu et élaboré dans le langage de quelques uns“. Cela dit, Maurice Mashaal par exemple donne un sens plus large au terme “vulgariser“, en autorisant “tous” à être assez restreint. Selon lui, “les chercheurs ne vulgarisent pas uniquement en direction du grand public. Tout d’abord, ils sont généralement curieux et aiment savoir ce qui se passe dans d’autres spécialités que la leur“. On voit alors une catégorie un peu distincte, entre le “grand public” et les lecteurs spécialistes, éventuellement dans d’autres spécialités que la personne qui transmet le message. Cette classification est retenue très clairement dans le rapport Kunth, qui retient trois formes de communication scientifique:

• la diffusion de l’information scientifique spécialisée (qui se fait dans les conférences spécialisées, les séminaires, et le publication académique pointues, traditionnelle)
• la diffusion scientifique interdisciplinaire
• la vulgarisation scientifique (diffusion de la culture scientifique et technique, s’adressant au grand public ou aux jeunes).

Avec cette classification, je ne classerais pas mon blog dans la catégorie des blogs de vulgarisation. J’aime l’idée de cette catégorie un peu intermédiaire. Cela dit, dans Communicating Science, Professional, Popular, Literary, Nicholas Russell rappelle que Nature et Science sont depuis toujours des revues de vulgarisation. On peut donc vulgariser tout en visant un public qui n’est pas forcément “tout le monde“, et qui autorise l’utilisation d’un jargon technique. Autrement dit, la seconde catégorie serait de la vulgarisation. Il va donc falloir creuser un peu plus pour comprendre si j’ai un blog de vulgarisation, ou pas.

• Un espace entre la recherche et l’enseignement

J’ai toujours revendiqué le fait d’être à la fois un chercheur et un enseignant, un “savant” et un “pédagogue” pour reprendre les termes de Dider Moreau (justement pour caractériser la vulgarisation). Et mon blog est précisément un espace entre la recherche et l’enseignement. Mais ce n’est ni l’un, ni l’autre. Et quelque part, l’activité de vulgarisation se situe aussi dans cet espace. Pour reprendre un point noté par Baudouin Jurdant, “la relation pédagogique instaurée par la vulgarisation se distingue en outre de celle en jeu dans l’enseignement par le fait qu’elle se développe dans un espace social ouvert, non marqué, abstrait. En bénéficiant de lieux institutionnellement marqués, clos par rapport à l’espace de la chose publique, et concrets, l’enseignement assure une certaine signifiance à l’économie de sa relation pédagogique. La vulgarisation au contraire, est un propos public. Son espace est celui de la société toute entière, abstrait puisqu’il est ouvert à tous et par là, impénétrable“. Daniel Jacobi, Bernard Schiele et Marie-France Cyr parlent d’ailleurs d'”éducation non-formelle” pour définir l’activité de vulgarisation. Mais le blog s’inscrit aussi dans cet dimension, avec un espace ouvert, et un discours plus informel.

Mais il ne s’agit pas non plus de recherche. Comme l’écrivaient Renaud Dulong et Werner Ackermann, “le scientifique et le vulgarisateur ne répondent pas en fait aux mêmes interrogations et le discours de vulgarisation ne peut donc être uniquement (que) le reflet appauvri des concepts élaborés par le scientifique.” Mais mon blog ne cherche pas à simplifier. Il n’y a aucun enrobage, au contraire. Il s’agit d’enlever l’enrobage des articles académiques, en expliquant peut être davantage des pistes qui n’aboutissent pas, des conjectures, en n’essayant pas d’impressionner la galerie. Une autre différence est que le langage utilisé dans un blog est beaucoup plus vulgaire que celui utilisé dans un article scientifique. Cela dit, comme le rappelle Florine Vital-Durant, en 1632, lorsque Galilée publie ses dialogues (dialogo sopra i due massimi sistemi del mondo), il les publie en italien, et pas en latin. L’idée d’écrire de la science dans un langage plus ordinaire, plus quotidien, sans pour autant faire de la vulgarisation ne date pas d’hier !

• Quelle différence entre un blog académique et la vulgarisation ?

Même si certains me reprochent l’utilisation du terme blog académique, je continue à revendiquer ce terme. Car mon statut d’académique, ou d’enseignant et de chercheur me donne une légitimité que n’a pas un vulgarisateur. D’ailleurs, dans carrière scientifique, morale scientifique et vulgarisation, Luc Boltanski et Pascale Maldidier remettraient presque en cause l’idée de journalistes faisant de la vulgarisation, lorsqu’ils notent que “la vulgarisation, devrait perdre toute légitimité lorsqu’elle est l’activité spécifique d’un groupe totalement extérieur à la cité savante“. Et on retrouve cette idée chez Daniel Jacobi, Bernard Schiele et Marie-France Cyr lorsqu’ils parlent des “usurpateurs“. En fait, Luc Boltanski et Pascale Maldidier vont même plus loin, lorsqu’ils affirment que “le monde des vulgarisateurs ne peut pas être considéré comme un univers autonome possédant un principe de légitimité propre et qu’il ne peut être saisi indépendamment du monde savant défini et dominé par l’institution universitaire. Il s’ensuit que les attitudes des vulgarisateurs et que les représentations qu’ils ont de leur propre pratique sont déterminées par la position marginale qu’ils occupent dans ce champ hiérarchisé et qu’elles se définissent essentiellement par rapport à la définition sociale de la vulgarisation comme activité illégitime qu’engendré, selon les mécanismes précédemment décrits, la communauté scientifique“. Mon blog est un blog de science tenu par un scientifique reconnu par ses pairs, et on retrouve ici la distinction que faisait Tom Roud sur son blog, commentant l’annonce de l’arrivée de blogs scientifiques hébergés par Le Monde.

Bref, je continue de croire que mon blog n’est pas un blog de vulgarisation, et que les blogs académiques, n’en déplaisent à certains, ont une réelle légitimité. Je revendique l’idée d’avoir un carnet de recherche ouvert à tous, comme le souhaitait hypotheses.org. Si un jour je fais un billet sur les modèles causaux, ce ne sera  pas par soucis de vulgariser un sujet pour un lectorat (qui me semble infiniment lointain), mais juste pour moi. Au lieu de noter dans un cahier, je posterais sur mon blog. Avec des illustrations, des codes, des données (bref, des choses qui tiennent mal dans un cahier peut être).

• Rendre public, ou viser un public ?

Bloguer permet d’être public (et éventuellement, d’avoir un lectorat régulier, ou pas) sans jamais rencontrer son public. Et pour moi qui suis un vieil ours, c’est idéal. Contrairement aux vulgarisateurs qui ont besoin d’un public, mon blog n’en a pas besoin. Lorsque je l’ai lancé, il existait avec 4 lecteurs, et il continuera d’exister même si la moitié ne me lit plus ! Peu importe, ce n’est pas sa finalité. Cela dit, la semaine passée, j’étais en Europe, j’ai fait plusieurs exposés, et plusieurs fois, j’ai entendu le “je suis un fidèle lecteur de votre blog”. Et c’est toujours une réflexion qui me surprend. J’ai toujours l’impression d’être un obscur professeur de statistique, ou d’économétrie, qui emploie tellement de jargon que personne ne devrait suivre (c’est en tous cas assez contradictoire avec certains propos tenus par des journalistes scientifiques, sur la vulgarisation). D’ailleurs, mardi dernier, j’ai vécu une aventure amusante. J’ai été contacté par un étudiant qui cherchait des données pour construire un modèle (en l’occurrence des informations climatiques, par département – en France – et sur plusieurs années). J’ai posé la question sur Twitter. En quelques secondes, le message a été relayé intensivement, et au bout de quelques minutes, j’avais des suggestions de trois sites. Comme me le disait ma femme, ce n’est pas la magie de Twitter, c’est que j’ai un réseau incroyable de followers actifs. Au delà des discours nauséabonds des CM qui expliquent que 6000 followers c’est fabuleux pour se faire connaître et assurer une notoriété en ligne, avoir 6000 followers, c’est surtout génial pour trouver de l’aide sur à peu près n’importe quel sujet. Au risque de paraître schizophrène, ami lecteur, je dois confesser que je ne blogue pas pour toi. Mais je te suis reconnaissant d’avoir accepté de m’intégrer dans une communauté aussi active, disponible et curieuse. Merci.

# 24 heures (de science)

Cette fin de semaine, ce sont les 24 heures de sciences, au Québec

Il paraît que c’est une fête. Coïncidence ou pas, j’ai fait mes 24 heures de science hier ! Après trois jours assez intense en début de semaine, en Belgique puis à Paris, avec un jury de thèse, un mini cours, un exposé, et plusieurs rendez-vous (pour des partenariats universitaires, avoir accès à des données, monter des formations, etc) je suis finalement rentré hier, jeudi. Compte tenu du prix des hôtels sur Paris, j’essaye autant que possible de squatter chez des collègues lorsque je dois y passer une nuit. Donc mercredi soir, j’étais chez un ami (et collègue), et on a parlé travail jusqu’à tard dans la nuit… A six heures, réveil matinal (enfants obligent) pour rejoindre l’aéroport. Comme le 8 mai est férié, j’avais essayé de prévoir un peu de temps car il y avait moins de transports en commun. Arriver trois heures en avance, c’est pénible, mais d’un autre côté, ça m’a laissé du temps pour avancer la rédaction d’un rapport d’arbitre que j’ai du mal à finir (faute de temps). Dans l’avion, j’ai pu me remettre à des papiers qu’il va falloir finir dans les semaines à venir (maintenant que la session d’hiver est finie, la recherche pourra reprendre !). J’ai pu arriver à l’heure, après presque 7 heures de travail non-stop (juste le temps de manger un plateau repas). A 16 heures (heure locale, soit 22 heures en France), j’arrivais enfin au bureau. Oui, les jeudis soirs, je donne un cours pour la session d’été. Bref, il fallait réviser un peu (relever aussi son courrier, récupérer les reçus pour les remboursements, etc), et j’ai pu donner mon cours, entre 17:30 et 21 heures. Le cours finissait à 21 heures (3 heures du matin, en France), et quand je suis rentré à la maison, les enfants dormaient depuis longtemps. J’ai alors du répondre aux courriels qui s’étaient accumulés depuis le matin, en suivant distraitement le match… et j’ai fini par m’endormir un peu avant 23 heures (5 heures du matin en France). Je n’ai pas pu tenir 24 heures d’affilée… je commence décidément à vieillir. Et je ne sais pas s’il y a des gens que ça amuse de faire de la science pendant 24 heures, mais personnellement, je n’ai qu’une envie, c’est d’aller dormir !

# Computational Actuarial Science

After some delay, the book Computational Actuarial Science with R is now annonced for July 2014. I don’t know if we will be able to get copies for the R in Insurance conference, in London, but I guess everyone is working on it. And kindly, CRC sent me the following flyer, with some reduction code to save 20% (when ordering on CRC’s website).

# Big Data, à Paris

Mercredi maitin, après une dernière réunion à Bruxelles, je prends le train pour être à Paris à l’heure du repas du midi pour faire une (rapide) intervention sur le big data, pour le groupe de travail big data de l’Institut des Actuaires. Je crois que l’intervention est prévue à 12:30, chez Optimind, rue de la Boëtie. J’ai essayé de mettre quelques éléments de réflexion dans des transparents, histoire de lancer des débats !

Je passerais la soirée à Paris (après quelques réunions l’après-midi), avant de redécoller jeudi matin.

# Graduate Crash Course on Risk Measures

Tomorrow morning, I will give a crash course on risk measures at Louvain-la-Neuve, in Belgium. This is a crash course of PhD students (and researchers) with a long introduction on the univariate static framework (and some mathematical tools that will be interesting later on, such as the Fenchel transform and more generally on convexity, as well as some results on optimal transport). I will also mention what was obtained in decision theory, inspired by Itzhak Gilboa‘s Theory of Decision under Uncertainty. Then I will mention extensions to derive multiple risk measures, based on Marc Henry and Alfred Galichon‘s work. Finally, I will conclude by introducing the difficulty to derive dynamic risk measures.

The slides are based on a document I am still working on. And unfortunately, the deeper I get to explain the roots of the axioms, or the assumptions, the more papers I discover (and I need to read, and understand). So I guess I will need some time to finalize my survey. Note that I decided to skip details on technical issues when working on $L^\infty$, and the weak topology on the dual of $L^\infty$. I will try to add additional references in the notes, but I wanted the slides to be as simple as possible. I also want to add more connections with statistical results, such as Neyman Pearson’s lemma, for instance (as mentioned in a paper by Alexander Schied). All my apologies for the typos, too.

# PhD defense, in Louvain

I was, this morning, at the (private) PhD defense of Mathieu Pigeon, at UCL, Université Cathololique de Louvain. Mathieu obtained some very interesting results on claims reserving based on individual observations. I will try to upload the slides of his defense, someday.

# Short (but Buzy) Trip in Europe

I will flight tomorrow evening to go to Belgium. I will attend a PhD defense this Monday, at Louvain-la-Neuve. On Tuesday, I will give a crash course on risk measures, for PhD Students, in Louvain-la-Neuve. Then, Wednesday morning, I have a meeting with some friends in an insurance company, in Brussels, and I will be around noon in Paris, to give a talk on Big Data, for the Big Data Task Force of the (French) Institute of Actuaries. I will then have some meetings in the afternoon. I should flight back on Thursday, and give my lecture on Predictive Modeling in the evening…

I guess that Saturday will be a mix between two very different activities, enjoying some time with the kids (and watch the game), and trying to finalize my slides…