Comment adapter le modèle Bradley – Terry – Luce en R, sans formule compliquée?

9

Le modèle Bradley – Terry – Luce (BTL) indique que , où est la probabilité que l'objet soit jugé "meilleur", plus lourd, etc., que l'objet , et et sont des paramètres.pjje=logjet-1(δj-δje)pjejjjeδjeδj

Cela semble être un candidat pour la fonction glm, avec family = binomial. Cependant, la formule serait quelque chose comme "Success ~ S1 + S2 + S3 + S4 + ...", où Sn est une variable fictive, c'est-à-dire 1 si l'objet n est le premier objet de la comparaison, -1 s'il l'est le second et 0 sinon. Alors le coefficient de Sn serait le correspondant .eltunen

Cela serait assez facile à gérer avec seulement quelques objets, mais pourrait conduire à une formule très longue et à la nécessité de créer une variable fictive pour chaque objet. Je me demande simplement s'il existe une méthode plus simple. Supposons que le nom ou le numéro des deux objets comparés sont des variables (facteurs?) Object1 et Object2, et Success est 1 si l'objet 1 est mieux jugé et 0 si l'objet 2 l'est.

Silverfish
la source
3
Il existe un package R pour le modèle Bradley-Terry. Regardez Rseek.
Cardinal
J'ai également fourni quelques liens sur une question connexe: stats.stackexchange.com/a/10741/930
chl
Le package @cardinal mentionné, btw: BradleyTerry2
conjugateprior

Réponses:

17

Je pense que le meilleur package pour les données de comparaison par paires (PC) dans R est le package prefmod , qui permet de préparer facilement les données pour s'adapter aux modèles BTL (log linéaires) dans R. Il utilise un GLM de Poisson (plus précisément, un logit multinomial dans Poisson voir par exemple cette discussion ).

La bonne chose est qu'il a une fonction prefmod::llbt.designqui convertit automatiquement vos données dans le format et la matrice de conception nécessaires.

Par exemple, supposons que vous ayez 6 objets tous par paires comparés. alors

R> library(prefmod)
R> des<-llbt.design(data, nitems=6)

va construire la matrice de conception à partir d'une matrice de données qui ressemble à ceci:

P1  0  0 NA  2  2  2  0  0  1   0   0   0   1   0   1   1   2
P2  0  0 NA  0  2  2  0  2  2   2   0   2   2   0   2   1   1
P3  1  0 NA  0  0  2  0  0  1   0   0   0   1   0   1   1   2
P4  0  0 NA  0  2  0  0  0  0   0   0   0   0   0   2   1   1
P5  0  0 NA  2  2  2  2  2  2   0   0   0   0   0   2   2   2
P6  2  2 NA  0  0  0  2  2  2   2   0   0   0   0   2   1   2

avec des lignes indiquant des personnes, des colonnes indiquant des comparaisons et 0 signifie indécis 1 signifie objet 1 préféré et 2 signifie objet 2 préféré. Les valeurs manquantes sont autorisées. Edit : Comme ce n'est probablement pas quelque chose à déduire simplement des données ci-dessus, je l'explique ici. Les comparaisons doivent être ordonnées de la manière suivante ((12) moyenne objet de comparaison 1 avec objet 2):

(12) (13) (23) (14) (24) (34) (15) (25) etc. 

L'ajustement est plus pratique que la gnm::gnmfonction, car il vous permet de faire une modélisation statistique. (Modifier: vous pouvez également utiliser la prefmod::llbt.fitfonction, qui est un peu plus simple car elle ne prend que les nombres et la matrice de conception.)

R> res<-gnm(y~o1+o2+o3+o4+o5+o6, eliminate=mu, family=poisson, data=des)
R> summary(res)
  Call:
gnm(formula = y ~ o1 + o2 + o3 + o4 + o5 + o6, eliminate = mu, 
    family = poisson, data = des)

Deviance Residuals: 
   Min      1Q  Median      3Q     Max  
-7.669  -4.484  -2.234   4.625  10.353  

Coefficients of interest:
   Estimate Std. Error z value Pr(>|z|)    
o1  1.05368    0.04665  22.586  < 2e-16 ***
o2  0.52833    0.04360  12.118  < 2e-16 ***
o3  0.13888    0.04297   3.232  0.00123 ** 
o4  0.24185    0.04238   5.707 1.15e-08 ***
o5  0.10699    0.04245   2.521  0.01171 *  
o6  0.00000         NA      NA       NA    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1 

(Dispersion parameter for poisson family taken to be 1)

Std. Error is NA where coefficient has been constrained or is unidentified

Residual deviance: 2212.7 on 70 degrees of freedom
AIC: 2735.3

Veuillez noter que le terme éliminer supprimera les paramètres de nuisance du résumé. Vous pouvez alors obtenir les paramètres de valeur (vos deltas) comme

## calculating and plotting worth parameters
R> wmat<-llbt.worth(res)
        worth
o1 0.50518407
o2 0.17666128
o3 0.08107183
o4 0.09961109
o5 0.07606193
o6 0.06140979

Et vous pouvez les tracer avec

R> plotworth(wmat)

Si vous avez de nombreux objets et que vous souhaitez écrire o1+o2+...+onrapidement un objet formule , vous pouvez utiliser

R> n<-30
R> objnam<-paste("o",1:n,sep="")
R> fmla<-as.formula(paste("y~",paste(objnam, collapse= "+")))
R> fmla
y ~ o1 + o2 + o3 + o4 + o5 + o6 + o7 + o8 + o9 + o10 + o11 + 
    o12 + o13 + o14 + o15 + o16 + o17 + o18 + o19 + o20 + o21 + 
    o22 + o23 + o24 + o25 + o26 + o27 + o28 + o29 + o30

pour générer la formule gnm(dont vous n'auriez pas besoin llbt.fit).

Il y a un article JSS , voir aussi https://r-forge.r-project.org/projects/prefmod/ et la documentation via ?llbt.design.

Momo
la source
1
C'est une réponse très approfondie. Je vous remercie. Il semble que prefmod soit un bon package à utiliser. Je suis intéressé par l'utilisation du modèle pour essayer de prédire les résultats des matchs de sport, soit dit en passant.
Silverfish
Pas de problème, content si cela a aidé. Je ne sais pas exactement comment vous entendez prédire, mais Leitner et al. ont utilisé ces modèles pour prédire les événements sportifs. Voir sa thèse epubdev.wu.ac.at/2925 . Bonne chance.
Momo
Peut-être que ce lien est meilleur epubdev.wu.ac.at/view/creators/…
Momo
Est-il possible de calculer la signification des différences entre les paires individuelles (par exemple o1 et o2) à partir de ces données? Ou devez-vous réorganiser la formule, utiliser o2 comme dernier facteur et vivre sans estimation d'erreur standard dans ce cas?
TNT
1
Cela fait un moment, donc je ne me souviens pas si vous pouvez facilement utiliser des restrictions linéaires mais ce que vous pouvez faire dans votre cas est d'utiliser l'un comme niveau de référence, par exemple o1, et d'utiliser la valeur t de l'autre, par exemple o2, du résumé - il constitue effectivement un test pour savoir si la différence entre o1 et o2 est nulle.
Momo