A Bayesian Look at the Rao-Kupper Model for Paired Comparisons

The method of paired comparisons may be regarded as a special rank order technique. It is a method long used in psychological experimentation and is well adapted to sensory difference testing. This study deals with estimation and testing of the parameters of the Rao-Kupper (1967) model for paired comparisons which allows for tied observations in paired comparison experiments using informative prior. An elicitation of the hyper parameters is carried out. Predictive probabilities that one treatment would be preferred to other treatment in a future single comparison are determined. Graphs for the marginal posterior densities are presented. All calculations are based on numerical integration technique and Gibbs sampling procedure. Appropriateness of the model is checked which declares it a good fit.


Introduction
In the method of paired comparisons (PCs), judges are presented with pairs of items and, for each pair, they are asked to choose the preferred item according to some criterion.Though these methods had their origin in psychophysics but were gradually incorporated into broader applications in experimental psychology.The value of paired comparison lies in the fact that it produces result that are statistically robust and conceptually meaningful.The modern Bayesian approach appeared as an alternative to the frequntist's approach in 1960s and1970s, and has grown into an increasingly substantial and effective methodology.The Bayesian paradigm allows us to incorporate prior information into statistical models for decision making.One feature of the certain type of paired comparison experiment is that prior knowledge may be available about the merits of the objects.Prior distribution which gives us specific and definite information about the variable is known as an informative prior.Paired comparison models due to availability of prior information have been studied by many authors through Bayesian approach.
A vast amount of literature exists on the subject.Aslam (2002) conducts the Bayesian analysis of the Bradley-Terry model and the Rao-Kupper model for paired comparisons.Aslam (2003) presents a method based on the prior predictive distribution to elicit the hyper parameters for the parameters of the Rao-Kupper model for paired comparisons.Altaf et.al. (2013) consider Bayesian analysis of the Van Baaren model using the informative and conjugate priors.Glickman and Jenson (2005) describe a method for designing paired comparison experiments, particular to tournament scheduling that incorporate prior information in a principled manner.Altaf et.al. (2013) address the Davidson model for paired comparisons introducing an amendment to the model and perform Bayesian analysis on the said model using informative and non-informative priors.Haung et.al. (2006) propose a new exponential model which helps to compare and estimate the individual's abilities while playing in different competitions.Francis et.al. (2006) extended the log linear form of the Bradely-Terry model to incorporate subject specific co-variates.Tsai and Bockenholt (2006) introduce new framework which simultaneously relaxes both the independence and invariance assumptions assumed in paired comparisons models.
This study deals with Bayesian analysis of the Rao-Kupper model for paired comparisons.Section 2 of the article presents the Rao-Kupper model and likelihood function along with necessary notations used in analysis of the model.In Section 3, Prior distributions for the parameters of the models are suggested and prior predictive distribution which is used to elicit the hyper parameters is presented.Section 4 and its subsequent parts, fully covers the Bayesian analysis of the said model considering a 5 treatment case using informative prior.Section 5, appropriateness of the model is checked.Finally Section 6 presents conclusion of the study.

The Rao-Kupper Model
Consider a paired comparison experiment with m treatments or objects in which each of the ( 1) / 2 m m  distinct treatment pairs are ranked r times.The rankings in the 'r' repetitions of each pair can be assumed to be independent.
The Bradley-Terry (1952) introduced a basic paired comparison model for such an experiment.The m treatments have "true" treatment ratings The model defined in (1) is very well known Bradley-Terry model.
In a paired-comparison experiment, situation often arises in which a judge may not be able to express any real preference in a number of pairs he judges.It is, therefore, important in such cases when judge is unable to make a difference between two treatments, he should declare a tie between treatments rather to make a definite preference.It is clear that any model which does not allow for the possibility of tie is not making full use of the information contained in the no-preference class.
The probability that treatment j T is preferred to treatment i T is denoted by The equations ( 2 1or 0 ijk n  as according to treatment i T is tied with treatment j T or not.Also 0.
= the number of times treatment i T and treatment j T are tied.ij r is the number of times treatment i T is compared with treatment j T and The following notations are useful for further simplification of likelihood function.
T is preferred to treatment j T and the number of times treatments i T and j T are tied.
T is preferred to any other treatment.And the number of times treatment i T and j T are tied.
= the total number of times treatment i T and j T are tied.
The probability of the observed result in the k th repetition of the pair ( , ) Hence, the likelihood function of the observed outcome x {where x represents the data . ., 0. , Here  is the threshold parameter and 1 2 , ,..., m    are the treatment parameters.

The Choice of an Informative Prior
Developing prior distribution is undoubtedly the most controversial aspect of any Bayesian analysis (Lindley, 1983;Walters and Ludwig, 1994).Considerable care should be taken when selecting priors and the process by which priors are selected must be documented carefully.This is because inappropriate choices for priors can lead to incorrect inferences.It is often easier to develop a prior for each parameter in turn (marginal priors) than to develop a joint prior for all the parameters simultaneously.In fact, most of the Bayesians today have made the assumption that the priors for the various parameters are independent (e.g.Punt and Walker, 1998; Punt and Butterworth, 2000).
Assuming the treatment parameters 1 2 , ,..., m    and threshold parameter  to be independent, prior distribution for the treatment parameters is considered to be Dirichlet distribution and threshold parameter is assumed to have Gamma prior distribution.
The Dirichlet distribution, used as a prior distribution for the treatment parameters 1 2 , ,..., m    has the form: a stands for a generalization of the beta function.
The prior distribution of  is assumed to have a Gamma distribution with hyper parameters 1 2 , 0 b b  , so  has density: Both the priors are considered to be independent so the joint prior Dirichlet-Gamma distribution has the following form.( 1) ( 1) ( , , ,..., ) ( ,...., ) ( )

Assessment of the Prior Distribution
Practical implementation of Bayesian methods requires the assessment of prior distribution i.e. the elicitation of hyperparameters.The elicitation of prior belief has received a little attention despite the excellent review of Kadane and Wolfson (1998) which contains dozens of references.There exist many new and workable methods for assessment of prior distribution in statistical literature.We base our elicitation method on the fact given by Geisser (1980) and Kadane (1980) that the natural elements for statistical inference are characteristics of the predictive distribution of an analyst.

The Prior Predictive Distribution of the Model
Consider the case of m treatments.The prior predictive distribution of the sufficient statistic: 0. ., n n using the likelihood function (4) and joint prior distribution (5) applying constraint (1 ) where , , 1,..., In practice the integral is more efficiently computed as function of  and  where (1 ) The evaluation of ( 7) is intractable so the integral is evaluated numerically.

Bayesian Analysis of the Model using Informative Prior for m=5
To conduct the Bayesian analysis for five treatments the data given by Glenn and David (1960) is used.These data are taken for a paired comparison experiment, and is described for five treatments and 30 replications of each of 10 distinct treatment pairs.The data are presented in Table 1.The joint posterior distribution of the parameters 1, 2 3 4 5 , , ,      and  is:

Table 1: Data for Comparison of 5 Treatments Allowing Ties
) The marginal posterior density of the parameter 1  is : Expressions for the marginal posterior densities of 2 3 4 5 , , , and

Elicitation of Hyperparameters for m=5
In the method of elicitation suggested by Aslam (2003) being used, the predictive distribution is compared with the expert's assessment of this distribution and hyperparameters are chosen in such a way that make the assessment agree closely with a member of the family.He provides two methods to elicit the hyperparameters of the Rao-Kupper model

Method (a) Via Eliciting the Prior Predictive Probabilities Method (b) Via Eliciting the Confidence Levels
In method (a), elicitation of prior predictive probabilities is difficult for large number of comparisons, so the method (b) is adopted in which the prior predictive confidence levels (confidence level is a probability for a given interval) are elicited.The following function is used to elicit where 'K' denotes the number of confidence levels, FCL is for fitted confidence levels and ECL stands for the elicited confidence level.
Elicitation of hyperparameters for 5 treatments is based on the estimated values of confidence levels for the prior predictive distribution defined in (7).Using expression (7) hyperparameters are elicited via prior predictive confidence level.Minimum value of function ( 10) is sorted.Prior predictive probabilities are calculated for vector

Marginal Posterior Densities for the parameters of the Rao-Kupper Model
When m=5


The posterior means are determined for the parameters of the Rao-Kupper model using Drichlet-Gamma informative prior by numerical integration as well as by Gibbs sampling for the data given in Table 1.A program is designed in SAS package to find the posterior means.The obtained posterior means are presented in a Table 2. Posterior means are also determined by using uniform prior for comparison purpose.The posterior estimates obtained for both priors are in coincidence up to 1 decimal place and show same preference order of treatments as 3 The results obtained by Gibbs sampling are also close enough.
Let the parameters i  and j  , , 1, 2, 3, 4, 5 i j i j   are compared using hypotheses: The posterior probability for ( )  and the posterior probability for 1 Using the data given in Table 1, the posterior probabilities obtained for ij p and ij q by running a SAS program are given in Table :3 The predictive probability (12) P that treatment 1 T would be preferred to treatment 2 T in a future single comparison may be obtained by using posterior distribution (8) and the model probability (2) as: Where K is a normalizing constant and 5 Similarly the predictive probability (0.12) P that there would be a tie between two treatments 1 2 and T T in a future single comparison can be derived by using ( 8) and (3) as under: For the data given in Table 1, the predictive probabilities are obtained using the quadrature method with the help of a program written in SAS and are given in Table 3. Again the results show that predictive probabilities are completely in agreement with the previous ordering of treatments.Predictive probabilities obtained for no preference category turn out to be very small.

Appropriateness of the Model
The observed number of preferences with that of expected number of preferences are compared to check the goodness of fit of the model.A 2  Statistic is used for testing purpose, as it is widely used by renowned statisticians such as Glenn and David (1960) and Rao-Kupper (1967).
Let us consider two hypotheses.

( ) ( ) ( )
The null hypothesis is rejected when calculated   The observed and expected number of preferences is given in Table 4.The value of 2  Statistic is obtained to be 10.65 and the p-value at 15 degree of freedom for the obtained value of 2   is found to be 0.78, which indicates that model is worthy of fit.

Conclusion
The Bayesian analysis of the Rao-Kupper model for paired comparisons has been conducted for five treatments.Posterior means, which are used to rank the treatments, are obtained by two methods, quadrature method and the gibbs sampling technique.To compare the results posterior means are obtained using informative prior and a uniform prior, though only numerical expressions under informative prior are given.It is observed that estimates for treatments are little smaller under uniform prior than informative prior but they conclude to same ranking i.e.
) and (3) are known as the Rao-Kupper model.If  =1 then the Rao- Kupper model yields the Bradley-Terry model.The following notations are used for the analysis of the model.treatment i T is preferred to treatment j T or not, in the k th repetition of comparison.0.
  and K is normalizing constant:

0H: 1 H
The model is true for some values 0 : The model is not true for any value of  Let us consider .ˆi ij n = the expected number of times treatment i T is preferred to treatment j T and 0. ˆij n = the expected number of times treatments i T and j T are tied.Thus 2 Statistic degrees of freedom (which is the level of significance at which the null hypothesis is rejected).The expected number of the preferences are obtained from the following expressions.

.
 .Same treatment ranking is obtained through graphical display of marginal posterior densities and hypothesis testing procedure.
, for the given intervals.Those values for hyperparameters are selected against which ( )   is minimum.A program designed in SAS package is used to get these hyperparameters.It is found that ( )   =0.075 is minimum against the vector  = (4.98,4.98, 4.98, 4.98, 1.49, 1.55, 5.48).So these set of values are the elicited hyper parameters.
k FCL k 