R/flatlizards.R

In BradleyTerry2: Bradley-Terry Models

#' Augrabies Male Flat Lizards: Contest Results and Predictor Variables
#' 
#' Data collected at Augrabies Falls National Park (South Africa) in
#' September-October 2002, on the contest performance and background attributes
#' of 77 male flat lizards (*Platysaurus broadleyi*).  The results of
#' exactly 100 contests were recorded, along with various measurements made on
#' each lizard.  Full details of the study are in Whiting et al. (2006).
#' 
#' There were no duplicate contests (no pair of lizards was seen fighting more
#' than once), and there were no tied contests (the result of each contest was
#' clear).
#' 
#' The variables `head.length`, `head.width`, `head.height` and
#' `condition` were all computed as residuals (of directly measured head
#' length, head width, head height and body mass index, respectively) from
#' simple least-squares regressions on `SVL`.
#' 
#' Values of some predictors are missing (`NA`) for some lizards,
#' \sQuote{at random}, because of instrument problems unconnected with the
#' value of the measurement being made.
#' 
#' @name flatlizards
#' @docType data
#' @format This dataset is a list containing two data frames:
#' `flatlizards$contests` and `flatlizards$predictors`.
#' 
#' The `flatlizards$contests` data frame has 100 observations on the
#' following 2 variables: \describe{ 
#' \item{winner}{a factor with 77
#' levels `lizard003` ... `lizard189`.} 
#' \item{loser}{a factor
#' with the same 77 levels `lizard003` ... `lizard189`.} }
#' 
#' The `flatlizards$predictors` data frame has 77 observations (one for
#' each of the 77 lizards) on the following 18 variables: \describe{
#' \item{id}{factor with 77 levels (3 5 6 ... 189), the lizard
#' identifiers.} 
#' \item{throat.PC1}{numeric, the first principal
#' component of the throat spectrum.} 
#' \item{throat.PC2}{numeric, the
#' second principal component of the throat spectrum.}
#' \item{throat.PC3}{numeric, the third principal component of the
#' throat spectrum.} 
#' \item{frontleg.PC1}{numeric, the first principal
#' component of the front-leg spectrum.} 
#' \item{frontleg.PC2}{numeric,
#' the second principal component of the front-leg spectrum.}
#' \item{frontleg.PC3}{numeric, the third principal component of the
#' front-leg spectrum.} 
#' \item{badge.PC1}{numeric, the first principal
#' component of the ventral colour patch spectrum.}
#' \item{badge.PC2}{numeric, the second principal component of the
#' ventral colour patch spectrum.} 
#' \item{badge.PC3}{numeric, the third
#' principal component of the ventral colour patch spectrum.}
#' \item{badge.size}{numeric, a measure of the area of the ventral
#' colour patch.} 
#' \item{testosterone}{numeric, a measure of blood
#' testosterone concentration.} 
#' \item{SVL}{numeric, the snout-vent
#' length of the lizard.} 
#' \item{head.length}{numeric, head length.}
#' \item{head.width}{numeric, head width.}
#' \item{head.height}{numeric, head height.}
#' \item{condition}{numeric, a measure of body condition.}
#' \item{repro.tactic}{a factor indicating reproductive tactic; levels
#' are `resident` and `floater`.} }
#' @seealso [BTm()]
#' @references Turner, H. and Firth, D. (2012) Bradley-Terry models in R: The
#' BradleyTerry2 package.  *Journal of Statistical Software*,
#' **48**(9), 1--21.
#' 
#' Whiting, M. J., Stuart-Fox, D. M., O'Connor, D., Firth, D., Bennett, N. C.
#' and Blomberg, S. P. (2006). Ultraviolet signals ultra-aggression in a
#' lizard. *Animal Behaviour* **72**, 353--363.
#' @source The data were collected by Dr Martin Whiting,
#' http://whitinglab.com/?page_id=3380, and they appear here with his
#' kind permission.
#' @keywords datasets
#' @examples
#' 
#' ##
#' ##  Fit the standard Bradley-Terry model, using the bias-reduced
#' ##  maximum likelihood method:
#' ##
#' result <- rep(1, nrow(flatlizards$contests))
#' BTmodel <- BTm(result, winner, loser, br = TRUE, data = flatlizards$contests)
#' summary(BTmodel)
#' ##
#' ##  That's fairly useless, though, because of the rather small
#' ##  amount of data on each lizard.  And really the scientific
#' ##  interest is not in the abilities of these particular 77
#' ##  lizards, but in the relationship between ability and the
#' ##  measured predictor variables.
#' ##
#' ##  So next fit (by maximum likelihood) a "structured" B-T model in
#' ##  which abilities are determined by a linear predictor.
#' ##
#' ##  This reproduces results reported in Table 1 of Whiting et al. (2006):
#' ##
#' Whiting.model <- BTm(result, winner, loser, ~ throat.PC1[..] + throat.PC3[..] +
#'                      head.length[..] + SVL[..],
#'                      data = flatlizards)
#' summary(Whiting.model)
#' ##
#' ##  Equivalently, fit the same model using glmmPQL:
#' ##
#' Whiting.model <- BTm(result, winner, loser, ~ throat.PC1[..] + throat.PC3[..] +
#'                      head.length[..] + SVL[..] + (1|..), sigma = 0,
#'                      sigma.fixed = TRUE, data = flatlizards)
#' summary(Whiting.model)
#' ##
#' ##  But that analysis assumes that the linear predictor formula for
#' ##  abilities is _perfect_, i.e., that there is no error in the linear
#' ##  predictor.  This will always be unrealistic.
#' ##
#' ##  So now fit the same predictor but with a normally distributed error
#' ##  term --- a generalized linear mixed model --- by using the BTm
#' ##  function instead of glm.
#' ##
#' Whiting.model2 <- BTm(result, winner, loser, ~ throat.PC1[..] + throat.PC3[..] +
#'                       head.length[..] + SVL[..] + (1|..), 
#'                       data = flatlizards, trace = TRUE)
#' summary(Whiting.model2)
#' ##
#' ##  The estimated coefficients (of throat.PC1, throat.PC3,
#' ##  head.length and SVL are not changed substantially by
#' ##  the recognition of an error term in the model; but the estimated
#' ##  standard errors are larger, as expected.  The main conclusions from
#' ##  Whiting et al. (2006) are unaffected.
#' ##
#' ##  With the normally distributed random error included, it is perhaps
#' ##  at least as natural to use probit rather than logit as the link
#' ##  function:
#' ##
#' require(stats)
#' Whiting.model3 <- BTm(result, winner, loser, ~ throat.PC1[..] + throat.PC3[..] +
#'                       head.length[..] + SVL[..] + (1|..),
#'                       family = binomial(link = "probit"),
#'                       data = flatlizards, trace = TRUE)
#' summary(Whiting.model3)
#' BTabilities(Whiting.model3)
#' ##  Note the "separate" attribute here, identifying two lizards with
#' ##  missing values of at least one predictor variable 
#' ##
#' ##  Modulo the usual scale change between logit and probit, the results
#' ##  are (as expected) very similar to Whiting.model2.
#' 
"flatlizards"