Tag Archives: multinomiale

Ex-æquo aux élections ?

La saison 1 de Mad Men se termine avec l’élection de 1960, Nixon contre Kennedy, qui s’est jouée dans un mouchoir de poche.

L’autre jour, Benoit Rittaud évoquait sur http://images.math.cnrs.fr/ la possibilité que des candidats soient ex-æquo dans une élection, et les aspects légaux (car visiblement la constitution française n’a pas vraiment prévu ce cas). Le plus intéressant aurait été – à mon avis (mais mon avis est fortement biaisé) – de calculer la probabilité qu’un tel évènement survienne. Par exemple, qu’à l’issu du premier tour, le deuxième et le troisième candidat soient ex-æquo, de telle sorte qu’on ne puisse déterminer qui ira au second tour. Justement, un des commentaire affirme “la probabilité qu’un tel événement arrive est environ l’inverse du nombre typique de voix séparant le deuxième homme (ou femme) du troisième”. N’ayant réussi à démontrer ce résultat, je me suis dit qu’on pourrait se lancer dans des simulations pour quantifier cette probabilité.

Commençons par supposer qu’il y a quatre candidats (et la possibilité de s’abstenir). On suppose que le taux d’abstention est de 20%, et que les intentions de vote pour chacun des candidats soit de 20%.

EXE=function(N=100,ns=1000000){
ExE=rep(NA,ns)
for(s in 1:ns){
M=sample(c("A","B","C","D","E"),size=N,prob=rep(.2,5),
replace=TRUE)
tb=table(M)
Ms=rev(sort(tb[-(names(tb)=="A")]))
ExE[s]=(Ms[2]==Ms[3])
}
return(mean(ExE))}

Le soucis est que faire plusieurs centaines de millions de tirages de plus de 40 millions de lois multinomiales pourrait prendre du temps (genre vraiment beaucoup). Une stratégie intuitive peut être de faire des simulations pour une population plus faible, puis d’extrapoler.

E=function(N) EXE(N)
taille=c(20,50,100,200,500,1000,2000)
proba=Vectorize(E)(taille)
plot(taille,proba,type="b",log="xy")

base=data.frame(y=proba,x=taille)
reg=lm(log(y)~I(log(x)),data=base)
s2=summary(reg)$sigma
m=exp(predict(reg,newdata=
data.frame(x=41194689)))*exp(s2^2/2)

L’extrapolation linéaire n’étant peut-être pas parfaite, on peut tenter un ajustement quadratique (comme cela a été fait sur la figure ci-dessus)

reg2=lm(log(y)~I(log(x))+I(log(x)^2),data=base)
s2=summary(reg2)$sigma
m2=exp(predict(reg2,newdata=data.frame(x=41194689)))*exp(s2^2/2)

La différence de prédiction avec les deux modèles est significative comme on le voit sur le graphique. Toutefois, avoir des ex-æquo avec une très très grande population est tellement rare qu’il faut réellement beaucoup de simulations pour avoir un ordre de grandeur valide. L’extrapolation semble être une alternative valide ici (mais je suis ouvert à toute critique, ou discussion).

Maintenant que l’on a une méthodologie (ou presque), utilisons des ordres de grandeurs réalistes. Par exemple en se basant sur l’élection présidentielle de 2002. 41.19 millions d’électeurs, 28.4% d’abstentions, et pour les votes exprimés, un candidat en tête avec 19.9% des votes, un second à 16.8% et le troisième à 16.2% (et les autres plus loin). Les intentions de votes pour les principaux candidats figurent dans la loi multinomiale ci-dessous

EXE=function(N=100,ns=1000){
ExE=rep(NA,ns)
for(s in 1:ns){
M=sample(c("abs","JC","JMLP","LJ","FB","AL","JPC"
,"NM","autres"),size=N,prob=c(.284,.1988*.716,
.1686*.716,.1618*.716,.0684*.716,.0572*.716,
.0533*.716,.0525*.716,.1714104),replace=TRUE)
tb=table(M);
Ms=rev(sort(tb[-(names(tb)%in%c("abs","autres"))]));
ExE[s]=(Ms[2]==Ms[3])
}
return(mean(ExE))}

Si on en croit l’affirmation du commentaire, la probabilité d’avoir eu des ex-æquo (entre le second et le troisième candidat) serait de l’ordre de

> 1/(41194689*(.1686*.716-.1618*.716))
[1] 4.985823e-06

Je dois avouer que les calculs n’ont rien donné (on verra tout à l’heure ce qui se passe si on change la définition d’ex-æquo). Numériquement, les calculs tournent encore (après 4 jours) sur le serveur… Par contre, si la population ne contenait que 10000 électeurs, la probabilité d’avoir des ex-æquo (entre le 2ème et le 3ème candidat) serait de l’ordre de 3 sur 10000 (avec un intervalle de confiance que je n’ai toujours pas calculé)

> EXE(10000,100000000)
[1] 3.245562e-05

A titre d’illustration, on peut regarder ce qui s’est passé, commune par commune en France (les données sont en ligne sur le site du Ministère de l’Intérieur). L’exercice est biaisé parce qu’il y a de l’hétérogénéité spatiale… mais ça permet toujours d’illustrer… Les nombres de voix obtenues, pour le second et le troisième candidat, par commune est le suivant

> election=read.table(
+ "/Users/UQAM/Documents/Pres2002Tour1.csv",sep=";")
+ election=election[-1,]
> VOIX=election[,seq(14,74,by=4)]
> DEUXIEME=TROISIEME=TOTAL=rep(NA,nrow(VOIX))
> for(i in 1:nrow(VOIX)){
+ u=rev(sort(as.numeric(as.character(
+ election[i,seq(14,74,by=4)]))))
+ TOTAL[i]=sum(u)
+ DEUXIEME[i]=u[2]
+ TROISIEME[i]=u[3]
+ }

Sur 37000 communes,

> nrow(election)
[1] 37513

ou plutôt sur les 24000 communes qui ont eu plus de 5000 votes exprimés

> sum((TOTAL>5000))
[1] 24230

il y a eu 56 cas avec de deuxième et troisième candidats ex-aequo,

> base=data.frame(TOTAL,DEUXIEME,TROISIEME,COMMMUNE=election[,3])
> selection=(TOTAL>5000)
> sousbase=base[selection,]
> sum(sousbase$DEUXIEME==sousbase$TROISIEME)
[1] 56

On a même le détail,

> sousbase[indice,]
TOTAL DEUXIEME TROISIEME                     COMMMUNE
1134   5819      699       699                        Crouy
1507   5386      717       717          Louroux-Bourbonnais
2603   6235     1151      1151                  Saint-Peray
2916   6346      935       935                        Coucy
2970   7529     1413      1413                   Haraucourt
5089   5567      717       717                      Vignats

(etc). On peut alors regarder, en fonction de la taille des communes (ici le nombre de votes exprimés) la probabilité d’avoir un ex-æquo, dans une commune donnée. On va subdivisionner en 50 classes (pour commencer), et compter à chaque fois le nombre d’ex-æquo,

> n=50
> niveau=seq(1/(2*n),1-1/(2*n),by=1/n)
> base$CLASSE=cut(TOTAL,quantile(TOTAL,seq(0,1,by=1/n)),
+ labels=niveau)
> base$EXAEQUO=(base$DEUXIEME==base$TROISIEME)*1
> m=by(base$EXAEQUO, base$CLASSE, sum)

Par exemple pour les 2% des plus grosses villes, on a eu 5 cas de communes avec des ex-aequos

------------------------------------------------------------
base$CLASSE: 0.99
[1] 5

Si on regarde plus en détails, ce sont les communes avec plus de 7500 votes exprimés,

> quantile(TOTAL,.98)
98%
7427.76
> selection=(TOTAL>quantile(TOTAL,.98))
> sousbase=base[selection,]
> sum(sousbase$DEUXIEME==sousbase$TROISIEME)
[1] 5
> indice=which(sousbase$DEUXIEME==sousbase$TROISIEME)
> sousbase[indice,]
TOTAL DEUXIEME TROISIEME       COMMMUNE CLASSE EXAEQUO
2970   7529     1413      1413     Haraucourt   0.99       1
14359  7496     1392      1392          Renac   0.99       1
25028  8104     1381      1381         Heloup   0.99       1
30014  8087     1340      1340         Barnay   0.99       1
36770  7435      764       764 Pont-sur-Yonne   0.99       1

On peut représenter ces nombres de cas d’ex-æquo et faire une régression de Poisson (car on fait du comptage ici),

> taille=length(TOTAL)/n
> plot(quantile(TOTAL,niveau),m/taille)
> sb=data.frame(niveau,exaequo=as.vector(m))
> library(splines)
> reg=glm(exaequo~bs(niveau),family=poisson,data=sb)
> u=seq(0,1,by=.001)
> lines(quantile(TOTAL,u),predict(reg,
+ newdata=data.frame(niveau=u),
+ type="response")/taille,lwd=2,col="red")

pour des facilités de lecture, on a en ordonnées des probabilités et en abscisse les tailles moyennes des villes

On peut faire la même chose avec 200 classes,

Bref, pour une commune de taille moyenne (disons moyenne supérieure, avec 5000 votes exprimés), on a 2 chances sur 1000 d’avoir un second et un troisième candidat ex-æquo. Ce qui laisserait entendre que pour 40 millions, la probabilité devrait être faible, voire très faible…

Bon, au delà du problème conceptuel des ex-æquo (ou légal, car je doute qu’il soit possible d’avoir une égalité parfaite, le Ministère de l’Intérieur ferait un recompte et trouverait forcément un résultat différent), il semble que la probabilité semble beaucoup plus faible que celle obtenue par la règle ad-hoc….

Mais peut-être serait-il réaliste de dire que l’élection n’a pas su dégager une majorité (c’est le but d’une élection il me semble) si le nombre de voix séparant deux candidats est trop faible. Par exemple, si moins de 100 voix séparent les candidats, on refait un vote ! Pour revenir sur mon exemple précédant, avec 10000 votants

> EXE100(10000,1000000)
[1] 0.012641

soit une chance sur 100 (je me suis autorisé à faire moins de simulations car l’évènement me semble moins improbable… donc il y a besoin de moins de simulations pour l’observer à l’occasion). Si on va plus loin, et que l’on fait des simulations que l’on projette, comme auparavant, pour estimer la probabilité d’avoir des ex-æquo non pas pour une ville, mais au niveau de la France entière, on est plutôt aux alentours de

> (m=exp(predict(reg,newdata=
+ data.frame(x=41194689)))*exp(s2^2/2))
1
1.505422e-12

pour un modèle linéaire, et

> (m=exp(predict(reg2,newdata=
+ data.frame(x=41194689)))*exp(s2^2/2))
1
4.074588e-98

pour un modèle quadratique (avec la confiance – relativement faible – que l’on peut accorder à une telle projection). Ce qui donne une différence énorme… Si on visualise tout ça,

lines(c(1,100000000),exp(predict(reg,newdata=
data.frame(x=c(1,100000000))))*exp
(s2^2/2),lty=2,col="light blue")

que l’on représente en bleu, en vert un ajustement quadratique

lines(10^(seq(1,8,by=.1)),exp(predict(reg,newdata=
data.frame(x=10^(seq(1,8,by=.1)))))*exp
(s2^2/2),lty=2,col="dark green")

alors que le calcul suggéré dans le commentaire donnerait la courbe rouge,

lines(c(1,100000000),1/(c(1,100000000)*
(.1686*.716-.1618*.716)),col="red",lty=2)

Peut-on parler d’évènement improbable ? A titre de comparaison, quand on parle de quelque chose d’improbable, ou qui arrive avec une probabilité infinitésimale, on pense au loto. A l’ancien loto (français), il fallait choisir 6 numéros parmi 49 (sans remise). La probabilité de gagner le gros lot était de

Pour l’Euromillions, il faut choisir 5 numéros parmi 50 numéros auxquels il faut ajouter les combinaisons possibles avec les Étoiles, soit 2 numéros à choisir parmi 9. La probabilité de gagner le gros lot est alors de

Bref, tomber sur des ex-æquo lors d’une élection nationale, en France, c’est a priori plus rare que gagner le gros lot en jouant une fois au loto… On peut donc dire qu’il ne devrait pas y avoir d’ex-æquo à l’issu du premier tour !

Sondages et probabilités en trois graphiques

Petite précision liminaire: je suis nul en sondage. Je n’ai jamais rien compris aux sondages. C’est la matière où j’ai eu ma pire note quand j’étais étudiant. Mais si on voit un sondage comme l’estimation d’une probabilité de succès, dans des tirages de lois binomiales (ou multinomiales), alors je me sens plus à l’aise pour commenter des choses lues ici ou là. Car en ces temps d’élections (qui approchent), on lit beaucoup de choses sur les sondages ou sur les intervalles de confiance des sondages, voire sur des comparaison de sondages… mais rares sont les justifications. Donc on va creuser un peu (histoire de voir ce qui peut être dit, dans un cadre largement simplifié).

  • “Il y a 1 chance sur 4 qu’un sondage donne une variation d’un point en + ou en – sur une proportion inchangée”

comme le disait @adelaigue (sur Twitter). Faire un sondage politique (si on résume simplement) c’est estimer la probabilité https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic01.gif de voter pour A, sachant que l’on a alors une probabilité https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic02.gif pour qu’une personne prise au hasard ne vote par pour A. Sur https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic03.gif personnes interrogées, le nombre de personnes qui déclarera voter pour A est https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic04.gif, qui suit une loi binomiale de paramètres https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic03.gif et https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic01.gif. Or par approximation (ou convergence) de la loi binomiale par une loi Gaussienne, on peut dire que

https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic05.gif

Or l’estimateur naturel de https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic01.gif est https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic06.gif, i.e.

https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic07.gif

Aussi, on en déduit l’intervalle de confiance classique, tel que https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic01.gif appartient à cet intervalle avec une probabilité https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic11.gif,

https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic08.gif

https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic12.gif désigne la fonction quantile de la loi normale centrée réduite; ou encore, en remplaçant les quantités inconnues par des estimateurs,

https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic09.gif

Les bornes de l’intervalle de confiance de niveau https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic11.gif sont alors

https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic10.gif

C’est ce qu’on peut visualiser sur la figure ci-dessous,

Comme le notait @adelaigue, si on souhaite avoir

https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic13.gif

(car on parle d’un point de base 100), il faudrait que

https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic14.gif

Classiquement, https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic03.gif est aux alentours de 1000. Et ce qu’on cherche c’est une probabilité d’appartenir à un intervalle, i.e. https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic11.gif,

https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adi15.gif

Graphiquement, en fonction de https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic01.gif on a alors

p=seq(0,1,by=.0025)
plot(p,2*(1-pnorm(.01/sqrt(p*(1-p)/1000))))

(ca sera notre premier graphique). Effectivement, entre 35% et 65%, il y a environ 50% de chances d’être autour de la fréquence observée, à un point près. Ou 50% de chances d’avoir plus d’un point de différence. Autrement dit, si on a observé avec 1000 personnes une fréquence de 40%, il y a une chance sur 4 que la probabilité soit entre 39% et 40%, une chance sur 4 qu’elle soit entre 40% et 41%, et surtout une chance sur 4 qu’elle soit inférieure à 39%, et 1 chance sur 4 qu’elle soit supérieure à 41%. Comme le notait @adelaigue, sur les bords, il semble que la formule ne soit pas valide. Mais si on prend justement les cas limites (probabilité très proche de 0, voire franchement nulle), est-ce qu’on peut dire quelque-chose ?

  • peut-on dire quelque chose (malgré tout) quand on ne peut rien dire ?

ou plus précisément, si en interrogeant https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic03.gif personnes et que personne n’a déclaré vouloir voter pour A, quelle serait l’intervalle de confiance de la probabilité de voter pour A ? J’avais évoqué ce point dans un très vieux billet. La première piste est la réponse apporté par Pearson-Clopper (Pearson & Clopper (1934)), basée sur l’utilisation de la loi binomiale, et non plus d’une approximation.

library(Hmisc)
bornesup=function(N){ binom.test(x=0, n=N, p = 0,
alternative="less")$conf.int[2]*100 }
BS=Vectorize(bornesup)
vn=seq(10,500,by=10)
plot(vn,BS(vn))

Une alternative est de faire du bayésien. On utilise le fait que la loi conjuguée de la loi binomiale est la loi beta (comme évoqué dans des précédents billets). Si on suppose que la loi a priori de https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic01.gifest une loi alors la loi a posteriori, sachant que sur https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic03.gif lancés, on n’a observé que des 0, devient

bornesup=function(N){ qbeta(.95,1,1+N)*100 }
BS=Vectorize(bornesup)
lines(vn,BS(vn),lwd=2,col="blue")

La courbe ci-dessus correspond au calcul dans le cas du modèle de Pearson-Clopper, en rouge, et avec le quantile de la loi Beta en bleu (qui sont confondues sur ce graphique). Autrement dit, avec 300 personnes interrogées, si aucune ne dit vouloir voter pour A, l’intervalle de confiance pour https://f.hypotheses.org/wp-content/blogs.dir/253/files/2015/12/adic01.gif est [0%;1%] (environ). Ca sera notre second graphique.

  • si on a trois alternatives, quelles sont les vraies probabilités de gagner ?

c’était déjà sur cette question que je concluais mon précédant billet, mais il est possible d’exploiter davantage cette idée de modèle bayésien avec un a priori suivant une loi beta (comme au dessus) ou de Dirichlet (si on a plus de deux alternatives). La première chose que l’on peut noter, c’est que s’il y a trois alternatives (A, B et C), et que la probabilité de voter pour A est 40%, la probabilité que A gagne… dépend des autres probabilités. On ne peut pas faire comme si A était seul, contre tous les autres:
– si les probabilités de voter pour B et C sont de l’ordre de 30%, on se dit que A a de fortes chances de gagner (oui, pour l’instant ce n’est pas très quantitatif).
– si les probabilités de voter pour B et C sont respectivement 55% et 5%, alors on se dit que A a de faibles chances de gagner.
Plus précisément, si on souhaite savoir si A va l’emporter ou pas, et que 200 personnes ont été interrogées, on peut essayer de faire une abaque, avec des courbes de niveau… dans le simplexe. Pour rappel, les scores obtenus par A, B et C se lisent sur le triangle ci-dessous (le point représenté correspond au triplet (40%,35%,25%))

Les probabilités sont les suivantes, avec n=200 (personnes qui se seraient exprimés)

où la région rouge correspond à une (quasi) certitude de gagner, pour A, et la région blanche, au contraire, une probabilité infinitésimale de gagner.

Ces résultats ont été obtenus par simulation,
vx=v[i]
vz=v[j]
vy=1-vz-vx
alpha=c(vx,vy,vz)
alpha=c(alpha,1-sum(alpha))*200
RD=rdirichlet(ns, alpha)
PROBA[i,j]=mean(RD[,1]>=apply(
cbind(RD[,3],RD[,2]),1,max))

Pour ceux que ne sont pas très à l’aise avec les excursions dans le simplexe, on peut se ramener dans un repère orthonormé. Ici, les scores obtenus par A et B se lisent sur le triangle ci-dessous (le point représenté correspond au triplet (40%,35%,25%), le même qu’auparavant

On a alors les probabilités de gagner suivantes, avec à droite la représentation sur le simplexe, en zoomant autour de notre triplet,

 

(ça sera le troisième et dernier graph). Si on regarde plus en détails, l’intervalle de confiance sur la probabilité de gagner dépend de la région où on se trouve dans le simplexe (un peu comme l’intervalle de confiance dépendait de la probabilité de l’emporter, dans la première partie de ce billet).

Maintenant, comme le rappelait @Christian dans un commentaire, tout cela est la version très simpliste de la pratique des sondages. Dans les vrais sondages (et pas des exercices de probabilités comme on vient d’en faire ici), on utilise des quotas, on fait du redessement, et c’est beaucoup plus complexe (enfin… c’est ce que j’essaye de me dire, pour me convaincre que la mauvaise note que j’avais obtenue lorsque j’étais étudiant était motivée).

Cochrane, Pearson et le test du chi-deux

En cours, nous avons poursuivi sur la loi multinomiale, et le test du chi-deux. Je voulais mettre un petit billet pour récapituler les différents points, et montrer une application numérique (nous reviendrons en détails mercredi sur des applications des outils vus jusqu’à présent).

  • Inférence avec la loi multinomiale

On suppose qu’une variable http://freakonometrics.blog.free.fr/public/maths/coch-01.gif peut prendre http://freakonometrics.blog.free.fr/public/maths/coch-02.gif modalités, notées http://freakonometrics.blog.free.fr/public/maths/coch-03.gif, avec http://freakonometrics.blog.free.fr/public/maths/coch-04.gif. On posera

http://freakonometrics.blog.free.fr/public/maths/coch-05.gif

en notant que http://freakonometrics.blog.free.fr/public/maths/coch-06.gif appartient au simplexe de http://freakonometrics.blog.free.fr/public/maths/coch-07.gif au sens où

http://freakonometrics.blog.free.fr/public/maths/coch-08.gif

On a vu que l’estimateur du maximum de vraisemblance s’obtenait en faisant un peu d’optimisation sous contrainte, et que

http://freakonometrics.blog.free.fr/public/maths/coch-10.gif

(en reprenant les notations du cours). On avait montré que

http://freakonometrics.blog.free.fr/public/maths/coch-11.gif

et on a vu

http://freakonometrics.blog.free.fr/public/maths/coch-13.gif

(ce qui peut se retrouver en introduisant la variable binomiale http://freakonometrics.blog.free.fr/public/maths/coch-16.gif). Mais plus généralement, on finira les calculs permettant d’établir que, pour http://freakonometrics.blog.free.fr/public/maths/coch-17.gif

http://freakonometrics.blog.free.fr/public/maths/coch-18.gif

(ce qui permet d’obtenir la matrice de variance covariance de http://freakonometrics.blog.free.fr/public/maths/coch-20.gif).
En utilisant le théorème central limite on peut montrer que

http://freakonometrics.blog.free.fr/public/maths/coch-23.gif

Sous une forme multivariée, on écrira

http://freakonometrics.blog.free.fr/public/maths/coch-25.gif
http://freakonometrics.blog.free.fr/public/maths/coch-26.gif avec http://freakonometrics.blog.free.fr/public/maths/coch-27.gif et pour http://freakonometrics.blog.free.fr/public/maths/coch-17.gif, http://freakonometrics.blog.free.fr/public/maths/coch-28.gif.

On notera que la somme de la ième colonne de http://freakonometrics.blog.free.fr/public/maths/coch-29.gif est http://freakonometrics.blog.free.fr/public/maths/coch-30.gif.
Aussi, http://freakonometrics.blog.free.fr/public/maths/coch-29.gif n’est pas inversible (c’est le fait que notre paramètre appartient au simplexe).
Pour s’en sortir, la première idée est d’utiliser un peu d’algèbre linéaire. Une matrice http://freakonometrics.blog.free.fr/public/maths/coch-31.gif est une matrice de projection si elle est idempotente, i.e. http://freakonometrics.blog.free.fr/public/maths/coch-32.gif. Ses valeurs propres sont alors 0 ou 1, et si http://freakonometrics.blog.free.fr/public/maths/coch-35.gif est le nombre de fois où 1 est valeur propre, et si http://freakonometrics.blog.free.fr/public/maths/coch-36.gif, alors http://freakonometrics.blog.free.fr/public/maths/coch-37.gif.
Posons http://freakonometrics.blog.free.fr/public/maths/coch-38.gif. Alors

http://freakonometrics.blog.free.fr/public/maths/coch-39.gif

Or compte tenu de la forme de http://freakonometrics.blog.free.fr/public/maths/coch-29.gif,

http://freakonometrics.blog.free.fr/public/maths/coch-40.gif

qui est une matrice de projection dont la trace est http://freakonometrics.blog.free.fr/public/maths/coch-41.gif (qui est aussi le nombre de fois où 1 est valeur propre). Donc

http://freakonometrics.blog.free.fr/public/maths/coch-42.gif

Le test du chi-deux sera basé sur le fait qu’asymmptotiquement

http://freakonometrics.blog.free.fr/public/maths/coch-44.gif

Une autre idée consiste à construire http://freakonometrics.blog.free.fr/public/maths/coch-41.gif variables aléatoires qui seront indépendantes. Mais on peut plutôt regarder les applications de ce test, en particulier comme test d’indépendance.
Pour information, Frank Yates a proposé un correction “pour continuité“, ici. La statistique considérée est alors

http://fre<br /><br /><br /><br />
akonometrics.blog.free.fr/public/maths/coch-46.gif
  • Retour sur le théorème de Cochrane

Soit http://freakonometrics.blog.free.fr/public/maths/coch-50.gif de dimension http://freakonometrics.blog.free.fr/public/maths/coch-51.gif. Posons http://freakonometrics.blog.free.fr/public/maths/coch-59.gif, pour http://freakonometrics.blog.free.fr/public/maths/coch-04.gif, où on notera http://freakonometrics.blog.free.fr/public/maths/coch-60.gif le rang de http://freakonometrics.blog.free.fr/public/maths/coch-62.gif, en supposant que les http://freakonometrics.blog.free.fr/public/maths/coch-62.gif sont positives semidefinies, alors on a équivalence entre

  •  http://freakonometrics.blog.free.fr/public/maths/coch-63.gif
  • http://freakonometrics.blog.free.fr/public/maths/coch-64.gif pour http://freakonometrics.blog.free.fr/public/maths/coch-04.gif,
  • les http://freakonometrics.blog.free.fr/public/maths/coch-65.gif sont des variables indépendantes.

Les http://freakonometrics.blog.free.fr/public/maths/coch-65.gif s’interprètent comme des longueurs (euclidienne) de projections d’un vecteur Gaussien sur des sous-espaces orthogonaux (de dimension respective http://freakonometrics.blog.free.fr/public/maths/coch-60.gif). Si ces vecteurs sont indépendants, et suivent des lois du chi-deux à http://freakonometrics.blog.free.fr/public/maths/coch-60.gif degrés de libertés, avec http://freakonometrics.blog.free.fr/public/maths/coch-63.gif, alors les sous-espaces sont orthogonaux, et supplémentaires. On peut y voir une espèce d’extension du théorème de Pythagore, en remplaçant la notion de vecteurs orthogonaux par des variables indépendantes suivant des lois du chi-deux, et la somme des carrés des longueurs devient la somme des degrés de liberté.

  • Application comme test d’indépendance

Considérons deux variables http://freakonometrics.blog.free.fr/public/maths/coch-66.gif pouvant prendre toutes les deux deux modalités (disons deux lois binomiales). On est alors face a une loi multinomiale à 4 modalités

  • http://freakonometrics.blog.free.fr/public/maths/coch-79.gif avec probabilité http://freakonometrics.blog.free.fr/public/maths/coch-70.gif
  • http://freakonometrics.blog.free.fr/public/maths/coch-78.gif avec probabilité http://freakonometrics.blog.free.fr/public/maths/coch-73.gif
  • http://freakonometrics.blog.free.fr/public/maths/10gif.gif avec probabilité http://freakonometrics.blog.free.fr/public/maths/coch-74.gif
  • http://freakonometrics.blog.free.fr/public/maths/coch-77.gif avec probabilité http://freakonometrics.blog.free.fr/public/maths/coch-75.gif

Un test d’indépendance revient à tester si la loi multinomiale peut s’écrire

http://freakonometrics.blog.free.fr/public/maths/chi2-ab.gif
http://freakonometrics.blog.free.fr/public/maths/chi_ab2.gif
http://freakonometrics.blog.free.fr/public/maths/chi_ab3.gif
http://freakonometrics.blog.free.fr/public/maths/chi_ab4.gif

pour des vecteurs http://freakonometrics.blog.free.fr/public/maths/chi-a.gif et http://freakonometrics.blog.free.fr/public/maths/chi-b.gif. tels que http://freakonometrics.blog.free.fr/public/maths/chi-a2.gif et http://freakonometrics.blog.free.fr/public/maths/chi-b2.gif. On a alors http://freakonometrics.blog.free.fr/public/maths/chi212121.gif contraintes sur les paramètres. Ces deux contraintes additionnelles font que la statistique de test s’écrit

http://freakonometrics.blog.free.fr/public/maths/CHI-INDEP.gif

qui va suivre asymptotiquement une loi http://freakonometrics.blog.free.fr/public/maths/CHI1.gif i.e. http://freakonometrics.blog.free.fr/public/maths/CHI12.gif d’après le théorème de Cochrane.

  • Peine de mort pour les condamnées pour meurtre en Floride 1976-1987

en fonction de la “race” du meurtrier et de la victime,

  • meurtrier de “race blanche” et victime de “race blanche“: 53 condamnés à mort, 414 non condamnés à mort
  • meurtrier de “race blanche” et victime de “race noire“: 0 condamné à mort, 16 non condamnés à mort
  • meurtrier de “race noire“et victime de “race blanche“: 11 condamnés à mort, 37 non condamnés à mort
  • meurtrier de “race noire“et victime de “race noire“: 4 condamnés à mort, 139 non condamnés à mort

On peut alors faire des tests d’indépendance, entre la “race” du meurtrier et le verdict par exemple.

MEURTRIER=matrix(c(53+0,11+4,414+16,139+37),2,2)
VICTIME  =matrix(c(53+11,0+4,414+37,139+16),2,2)
n=sum(MEURTRIER)
(PROBMEUTR=MEURTRIER/n)
           [,1]      [,2]
[1,] 0.07863501 0.6379822
[2,] 0.02225519 0.2611276

SL=rowSums(PROBMEUTR)
SC=colSums(PROBMEUTR)
(MEUTRINDEP=outer(SL, SC, "*"))
           [,1]      [,2]
[1,] 0.07229966 0.6443176
[2,] 0.02859055 0.2547922

(Q=n*sum((PROBMEUTR - MEUTRINDEP)^2/MEUTRINDEP))
[1] 1.468519

(Qcorrec=n*sum((abs(PROBMEUTR - MEUTRINDEP)-.5/n)^2/MEUTRINDEP))
[1] 1.144741

pchisq(Qcorrec, (2-1)*(2-1), lower.tail
 = FALSE)
[1] 0.2846528

qchisq(.95, (2-1)*(2-1))
[1] 3.841459

chisq.test(MEURTRIER)

Pearson's Chi-squared test with Yates' continuity correction

data:  MEURTRIER 
X-squared = 1.1447, df = 1, p-value = 0.2847

On rejette donc l’hypothèse d’indépendance.

et c’est parti pour les élections…

Ce soir, @imparibus reprenait  une information du Point (ici). “Un sondage de l’institut Sensus diffusé jeudi crédite Mme Rousseff, 62 ans, de 46,8% des intentions de vote et José Serra, 68 ans, de 42,7% […]  la marge d’erreur du sondage est de 2,2%, a expliqué un reponsable de l’enquête [… ce qui donne une ] grande probabilité que les 2 candidats soient à égalité”. @imparibus affirmait que c’était faux,@adelaigue me demandant mon avis…
Bon, je suis flâté qu’Alexandre m’appelle à l’aide, mais je suis plutôt nul en sondage… Bref, pour faire une réponse intelligente et détaillée, il me faudra un peu plus de place que ce que ne me propose Twitter…
Soit https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial01.png la probabilité de voter pour R, et https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial02.png la probabilité de voter pour S. On note https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial03.png la probabilité de voter pour un autre candidat, https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial04.png.
Si https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial05.png et https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial06.png sont les fréquences empiriques obtenues à partir de https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial15.png personnes sondées, rappelons que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial07.png
https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial08.png

et

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial09.png

Aussi,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial11.png

Si on veut tester https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial12.png (hypothèse https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial13.png) car c’est finalement de ça qu’on parle ici, l’idée la plus naturelle serait d’introduire la statistique de test d’égalité des proportions,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial14.png

qui suit – à la louche – une loi normale sous https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial13.png.
C’est bien joli tout ça, mais je ne connais pas https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial15.png ! En fait, on va essayer de s’aider de la petite phrase “la marge d’erreur du sondage est de 2,2%”
Alors si j’en crois le site d’ipsos (ici),  cette marge d’erreur correspond à 2000 personnes interrogées…
Or ici https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial05.png =46% et et https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial06.png = 42.7%, avec

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial16.png

soit,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial18.png

Si on compare à une loi normale, on a 2.6% de chances d’atteindre une telle valeur…
On peut aussi faire un test du chi-deux. En cas d’égalité des proportions, on devrait avoir https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial19.png=44.75%
et si on utilise un test du chi-deux pour mesurer la distance entre les deux lois multinomiale (sous https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial13.png et sous H1), on a comme distance du chi-deux

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial21.png

soit ici https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial22.png=3.71 avec 2000 personnes interrogées. Sous https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial13.pnghttps://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial22.png devrait suivre une loi du chi-deux à 1 degré de liberté, ce qui donne une p-value de 5.26%.
On peut aussi tenter un test de rapport de vraisemblance. Sous https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial13.png, on a

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial24.png

alors que sous https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial23.png

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial25.png

(en utilisant la forme des lois multinomiales). Le rapport des vraisemblances s’écrit alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial26.png

soit ici https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial27.png= 10.436. La statistique de test est icihttps://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial27.png qui doit suivre, sous https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial13.png une loi du chi-deux à un degré de liberté. La p-value associée est alors de l’ordre de 0,1%….
Bon, il semble qu’on puisse raffiner un peu les tests, en reprenant les idées de Nass (1959) et de Hoel (1938) – mais qui sont intéressant si on a peu de données (ce qui n’est pas vraiment le cas ici).
Bref, j’aurais tendance à être plutôt d’accord avec @imparibus. Maintenant, je suis nul en sondage, donc si quelqu’un a d’autres idées, je suis preneur ! Surtout que j’ai peut-être fait quelques erreurs de calculs ici ou là….