Wie passt das Bradley-Terry-Luce-Modell ohne komplizierte Formel in R?

9

Das Bradley-Terry-Luce (BTL) -Modell besagt, dass , wobei die Wahrscheinlichkeit ist, dass das Objekt als "besser" beurteilt wird. schwerer usw. als Objekt und und sind Parameter.pjich=lÖGicht- -1(δj- -δich)pichjjichδichδj

Dies scheint ein Kandidat für die glm-Funktion mit family = binomial zu sein. Die Formel wäre jedoch so etwas wie "Erfolg ~ S1 + S2 + S3 + S4 + ...", wobei Sn eine Dummy-Variable ist, dh 1, wenn Objekt n das erste Objekt im Vergleich ist, -1, wenn dies der Fall ist die zweite und 0 sonst. Dann wäre der Koeffizient von Sn das entsprechende .delteinn

Dies wäre mit nur wenigen Objekten recht einfach zu handhaben, könnte jedoch zu einer sehr langen Formel und der Notwendigkeit führen, für jedes Objekt eine Dummy-Variable zu erstellen. Ich frage mich nur, ob es eine einfachere Methode gibt. Angenommen, der Name oder die Nummer der beiden verglichenen Objekte sind Variablen (Faktoren?) Objekt1 und Objekt2, und Erfolg ist 1, wenn Objekt 1 besser beurteilt wird, und 0, wenn Objekt 2 ist.

Silberfisch
quelle
3
Für das Bradley-Terry-Modell gibt es ein R-Paket. Schau auf Rseek.
Kardinal
Ich habe auch einige Links zu einer verwandten Frage bereitgestellt
chl
Das Paket @cardinal erwähnt, übrigens: BradleyTerry2
Conjugateprior

Antworten:

17

Ich denke, das beste Paket für PC-Daten (Paired Comparison) in R ist das Prefmod-Paket , mit dem Daten bequem für (logarithmisch lineare) BTL-Modelle in R vorbereitet werden können. Es verwendet ein Poisson GLM (genauer gesagt ein multinomiales Logit in Poisson) Formulierung siehe zB diese Diskussion ).

Das Schöne ist, dass es eine Funktion hat prefmod::llbt.design, die Ihre Daten automatisch in das erforderliche Format und die erforderliche Designmatrix konvertiert.

Angenommen, Sie haben 6 Objekte, die alle paarweise verglichen werden. Dann

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

erstellt die Entwurfsmatrix aus einer Datenmatrix, die folgendermaßen aussieht:

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

mit Zeilen, die Personen bezeichnen, Spalten, die Vergleiche bezeichnen, und 0 bedeutet unentschlossen 1 bedeutet Objekt 1 bevorzugt und 2 bedeutet Objekt 2 bevorzugt. Fehlende Werte sind zulässig. Bearbeiten : Da dies wahrscheinlich nicht einfach aus den obigen Daten abzuleiten ist, schreibe ich es hier. Die Vergleiche müssen folgendermaßen angeordnet werden ((12) bedeutet Vergleichsobjekt 1 mit Objekt 2):

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

Die Anpassung erfolgt am bequemsten mit der gnm::gnmFunktion, da Sie damit statistische Modelle erstellen können. (Bearbeiten: Sie können auch die prefmod::llbt.fitFunktion verwenden, die etwas einfacher ist, da nur die Anzahl und die Entwurfsmatrix berücksichtigt werden.)

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

Bitte beachten Sie, dass der Eliminierungsbegriff die Störparameter in der Zusammenfassung weglässt. Sie können dann die Wertparameter (Ihre Deltas) als erhalten

## 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

Und Sie können sie mit zeichnen

R> plotworth(wmat)

Wenn Sie viele Objekte haben und o1+o2+...+onschnell ein Formelobjekt schreiben möchten , können Sie verwenden

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

um die Formel zu generieren gnm( für die Sie nicht benötigen würden llbt.fit).

Es gibt einen JSS-Artikel , siehe auch https://r-forge.r-project.org/projects/prefmod/ und die Dokumentation über ?llbt.design.

Momo
quelle
1
Das ist eine sehr gründliche Antwort. Vielen Dank. Es scheint, als wäre prefmod ein gutes Paket. Ich bin übrigens daran interessiert, mit dem Modell zu versuchen, die Ergebnisse von Sportspielen vorherzusagen.
Silberfischchen
Kein Problem, froh, wenn es geholfen hat. Ich weiß nicht genau, wie Sie vorhersagen wollen, aber Leitner et al. haben diese Modelle verwendet, um Sportereignisse vorherzusagen. Siehe seine Dissertation epubdev.wu.ac.at/2925 . Viel Glück.
Momo
Vielleicht ist dieser Link besser epubdev.wu.ac.at/view/creators/…
Momo
Ist es möglich, aus diesen Daten Signifikanzen für die Unterschiede zwischen einzelnen Paaren (zB o1 und o2) zu berechnen? Oder müssen Sie die Formel neu ordnen, o2 als letzten Faktor verwenden und in diesem Fall ohne Std.error-Schätzung leben?
TNT
1
Es ist eine Weile her, daher erinnere ich mich nicht, ob Sie bequem lineare Einschränkungen verwenden können, aber was Sie in Ihrem Fall tun können, ist, eine als Referenzpegel zu verwenden, sagen wir o1, und den t-Wert der anderen zu verwenden, sagen wir o2, aus der Zusammenfassung - es stellt effektiv einen Test dar, ob die Differenz zwischen o1 und o2 Null ist.
Momo