Jak dopasować model Bradley – Terry – Luce w R, bez skomplikowanej formuły?

9

Model Bradleya-Terry'ego-Luce (BTL) stwierdza, że , gdzie to prawdopodobieństwo, że obiekt zostanie oceniony jako „lepszy”, cięższe itp., niż obiekt oraz i są parametrami.pjotja=losoljat-1(δjot-δja)pjajotjotjaδjaδjot

Wydaje się, że jest to kandydat na funkcję glm, z rodziną = dwumianową. Jednak wzór byłby podobny do „Sukces ~ S1 + S2 + S3 + S4 + ...”, gdzie Sn jest zmienną fikcyjną, to jest 1, jeśli obiekt n jest pierwszym obiektem w porównaniu, -1 jeśli jest drugi i 0 w przeciwnym razie. Wówczas współczynnik Sn będzie odpowiadał .remiltzan

Byłoby to dość łatwe do zarządzania przy użyciu tylko kilku obiektów, ale mogłoby to prowadzić do bardzo długiej formuły i potrzeby utworzenia zmiennej zastępczej dla każdego obiektu. Zastanawiam się tylko, czy istnieje prostsza metoda. Załóżmy, że nazwa lub liczba dwóch porównywanych obiektów to zmienne (czynniki?) Obiekt1 i Obiekt2, a Sukces to 1, jeśli obiekt 1 zostanie oceniony lepiej, a 0, jeśli obiekt 2 to.

Silverfish
źródło
3
Istnieje pakiet R dla modelu Bradley-Terry. Spójrz na Rseek.
kardynał
Udostępniłem również kilka linków do powiązanego pytania: stats.stackexchange.com/a/10741/930
chl.
Wspomniany pakiet @cardinal, btw: BradleyTerry2
conjugateprior

Odpowiedzi:

17

Myślę, że najlepszym pakietem dla danych porównywania w parze (PC) w R jest pakiet prefmod , który pozwala wygodnie przygotowywać dane do dopasowania (log liniowy) modeli BTL w R. Używa Poissona GLM (dokładniej wielomianowego logitu w Poissonie sformułowanie patrz np. ta dyskusja ).

Zaletą jest to, że ma funkcję, prefmod::llbt.designktóra automatycznie konwertuje dane do niezbędnego formatu i niezbędnej matrycy projektowej.

Załóżmy na przykład, że masz 6 obiektów porównanych parami. Następnie

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

zbuduje macierz projektową z macierzy danych, która wygląda następująco:

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

z wierszami oznaczającymi ludzi, kolumnami oznaczającymi porównania, a 0 oznacza niezdecydowany 1 oznacza preferowany obiekt 1, a 2 oznacza preferowany obiekt 2. Brakujące wartości są dozwolone. Edycja : Ponieważ prawdopodobnie nie jest to coś, co można wywnioskować na podstawie powyższych danych, przeliterowałem to tutaj. Porównania należy uporządkować w następujący sposób ((12) oznacza średni obiekt porównawczy 1 z obiektem 2):

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

Dopasowanie jest najwygodniej przeprowadzane za pomocą gnm::gnmfunkcji, ponieważ pozwala na modelowanie statystyczne. (Edycja: Można również użyć prefmod::llbt.fitfunkcji, która jest nieco prostsza, ponieważ wymaga tylko zliczeń i macierzy projektu.)

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

Należy pamiętać, że termin eliminacji pomija parametry uciążliwe w podsumowaniu. Następnie możesz uzyskać wartościowe parametry (swoje delty) jako

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

I możesz je wykreślić

R> plotworth(wmat)

Jeśli masz wiele obiektów i chcesz szybko napisać obiekt formuły o1+o2+...+on, możesz użyć

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

wygenerować formułę gnm( dla której nie byłbyś potrzebny llbt.fit).

Istnieje artykuł JSS , patrz także https://r-forge.r-project.org/projects/prefmod/ i dokumentacja przez ?llbt.design.

Momo
źródło
1
To bardzo dokładna odpowiedź. Dziękuję Ci. Wygląda na to, że prefmod byłby dobrym pakietem do użycia. Nawiasem mówiąc, jestem zainteresowany wykorzystaniem modelu do próby przewidzenia wyników meczów sportowych.
Silverfish,
Nie ma problemu, cieszę się, że to pomogło. Nie wiem dokładnie, co masz zamiar przewidzieć, ale Leitner i in. wykorzystali te modele do przewidywania wydarzeń sportowych. Zobacz jego pracę dyplomową epubdev.wu.ac.at/2925 . Powodzenia.
Momo
Może ten link jest lepszy epubdev.wu.ac.at/view/creators/…
Momo
Czy na podstawie tych danych można obliczyć istotność różnic między poszczególnymi parami (np. O1 i o2)? A może musisz zmienić formułę, użyć o2 jako ostatniego czynnika i żyć w takim przypadku bez szacowania błędu standardowego?
TNT,
1
Minęło trochę czasu, więc nie pamiętam, czy możesz wygodnie użyć ograniczeń liniowych, ale w twoim przypadku możesz użyć jednego jako poziomu odniesienia, powiedzmy o1, i użyj wartości t drugiego, powiedzmy o2, z podsumowania - skutecznie stanowi test, czy różnica między o1 a o2 wynosi zero.
Momo,