Qui dit fêtes de fin d’année dit correction de copies… Avant les vacances, il y avait l’examen final du cours de statistique, de Master. L’énoncé est en ligne, ainsi que des éléments de correction.
Monthly Archives: December 2016
L’esprit critique, c’est tellement old-school
“I’m sorry. I did not start it, but I may have contributed to it” disait S.M. Ong au début du mois en évoquant les “fake-news”.
Ces dernières heures, j’ai pas mal échangé sur twitter – ce qui n’est pourtant pas dans mes habitudes – suite à un tweet que j’avais posté, en début de semaine,
“Un touriste suisse a eu la mauvaise idée de brancher son appareil à raclette sur la prise électrique d’un TGV” https://t.co/um6pBo91nk
— Arthur Charpentier (@freakonometrics) 26 décembre 2016
J’avais trouvé l’article amusant, vraiment bien écrit, donc je l’ai mentionné sur Twitter. Tout comme je peux mentionner une carte intéressante sur le site du New York Times ou une statistique intéressante sur le site du Guardian. Après à peine une heure, j’ai été surpris du nombre de personnes qui semblaient avoir pris ça au premier degré. Et ça ne s’est pas démenti avec le temps, même si j’ai posté rapidement un tweet rappelant que scienceinfo.fr était un site parodique…
désolé, après avoir lu certaines réactions à mon tweet, je crois j’aurais du préciser que https://t.co/Hvl9t9WcNk était un site parodique… pic.twitter.com/TzM60SbSFz
— Arthur Charpentier (@freakonometrics) 26 décembre 2016
J’ai reçu un tweet cet après-midi où on tentait de m’apprendre que c’était un “hoax”
@freakonometrics c’est un hoax. Cf https://t.co/j0N9DFUozd …
— Tristan Nitot (@nitot) 28 décembre 2016
ce qui m’a interpelé car pour moi, un “hoax” et un article sur un site parodique, ce sont des choses différentes. Puis on m’a dit que je n’avais pas à twitter des choses pareils, des “fakes”, ayant un compte twitter certifié (tout comme Jean-Marc Morandini, ou…. Le Gorafi)
@freakonometrics Vous avez un compte certifié et vous tombez dans ces fakes … L’histoire de l’appareil a raclette 😂
— Coco (@coco39100) 28 décembre 2016
J’avoue que tous ces échanges m’ont interpelé. Surtout dans le contexte des débats sur les “fake news” qui fleurissent ces jours-ci. Mais comme s’interrogeait S.M. Ong, “what’s fake news and what’s satirical fake news these days“, la distinction est-elle simple à faire ?
Hoax et “fake news”
Les “hoax”, parfois traduit “infaux“, est un mot anglais qui désigne un canular créé notamment à des fins malveillantes, comme le rappelle wikipedia. Ces “hoax”, ou “fake news” nous ont fait entrer dans une ère de post-vérité, “post-truth”, désigné mot de l’année 2016 par l’Oxford Dictionary, comme le notait Alison Flood.
Benoît avait une réflexion intéressante à ce sujet
@freakonometrics Et le mobile accentue cette “réduction” du champ visible (url quasi invisible, environnement très standardisé, rebond fort)
— Benoit F. (@benoitfuric) 28 décembre 2016
Quand on lit un article sur son téléphone cellulaire, on en voit pas le site, le contexte général, on ne voit plus que l’article. On ne voit pas que le site a tendance à publier des articles satiriques. On ne voit même pas le site sur lequel on est arrivé….
Les sites parodiques (ou satiriques)
Ce qui m’a surpris, suite à mon tweet, c’est que si je mentionne un article publié dans Le Gorafi, personne ne crie au scandale, personne ne me rappelle que “c’est un hoax, attention”. Ici, nous ne sommes pas dans le “hoax”, nous sommes dans l’humour, comme beaucoup de journaux le pratiquent, ainsi que le rappelle Isabel Ermida dans “Language and Humour in the Media” (qui revient d’ailleurs sur l’origine du mot “hoax”). L’exemple historique le plus connu est The Onion, en ligne depuis 20 ans (mais qui existait avant sous format papier, depuis 1988). Elle Hunt mentionne d’ailleurs que The Onion n’est pas un site de “fake news”. En France, la référence depuis 2012 est Le Gorafi. Comme le rappelle le site “Tous les articles relatés ici sont faux (jusqu’à preuve du contraire) et rédigés dans un but humoristique”. Quand on lit “Toulouse: il se fait abattre de 46 balles dans le corps pour avoir demandé un ‘pain au chocolat’” on se doute que ce n’est pas sérieux (ou sinon, ça mérite réflexion…) ou encore “il parvient enfin à utiliser le théorème de Pythagore dans sa vie professionnelle“. Peut-on parler de “hoax” ?
En 1998, Lyle Zapato avait fait un test auprès de collégiens qui avaient cru, à l’unanimité (certes, il n’y avait que 25 enfants), un article sur internet qui mentionnait une espèce animale malheurement menacée, les “octopus paxarbolis”, c’est à dire des pieuvres qui vivent dans les arbres, raconté dans “Help Save The Endangered Pacific Northwest Tree Octopus From Extinction!“.
Il est étonnant qu’il soit aussi simple de croire vraie une histoire manifestement fausse. Mais comme le disait Jean-Noël, oui, ça fait sérieux, il y a même une photo de train, et de raclette !
@freakonometrics mais pourtant la photo prouvait que c’était vrai, on voyait un train et un appareil à raclette ! ;-D cc @nitot
— Jean-noël Lafargue (@Jean_no) 28 décembre 2016
Les sites “sérieux”… et l’esprit critique
C’est marrant parce que l’an passé, mon fils avait justement une introduction critique aux médias, au collège. On lui expliquait qu’il y avait des sites “sérieux”. Mais force est de constater que publier dans un journal renommé ne veut rien dire… Par exemple, le même jour que mon tweet sur la raclette, lefigaro.fr publiait “le français se hisserait au 3e rang des langues les plus parlées au monde”
ouch, https://t.co/OZLAtJAKFk je crois qu’il y a une mauvais interprétation du rang (pas “langues les plus parlées” mais en “puissance”)
— Arthur Charpentier (@freakonometrics) 26 décembre 2016
L’article en question parlait d’une “étude” (non citée) mais qui renvoyait vers un article paru dans quebec.huffingtonpost.ca, ce dernier étant un peu plus prudent en annonçant “le français serait la troisième langue du monde selon une nouvelle étude“. Tout en restant assez flou. Mais le titre de l’article de lefigaro.fr est complètement faux ! Le Français reste (relativement) peu parlé, mais reste une langue “puissante”, selon un terme qui n’est pas expliqué (il faut chercher un peu pour trouver plus de détails, sur kailchan.ca). Au même moment, liberation.fr publiait “on s’est aperçu que Google Maps mentait“, repris sur lci.fr. Comme le notait Philippe Reka, dire “on a demandé à des chercheurs chinois, japonais, indiens de faire des tests, et on a pu voir que si vous êtes à Pékin, vous avez une certaine frontière dans l’Himalaya et qu’à Delhi, vous en avez une autre“, c’est un peu fort de café… Il suffit en effet d’aller (depuis son canapé) sur http://google.cn/maps/ ou http://www.google.co.in/maps/. Pour reprendre l’animation de Philippe Rivière, on voit assez simplement que les tracés des frontières diffèrent (comme tout ceux qui jouent avec google maps le savent depuis longtemps)
Mais je m’arrête là, acrimed.org revient très régulièrement sur les “vraies news”, qui n’en sont souvent pas…
Comme le rappellait Benjamin Herold le meilleur antidote à la gratuité de l’information en ligne (mais je rajouterais à la course au buzz que l’on retrouve partout, y compris dans les sites – et les journaux – “sérieux”) est une culture de la pensée critique. C’est également ce que rappelle “Evaluating Information: The Cornerstone of Civic Online Reasoning” paru voilà déjà 5 ans.
Je vais continuer à publier sur mon compte twitter des liens vers des informations sérieuses, mais aussi des articles qui m’ont amusés. Tout comme sur mon blog, lorsque je publie “coffee and productivity” montrant que le café augmente la productivité,
ou “cigarette and life expectancy” montrant qu’avec le même raisonnement, fumer fait vivre plus longtemps…
J’ai envie de croire que les personnes qui suivent mon fil twitter ou qui viennent sur mon blog savent faire preuve d’esprit critique. Et apprécient que l’on soit sérieux, mais aussi satirique, à l’occasion, voie que l’on fasse preuve d’humour… Cela dit, la frontière entre “humour” et “sérieux” est complexe à définir… spurious-correlations c’est sérieux ?
“Divorce rate in Maine correlates with Per capita consumption of margarine (U.S.)” (via http://t.co/tyj58ay5L5) pic.twitter.com/UhFt1zVJeG
— Arthur Charpentier (@freakonometrics) 10 mai 2014
Quantile and Expectile Regression for random effects model
Our joint paper with Amadou Diogo Barry and Karim Oualkacha is now available on https://hal.archives-ouvertes.fr/hal-01421752
Quantile and expectile regression models pertain to the estimation of unknown quantiles/expectiles of the cumulative distribution function of a dependent variable as a function of a set of covariates and a vector of regression coefficients. Both approaches make no assumption on the shape of the distribution of the response variable, allowing for investigation of a comprehensive class of covariate effects. This paper fits both quantile and expectile regression models within a random effects framework for dependent/panel data. It provides asymptotic properties of the underlying model parameter estimators and suggests appropriate estimators of their variances-covariances matrices. The performance of the proposed estimators is evaluated through exhaustive simulation studies and the proposed methodology is illustrated using real data. The simulation results show that expectile regression is comparable to quantile regression, easily computable and has relevant statistical properties. In conclusion, expectiles are to the mean what quantiles are to the median, and they should be used and interpreted as quantilized mean.
Evolution du nombre de demandeurs d’emploi
Allez, on continue avec les projets rendus dans le cadre du projet R de la formation en Data Science pour l’Actuariat, Meryem Schalck a proposé une jolie animation des statistiques des demandeurs d’emploi. Sur http://dares.travail-emploi.gouv.fr/, on récupère les séries d’évolution du nombre de demandeurs d’emploi, par département (entre autres), sous forme de fichier Excel. Histoire de simplifier, je mets en ligne une copie en csv,
ufc<- read.csv("http://freakonometrics.free.fr/DemandeV2.csv", header=TRUE, sep = ";")
Pour charger les packages, Meryem propose une jolie petite fonction
Chargepackages <- function(x){
x <- as.character(x)
if(!require(x, character.only = TRUE)){
install.packages(x, dependencies = TRUE)
require(x, character.only = TRUE)
}
}
packages <- list("sp", "rgdal", "animation", "rgeos", "mapproj", "maptools", "dplyr", "ggplot2", "RColorBrewer", "classInt", "PBSmapping", "ggmap", "splancs", "osmar", "maps", "scales", "osmar", "geosphere", "RCurl", "bitops", "XML", "splancs", "PBSmapping", "htmltools"
)
lapply(packages, Chargepackages)
On va retravailler un peu les données pour commencer
x<-ufc$nombre
y<-ufc$demande
dt<-as.Date(ufc$date,"%d/%m/%Y")
n<-length(x)
et ensuite on fait notre animation,
colfunc <- colorRampPalette(c("green","yellow", "red"))
YlOrBr <- data.frame(col=colfunc(length(unique(y))),num=unique(y)[order(unique(y))])
ufc %>%
left_join(YlOrBr, by=c("demande"="num")) -> YlOrBr
oopt = ani.options(interval = 0.5)
library(scales)
YlOrBr$Date<-as.Date(YlOrBr$date,"%d/%m/%Y")
ani.record(reset = TRUE) # clear history before recording
for (i in 1:length(unique(dt))) {
YlOrBr[i>=x,] %>%
ggplot() + geom_point(aes(x=Date,y=demande),colour=as.character(YlOrBr$col[i>=x]))+
scale_x_date(breaks=date_breaks("6 month"),limits = c(min(dt),max(dt))) +
ylim(c(1000,max(y))) +
ggtitle ("Nombre de demandeurs d'emploi par mois entre 01/2000 et 04/2016") +
xlab("Date") +
ylab("Nombre de demandeurs d'emploi cat A ensemble") +
theme(axis.text.x = element_text(angle=45, vjust=0.5, size=5)) -> p
print(p)
ani.record()
}
oopts = ani.options(interval = 0.2)
ani.record()
saveHTML(ani.replay(),img.name = "record_plot")
et voilà
Localisation des accidents corporels
Encore un billet dans le cadre du projet R de la formation en Data Science pour l’Actuariat, avec quelques lignes de codes pour visualiser les accidents corporels, en France. Il y aura d’autres billets sur ces données dans les jours à venir. Le premier s’inspire du projet de Guillaume Gerber. L’idée est de normaliser par la population d’un département, quand on va visualiser le nombre d’accidents.
On va commencer par récupérer notre fond de carte
library(rgdal)
library(sp)
library(data.table)
library(dplyr)
library(plyr)
download.file(url = "http://professionnels.ign.fr/sites/default/files/GEOFLA_1-1_SHP_LAMB93_FR-ED111.tar.gz",
destfile="GEOFLA.tar.gz")
untar("GEOFLA.tar.gz")
Ces données contiennent la population. Mais par commune. Donc on va aggréger pour l’avoir par département
com <- readOGR(dsn = "./GEOFLA_1-1_SHP_LAMB93_FR-ED111/COMMUNES/", layer="COMMUNE", stringsAsFactors=FALSE,encoding = "UTF-8")
com <- readOGR(dsn = "/home/arthur/Téléchargements/projet-R-DSA-2016/GEOFLA_1-1_SHP_LAMB93_FR-ED111/COMMUNES/", layer="COMMUNE")
dep <- readOGR(dsn = "/home/arthur/Téléchargements/projet-R-DSA-2016/GEOFLA_1-1_SHP_LAMB93_FR-ED111/DEPARTEMENTS/", layer="DEPARTEMENT", stringsAsFactors=FALSE,encoding = "UTF-8")
pop <- as.data.table(tapply(X=com$POPULATION, INDEX = com$CODE_DEPT, FUN = sum), keep.rownames=TRUE)
colnames(pop) <- c("CODE_DEPT","POPULATION")
superficie <- as.data.table(tapply(X=com$SUPERFICIE, INDEX = com$CODE_DEPT, FUN = sum), keep.rownames=TRUE)
colnames(superficie) <- c("CODE_DEPT","SUPERFICIE")
dep@data <- inner_join(dep@data, pop)
dep@data <- inner_join(dep@data, superficie)
dep@data$POPULATION <- dep@data$POPULATION * 1000
On va ensuite récupérer les données d’accidents de la route
debut_url_base_accident = "https://www.data.gouv.fr/s/resources/"
acc_caract <- read.csv(file = paste(debut_url_base_accident,"base-de-donnees-accidents-corporels-de-la-circulation-sur-6-annees/20150806-155035/caracteristiques_2010.csv",sep=''),colClasses=c("com"="character","dep"="character"))
acc_caract <- rbind(acc_caract,read.csv(file = paste(debut_url_base_accident,"base-de-donnees-accidents-corporels-de-la-circulation-sur-6-annees/20150806-154723/caracteristiques_2011.csv",sep=""), colClasses=c("com"="character","dep"="character")))
acc_caract <- rbind(acc_caract,read.csv(file = paste(debut_url_base_accident,"base-de-donnees-accidents-corporels-de-la-circulation-sur-6-annees/20150806-154431/caracteristiques_2012.csv",sep=""), colClasses=c("com"="character","dep"="character")))
acc_caract <- rbind(acc_caract,read.csv(file = paste(debut_url_base_accident,"base-de-donnees-accidents-corporels-de-la-circulation-sur-6-annees/20150806-154105/caracteristiques_2013.csv",sep=""), colClasses=c("com"="character","dep"="character")))
acc_caract <- rbind(acc_caract,read.csv(file = paste(debut_url_base_accident,"base-de-donnees-accidents-corporels-de-la-circulation-sur-6-annees/20150806-153701/caracteristiques_2014.csv",sep=""), colClasses=c("com"="character","dep"="character")))
acc_caract <- rbind(acc_caract,read.csv(file = paste(debut_url_base_accident,"base-de-donnees-accidents-corporels-de-la-circulation/20160909-181230/caracteristiques_2015.csv",sep=""), colClasses=c("com"="character","dep"="character")))
acc_caract$dep[which(acc_caract$dep %in% "201")] <- "2A0"
acc_caract$dep[which(acc_caract$dep %in% "202")] <- "2B0"
acc_caract$dep <- substr(acc_caract$dep, 1, 2)
Maintenant, on peut compter, par année, par département (ou en aggrégeant, temporellement)
dep_with_nb_acc <- function(acc_caract, dep,nb_an=1,normalize=FALSE){
acc_nb_par_dep <- count(acc_caract,"dep")
acc_nb_par_dep$dep <- substr(acc_nb_par_dep$dep, 1, 2)
names(acc_nb_par_dep)[names(acc_nb_par_dep) == "dep"] <- "CODE_DEPT"
dep@data <- inner_join(dep@data, acc_nb_par_dep)
dep@data$SUPERFICIE <- dep@data$SUPERFICIE/100 # en km^2
dep@data$freq_par_hab <- (dep@data$freq/nb_an)/dep@data$POPULATION
m <- sum(dep@data$freq/nb_an)/sum(dep@data$POPULATION)
if(normalize) dep@data$freq_par_hab <- log(dep@data$freq_par_hab/m)
return(dep)
}
On va alors constituer nos bases,
data_plot <- c(
"2010_2015" = dep_with_nb_acc(acc_caract, dep,nb_an = 6),
"2010_2015_n" = dep_with_nb_acc(acc_caract, dep,nb_an = 6,normalize=TRUE))
La première correspond juste au nombre d’accident, sur 6 ans, normalisé par la population (ce qui pourraît être vu comme une fréquence d’accident corporel)
zmax = max(data_plot[[1]]@data$freq_par_hab)
spplot(obj = data_plot$'2010_2015',"freq_par_hab",at = seq(0, zmax, by = zmax/10),main = "")
mais on peut aussi normaliser par la fréquence nationale, afin de faire ressortir les départements les plus dangereux. On est également passé sur une échelle logarithmique,
zmin = min(data_plot[[8]]@data$freq_par_hab)
zmax = max(data_plot[[8]]@data$freq_par_hab)
spplot(obj = data_plot$'2010_2015_n',"freq_par_hab",at = seq(zmin, zmax, by = (zmax-zmin)/10),main = "")
Visualiser l’évolution de la taxe d’habitation
Dans le cadre du projet R de la formation en Data Science pour l’Actuariat, on continue à explorer les données sur data.gouv, par exemple https://data.gouv.fr/fr/datasets/impots-locaux/, sur les impôts locaux. C’est ce que Nicolas Jaudel avait proposé,
library(sp)
tax_hab=read.csv("http://freakonometrics.free.fr/F7815_taxe_hab.csv", col.names=c("CODDPT", "LIBDPT","E2001", "E2002","E2003","E2004","E2005", "E2006","E2007","E2008","E2009", "E2010"), header=FALSE, dec=",", sep=";", stringsAsFactors = FALSE)
On va extraire une information intéressante, par exemple la variation du taux moyen pour la taxe d’habitation, sur 10 ans,
tax_hab$taux_evo=tax_hab$E2010/tax_hab$E2001-1
Pour la carte des départements, on peut aller récupérer un fond de carte sur biogeo.ucdavis.edu.
download.file("http://freakonometrics.free.fr/FRA_adm2.rds","FRA_adm2.rds")
FR=readRDS("FRA_adm2.rds")
donnees_carte=data.frame(FR@data)
donnees_carte=cbind(donnees_carte,"taux_evo"=0)
for (i in 1:96){
donnees_carte[i,"taux_evo"]=tax_hab[tax_hab$CODDPT==donnees_carte$CCA_2[i],'taux_evo']*100
}
On a ainsi toutes les données, et on peut faire une carte
library(cartography)
cols <- rev(carto.pal(pal1 = "red.pal",n1 = 10,pal2="green.pal",n2=10))
plot(FR, col = "grey", border = "gray1", bg="#A6CAE0",xlim=c(-5.2,12))
choroLayer(spdf = FR,
df = donnees_carte,
var = "taux_evo",
breaks = seq(-30,70,5),
col = cols,
legend.pos = "topright",
legend.title.txt = "Evolution",
legend.values.rnd = 2,
add = TRUE)
Le sport, une activité de riches ?
Allez, un petit billet léger pour bien commencer les vacances de Noël. Toujours sur le projet R de la formation en Data Science pour l’Actuariat, Cyril Legrand a proposé de fusionner deux bases, une sur le revenu moyen par foyer (fiscal) par commune, et une sur le nom de licenciés dans les clubs sportifs. Il est faire un peu de retraitement des codes INSEE pour la jointure, car par exemple Marseille est codée sur les arrondissements dans une des bases, sur la ville dans l’autre.
De la Statistique Descriptive à la Statistique Mathématique
J’ai l’impression que la principale difficulté pour rentrer dans un cours de stat(istique) math(ématique) est de comprendre de quoi on parle. Que les grandeurs qui étaient évoquées dans le cours de stat(istique) des(criptive) deviennent aléatoires. On peut alors calculer l’espérance d’une moyenne, ou sa variance, ou la variance d’une variance. Dans mon cours, j’ai toujours pris soin de distinguer les minuscules (correspondant aux observations, aux réalisations de variables aléatoires) des majuscules (correspondant aux variables aléatoires sous-jacentes). Il y a quelques mois, Terry Tao suggérait d’aller un cran plus loin, en notant qu’une variable aléatoire était justement aléatoire. Par exemple, si on considère un tirage de Bernoulli, on a l’échantillon suivant
que l’on peut voir comme une réalisation des variables aléatoires suivantes
Pareil, une statistique peut être soit une fonction construit sur notre échantillon,
ou alors être à son tour une variable aléatoire,
Et comme la moyenne est maintenant vue comme une variable aléatoire, on peut regarder sa loi, sa moyenne, sa variance, etc.
Depuis quelques années, je m’efforçais de faire le cours au tableau, mais je me demande si des slides – avec des animations dès qu’on parle de grandeurs aléatoires – n’auraient pas plus de sens…. A suivre donc (on verra ma motivation à la rentrée).
A quelle distance d’une banque habite-t-on ?
Dans le cadre du projet de R de la formation en Data Science pour l’Actuariat, je vais continuer à mettre en ligne des morceaux de codes qui peuvent être utiles, dans un contexte spatial. Le dernier billet, sur cartographier le vote pour le Brexit, avait été repris (et bien amélioré) sur le site des voisins de rgeomatic. Aujourd’hui, je vais m’inspirer du travail d’Etienne Flichy qui mixe répartition de la population sur le territoire, et localisation des agences bancaires.
On parle des banques ici, mais si on a une base avec les coiffeurs, les boulangeries, etc, on peut faire la même chose ! (autant dire qu’on va pouvoir s’amuser quand la base sirene sera rendue ouverte – dans les semaines à venir). On va supposer que l’on a une base avec toutes les banques géocodées. Bon, pour l’exercice, on va utiliser la localisation des agences bancaires, en utilisant les données de cbanque.com. C’est assez facile d’aller scraper le site, quand on regarde la façon dont sont faites les pages, e.g. http://cbanque.com/pratique/agences/credit-cooperatif/35/. Là on récupère les adresses (postales) et on peut utiliser https://adresse.data.gouv.fr/csv/ (ou différents outils) pour géolocaliser les adresses.
Continue reading A quelle distance d’une banque habite-t-on ?
Rupture Detection
There are some graphs that you cannot forget. One graph that I found puzzling was mentioned on Andrew Gelman’s blog, a few years back, and was related to rupture detection
What I remember from this graph is that if you want to get a rupture, you can easily find one…
ENSAE, Devoir 3
Pour le dernier devoir du cours d’actuariat de l’assurance dommage, comme annoncé, je demande de faire un calcul de provisions, sur les données du fichier triangle_ensae.ods.
Je veux un rapport, de 5 pages, proposant un montant de provisions pour les 14 années passées, mais aussi pour la dernière année en cours (pour laquelle 103 733 € ont déjà été payés), avec deux méthodes différentes (montant dit best estimate). Dans un second temps, il faudra utiliser une méthode de simulations pour estimer le quantile à 99.5% des paiements à venir, toutes années confondues. Autrement dit, on veut savoir combien de provisions il convient de constituer pour pouvoir indemniser intégralement les assurés dans 99.5% des cas.
Cartographier le vote pour le Brexit
Je suis, ces temps-ci plongé dans les projets de R que j’avais donnés pour la formation Data Science pour l’Actuariat. Comme j’ai eu plein de choses intéressantes, je me suis dit que je pourrais tenter de reprendre dans des billets sur le blog (en plus, ça me permet de vérifier au passage que les codes tournent).
Pour commencer la série (oui, il y en aura d’autres), Flavien Thery a proposé une cartographie du vote pour le Brexit. La première étape était de récupérer un fond de carte (on va utiliser la seconde zone administrative ici)
library(sp)
library(raster)
download.file("http://biogeo.ucdavis.edu/data/gadm2.8/rds/GBR_adm2.rds","GBR_adm2.rds")
UK=readRDS("GBR_adm2.rds")
UK@data[159,"HASC_2"]="GB.NR"
plot(UK, xlim = c(-4,-2), ylim = c(50, 59), main="UK areas")
Cette carte est un peu surprenante… parce qu’on est sur le Royaume uni, ce qui inclu l’Angleterre, le Pays de Galles, l’Ecosse et l’Irlande du Nord. Et l’Irlande du Sud n’est pas sur la carte… Autant la rajouter (ça sera plus joli si on met du bleu pour l’eau)
download.file("http://biogeo.ucdavis.edu/data/gadm2.8/rds/IRL_adm0.rds","IRL_adm0.rds")
IRL=readRDS("IRL_adm0.rds")
plot(IRL,add=TRUE)
Mais on commence à rajouter l’Irlande, pourquoi pas la France, qui est juste en bas à droite, et qu’on devrait voir un peu…
download.file("http://biogeo.ucdavis.edu/data/gadm2.8/rds/FRA_adm0.rds","FRA_adm0.rds")
FR=readRDS("FRA_adm0.rds")
plot(FR,add=TRUE)
Bon, arrêtons là….
On peut ensuite récupérer les données du référendum sur le Brexit (que je stocke sur le blog histoire de gagner un peu de temps)
loc="http://freakonometrics.hypotheses.org/files/2016/12/EU-referendum-result-data.csv"
referendum=read.csv(loc,header=TRUE,dec=".",sep=",",stringsAsFactors = FALSE)
referendum=referendum[c(3,6,13,14)]
library(plyr)
referendum=ddply(referendum,.(Region,HASC_code),summarise,Remain=sum(Remain),Leave=sum(Leave))
On peut retrouver que le Brexit a gagné, avec 51,89% des votes exprimés (ce qui est conforme à ce que dit wikipedia)
> sum(referendum$Leave)/(sum(referendum$Leave)+sum(referendum$Remain))
[1] 0.5189184
On peut regarder, région par région, si c’est le leave ou le remain qui a gagné, en utilisant
referendum=referendum[c(referendum$Region!="Northern Ireland"),]
referendum=referendum[c(referendum$HASC_code!="Gibraltar"),]
row.names(referendum)=seq(1,nrow(referendum),1)
leave_or_remain=cbind(referendum,"Brexit?"=0)
leave_or_remain[,"Brexit?"]=ifelse(leave_or_remain$Remain<leave_or_remain$Leave,rgb(1,0,0,.7),rgb(0,0,1,.7))
map_data=data.frame(UK@data)
map_data=cbind(map_data,"Brexit"=0)
for (i in 1:nrow(map_data)){
if(map_data[i,"NAME_1"]=="Northern Ireland"){
map_data[i,"Brexit"]="blue"}else{
map_data[i,"Brexit"]=as.character(leave_or_remain[leave_or_remain$HASC_code==map_data$HASC_2[i],'Brexit?'])
}
}
plot(UK, col = map_data$Brexit, border = "gray1", xlim = c(-4,-2), ylim = c(50, 59), main="How the UK has voted?", bg="#A6CAE0")
plot(IRL, col = "lightgrey", border = "gray1",add=TRUE)
plot(FR, col = "lightgrey", border = "gray1",add=TRUE)
legend(-1,59,c("Leave","Remain"),fill=c("red","blue"),bty="n")
(on rajoute une petite légende pour être plus clair). Mais on peut aller plus loin en représentant le pourcentage obtenu, par région, par le camp du leave. Pour ça, on peut utiliser le package cartography
library(cartography)
cols <- carto.pal(pal1 = "red.pal",n1 = 5,pal2="green.pal",n2=5)
map_data=data.frame(UK@data)
map_data=cbind(map_data,"Percentage_Remain"=0)
for (i in 1:nrow(map_data)){
if(map_data[i,"NAME_1"]=="Northern Ireland"){
map_data[i,"Percentage_Remain"]=55.78}else{ map_data[i,"Percentage_Remain"]=100*leave_or_remain[leave_or_remain$HASC_code==map_data$HASC_2[i],"Remain"]/(leave_or_remain[leave_or_remain$HASC_code==map_data$HASC_2[i],"Remain"]+leave_or_remain[leave_or_remain$HASC_code==map_data$HASC_2[i],"Leave"])
}
}
plot(UK, col = "grey", border = "gray1", xlim = c(-4,-2), ylim = c(50, 59),bg="#A6CAE0")
plot(IRL, col = "lightgrey", border = "gray1",add=TRUE)
plot(FR, col = "lightgrey", border = "gray1",add=TRUE)
choroLayer(spdf = UK,
df = map_data,
var = "Percentage_Remain",
breaks = seq(0,100,10),
col = cols,
legend.pos = "topright",
legend.title.txt = "",
legend.values.rnd = 2,
add = TRUE)
(là encore, une légende permet d’aider à lire la carte). Amusant, n’est ce pas ?
Testing for Extreme Volatility Transmission
This morning, I will be at the 3rd International MACroeconomics workshop (IMAC). I will discuss Testing for Extreme Volatility Transmission with Realized Volatility Measures by Sessi Tokpavi (University of Orléans). Slides of the discussion are online.
Habilitation à Diriger des Recherches
« L’habilitation à diriger des recherches sanctionne la reconnaissance du haut niveau scientifique du candidat, du caractère original de sa démarche dans un domaine de la science, de son aptitude à maîtriser une stratégie de recherche dans un domaine scientifique ou technologique suffisamment large et de sa capacité à encadrer de jeunes chercheurs. Elle permet notamment d’être candidat à l’accès au corps des professeurs des universités. »
Ma défense était ce soir, ça y est, c’est fait… on va pouvoir faire avancer la recherche maintenant !
Actuariat de l’Assurance non-Vie #11
Dernier jeu de slides, pour compléter la section sur le provisionnement, avec des méthodes plus économétriques. Les slides sont enfin en ligne.
Je complèterais éventuellement s’il y a des questions spécifiques, sur un point, ou un autre.