r - 计算 BradleyTerry 能力排名

Question

我正在尝试使用 R 中的 BradleyTerry2 包输出 Bradley Terry 模型的能力估计值。我不断收到一个非常神秘的错误。他们的文档工作中的一些示例和其他示例返回了我使用数据时遇到的相同错误。此代码使用文档中的示例分析之一。因此，如果您加载库，“变色龙”数据应该已经存在

install.packages("BradleyTerry2")
library (BradleyTerry2)
summary(chameleon.model <- BTm(player1 = winner, player2 = loser,formula = ~ prev.wins.2 + ch.res[ID] + prop.main[ID] + (1|ID), id = "ID",data = chameleons))
BTabilities(chameleon.model)

我得到的错误是

Error in X[, est, drop = FALSE] : (subscript) logical subscript too long

有人知道怎么做吗？

Answer 1

我维护 BradleyTerry2。当能力公式中有特定于比赛的预测变量时，就会出现此错误。它应该按照 ?BTability 中记录的方式工作:

... the abilities are computed from the terms of the fitted model that involve player covariates only (those indexed by ‘model$id’ in the model formula). Thus parameters in any other terms are assumed to be zero.

我们不知道这不起作用，所以感谢您的错误报告。在修复之前，可以直接算出能力和标准误:

## get names of relevant coefficients
> nm <- grep("[ID]", names(coef(chameleon.model)),
>            fixed = TRUE, value = TRUE)
> nm
[1] "ch.res[ID]"    "prop.main[ID]"

## get names of corresponding predictors
> IDvar <- gsub("[ID]", "", nm, fixed = TRUE)
> IDvar
[1] "ch.res"    "prop.main"

## create coefficient vector and design matrix
> cf <- coef(chameleon.model)[nm]
> X <- as.matrix(chameleons$predictors[, IDvar])

## compute abilities
> abilities <- X %*% cf
> colnames(abilities) <- "abilities"

## compute standard errors
> V <- vcov(chameleon.model)[nm, nm]
> res <- cbind(abilities = abilities,
>               se = sqrt(diag(X %*% V %*% t(X))))
> head(res)
    abilities       se
C01  3.757229 1.655205
C02  2.404778 1.017782
C03  2.319346 1.133959
C04  1.892671 1.399391
C05  2.253472 1.101628
C06  2.015840 1.075806

对于预测变量中缺少值的个体，这将为 NA 提供 NA，但他们的能力是单独建模的，并在模型摘要中返回。上面的代码还假设连续协变量，如果您的某些预测变量是因素，您将需要更正式地创建一个模型矩阵，例如

X <- model.matrix(reformulate(IDvar, intercept = FALSE),
                  c(chameleons$predictors))