2.3 Ranking Responses with Bradley-Terry Models
2.3.2 Minorizing functions and MM algorithm
The strict concavity of the logarithm function implies for positive x and y that
− ln x ≥ 1 − ln y − (x/y) with equality if and only if x = y (16)
Therefore, as shown in Lange, Hunter and Yang (2000), if we fix γ(k)and define the function
Qk(γ) =
we may conclude that
Qk(γ) ≤ `(γ) with equality if γ = γ(k) (18)
where `(γ) is the log-likehood of (14). A function Qk(γ) satisfying conditions (18) is said minorize `(γ) at the point γ(k). It is easy to verify that, for any Qk(γ) satisfying the minorizing conditions (18),
Qk(γ) ≥ Qk(γ(k)) implies `(γ) ≥ `(γ(k)) (19)
Property (19) suggests an iterative algorithm in which we let γ(k) denote the value of the parameter vector before the k th iteration and define γ(k+1) to be the maximizer of Qk(γ);
thus γ(k+1) of (15) maximizes Qk(γ). Since this algorithm consists of alternately creating a minorizing function Qk(γ) and then maximizing it, Hunter and Lange (2000) call it an MM algorithm.
3 Ranking Rule
In this section we introduce a criterion to rank the responses using the methods intro-duced in Section 2. The ranking rule for the Wald test, the Generalized Score test and the Bayesian ranking methods are similar, but they are different from the rank method with Bradley-Terry model.
Now we first illustrate the ranking rule of the Wald test, the Generalized Score test and the Bayesian ranking methods. Assume that we have k responses and corresponding mj
value, j = 1, . . . , k. Let m(j) be the order statistics, that is, m(1) ≤ m(2) ≤, . . . , m(k). Let v(j) be the response corresponding to m(j). It is natural to rank the importance of responses in order of m(j). That is, the most important response is v(k), and the second important response is v(k−1). The proposed testing methods in Section 2.1 and Section 2.2 can be used to rank the responses. If the hypothesis π(k) = π(k−1) is rejected, we may claim that v(k) is the most important response. If it is not rejected, the two responses have same rank, that is, the response v(k) is as important as the response v(k−1), and next to test the hypothesis π(k−1) = π(k−2). Similarly, if it is rejected, we rank response v(k−1) first and response v(k−2) second. If it is not rejected, response v(k−1) and v(k−2) have same rank.
We use Question 1 as an example to illustrate the rule. For example, let m1=54, m2=49, m3=28, m4=71, m5=23, and then we have m(1)=23, m(2)=28, m(3)=49, m(4)=54, m(5)=71.
It is natural to rank the importance of responses in order of m(j). We may claim that price is more important than taste for consumers to purchase. However, this rank method based on the order of m(j) is not statistically significant. Hence, we follow the proposed rule and use one of these methods to rank all response. First, we rank response v(5) and response v(4), i.e. rank the response ”taste” and response ”price”. If H01 : π(5) = π(4) is not reject, it means that the response ”taste” and the response ”price” are equally important. Then we test H02: π(4) = π(3) versus H12: π(4) 6= π(3). If H02 is rejected ,we rank responses v(5)andv(4)
first and response v(3) third. And when hypothesis H03 : π(3) = π(2) and H04 : π(2) = π(1), we reject H03 and do not reject H04. Then we denote the ranking notations for the above result as ”taste” 1, ”capacity” 3, ”packaging” 4, ”price” 1, ”other” 4. Hence, we know that the response ”taste” and the response ”price” are top priorities when consumers purchase a drink, the response ”capacity” is third important factor, and the response ”packaging” and
the response ”other” are relatively unimportant for consumers.
Next, we illustrate the ranking rule of the method with Bradley-Terry model. It is easier than above rule. In this method, we compute all γ values and according to the size of values, we can obtain an order with descending. Hence, the order is the rank of these responses.
4 Simulation
In this section, a simulation study is conducted to evaluate the performance of these methods in this section. Because the Bayesian ranking method has prior distribution as-sumption, it is different from other methods. Hence, we do not discuss Bayesian rank-ing method in the simulation study. The true rank of these responses is according to the order of π(j). In this study, we regard two responses have the same rank if |π(i) − π(j)| ≤ where is a constant in a tolerance region. We compare the ranks of the 5 methods in terms of consistent rate, which is defined as the proportion that the rank of these methods and the true rank are consistent for n respondents in 1000 replicates.
For example, let π10000 = 0.032, π01000 = 0.015, π00100 = 0.087, π00010 = 0.061, π00001 = 0.009, π11000 = 0.008, π10100= 0.0082, π10010= 0.068, π10001= 0.002, π01100= 0.002, π01010= 0.00005, π01001= 0.0005, π00110= 0.0005, π00101 = 0.00006, π00011= 0.00319 ,others equal to 0.044 and = 0.01, resulting π1 = 0.626, π2 = 0.501, π3 = 0.585, π4 = 0.617, π5 = 0.479.
Thus the true is 1 4 3 1 5. We obtain a sample, and use the Wald test to obtain the rank 1 3 3 1 5. If the rank of each response derived from the the Wald test is smaller
than the true rank, we call this phenomenon is consistent. Here, the result of the ex-ample is consistent. The codes to run γ in the Bradley-Terry models can be download in http://sites.stat.psu.edu/~dhunter/code/btmatlab/. We apply these programs to rank responses in a multiple response question. Since there are three codes for this method in http://sites.stat.psu.edu/~dhunter/code/btmatlab/. The first code of using Bradly-Terry with MM method is denoted as btmm in Tables 1-6 and in the R code section. The second code of using Bradly-Terry with quasi-Newton accelerated MM method is denoted as btqn in Tables 1-6 and in the R code section. The third code of using Bradly-Terry with a Newton-Raphson method is denoted as btnr in Tables 1-6 and in the R code section.
Then the following table is the consistent rate for different k and n:
Table 1: The consistent rates of the 5 methods when π1 = 0.77, π2 = 0.28, π3 = 0.56, π4 = 0.21, π5 = 0.33, k =5 and =0.05 for all methods and α=0.05 for the Wald test and the Generalized Score test:
Sample size
Method
Wald test G.S. test btmm btqn btnr
n=100 0.997 0.996 0.682 0.682 0.682
n=200 0.999 0.999 0.814 0.814 0.814
n=300 0.999 0.999 0.879 0.879 0.879
n=500 1 1 0.951 0.951 0.951
n=800 1 1 0.968 0.968 0.968
n=1000 1 1 0.992 0.992 0.992
Table 2: The consistent rates of the 5 methods when π1 = 0.77, π2 = 0.28, π3 = 0.56, π4 = 0.21, π5 = 0.34, π6 = 0.43, k =6 and =0.05 for all methods and α=0.05 for the Wald test and the Generalized Score test:
Sample size
Method
Wald test G.S. test btmm btqn btnr
n=100 0.998 0.998 0.578 0.578 0.578
n=200 0.999 0.999 0.796 0.796 0.796
n=300 1 1 0.869 0.869 0.869
n=500 1 1 0.95 0.95 0.95
n=800 1 1 0.99 0.99 0.99
n=1000 1 1 0.992 0.992 0.992
Table 3: The consistent rates of the 5 methods when π1 = 0.77, π2 = 0.28, π3 = 0.56, π4 = 0.21, π5 = 0.34, π6 = 0.43, π7 = 0.12, k =7 and =0.05 for all methods and α=0.05 for the Wald test and the Generalized Score test:
Sample size
Method
Wald test G.S. test btmm btqn btnr
n=100 0.999 0.999 0.528 0.528 0.528
n=200 0.999 0.999 0.772 0.772 0.772
n=300 1 1 0.865 0.865 0.865
n=500 1 1 0.946 0.946 0.946
n=800 1 1 0.979 0.979 0.979
n=1000 1 1 0.99 0.99 0.99
Table 4: The consistent rates of the 5 methods when π1 = 0.77, π2 = 0.28, π3 = 0.56, π4 = 0.21, π5 = 0.34, π6 = 0.43, π7 = 0.12, π8 = 0.5, k =8 and =0.05 for all methods and α=0.05 for the Wald test and the Generalized Score test:
Sample size
Method
Wald test G.S. test btmm btqn btnr
n=100 0.999 0.998 0.355 0.355 0.355
n=200 0.999 0.999 0.626 0.626 0.626
n=300 1 1 0.796 0.796 0.796
n=500 1 1 0.91 0.91 0.91
n=800 1 1 0.972 0.972 0.972
n=1000 1 1 0.985 0.985 0.985
Table 5: The consistent rates of the 5 methods when π1 = 0.77, π2 = 0.28, π3 = 0.56, π4 = 0.21, π5 = 0.34, π6 = 0.43, π7 = 0.12, π8 = 0.5, π9 = 0.9, π10 = 0.62, k =10 and =0.05 for all methods and α=0.05 for the Wald test and the Generalized Score test:
Sample size
Method
Wald test G.S. test btmm btqn btnr
n=100 0.998 0.998 0.251 0.251 0.251
n=200 0.999 0.999 0.524 0.524 0.524
n=300 1 1 0.684 0.684 0.684
n=500 1 1 0.879 0.879 0.879
n=800 1 1 0.957 0.957 0.957
n=1000 1 1 0.985 0.985 0.985
Next, we compare the Wald test and the Generalized Score test for different α. Then the following tables are the consistent rates of the Wald test and the Generalized Score test for different α:
Table 6: The consistent rates of the Wald test and the Generalized Score test when π1 = 0.77, π2 = 0.28, π3 = 0.56, π4 = 0.21, π5 = 0.34, k =5 and =0.05:
Sample size
n=100 n=200
Significant level α
Method
Wald test G.S. test Wald test G.S. test
α=0.15 0.984 0.984 0.996 0.996
α=0.1 0.99 0.991 0.997 0.997
α=0.05 0.997 0.996 0.999 0.999
α=0.01 0.999 0.999 1 1
Table 7: The consistent rates of the Wald test and the Generalized Score test when π1 = 0.77, π2 = 0.28, π3 = 0.56, π4 = 0.21, π5 = 0.34, π6 = 0.43, π7 = 0.12, k =7 and =0.05:
Sample size
n=100 n=200
Significant level α
Method
Wald test G.S. test Wald test G.S. test
α=0.15 0.992 0.992 0.997 0.997
α=0.1 0.994 0.994 0.997 0.997
α=0.05 0.999 0.999 0.999 0.999
α=0.01 1 1 1 1
Table 8: The consistent rates of the Wald test and the Generalized Score test when π1 = 0.77, π2 = 0.28, π3 = 0.56, π4 = 0.21, π5 = 0.34, π6 = 0.43, π7 = 0.12, π8 = 0.5, π9 = 0.9, π10 = 0.62, k =10 and =0.05:
Sample size
n=100 n=200
Significant level α
Method
Wald test G.S. test Wald test G.S. test
α=0.15 0.992 0.992 0.998 0.998
α=0.1 0.997 0.998 0.999 0.999
α=0.05 0.998 0.998 0.999 0.999
α=0.01 1 1 1 1
According to above results, we find that consistent rate decreases as the number of responses increase for the methods with Bradley-Terry model when n is not enough large.
When the sample size is large, the results of these methods are almost consistent. In com-paring the Wald test and the Generalized Score test for different α, consistent rate increases when α decreases. Although the consistent rate is not high for small sample size case, it still has good result for large sample size case. It reveals that these methods are feasible in ranking responses when the sample size is not small.
5 R code
These ranking procedures has been written as a package RankResponse for R. RankRe-sponse is available from the Comprehensive R Archive Network at http://CRAN.R-project.
org/package=RankResponse, which include code function rank.wald, rank.gs, rank.LR, rank.LN, rank.L2R, rank.btmm, rank.btqn and rank.btnr.
rank.wald Rank responses based on the Wald test
Description
Rank responses of a single response question or a multiple response question by the Wald test procedure.
Usage
rank.wald(data,alpha,type=2) Argument
data A m × n matrix (dij), where dij = 0 or 1. If the i th respondent selects the j th response, then dij = 1, otherwise dij = 0.
alpha The significance level used in the Wald test.
type type=1 for a single response question;
type=2 for a multiple response question.
Value
The rank.wald returns the estimated probabilities of the responses being selected and the ranks of the responses by the Wald test procedure.
References
Wang, H. (2008). Ranking Responses in Multiple-Choice Questions. Journal of Applied Statistics, 35, 465-474.
Examples
## This is an example to rank three responses in a multiple response
## question when the number of respondents is 1000 and the
signifi-## cance level is 0.05.In this example,we do not use a real data, but
## generate data in the first three lines.
A <-sample(c(0,1),1000,p=c(0.21,0.79),replace=T) B <-sample(c(0,1),1000,p=c(0.86,0.14),replace=T) C <-sample(c(0,1),1000,p=c(0.42,0.58),replace=T) D <-cbind(A,B,C)
data <-matrix(D,nrow=1000,ncol=3)
# or upload the true data alpha<-0.05
rank.wald(data,alpha,type=2)
rank.gs Rank responses based on the Generalized score test
Description
Rank responses of a single response question or a multiple response question by the generalized score test procedure.
Usage
rank.gs(data,alpha,type=2) Argument
data A m × n matrix (dij), where dij = 0 or 1. If the ith respondent selects the jth response, then dij = 1, otherwise dij = 0.
alpha The significance level used in the Generalized score test.
type type=1 for a single response question ; type=2 for a multiple response question .
Value
The rank.gs returns the estimated probabilities of the responses being selected and the ranks of the responses by the Generalized score procedure.
References
Wang, H. (2008). Ranking Responses in Multiple-Choice Questions. Journal of Applied Statistics, 35, 465-474.
Examples
## This is an example to rank three responses in a multiple response
## question when the number of respondents is 1000 and the
signifi-## cance level is 0.05.In this example,we do not use a real data, but
## generate data in the first three lines.
A <-sample(c(0,1),1000,p=c(0.21,0.79),replace=T) B <-sample(c(0,1),1000,p=c(0.86,0.14),replace=T) C <-sample(c(0,1),1000,p=c(0.42,0.58),replace=T) D <-cbind(A,B,C)
data <-matrix(D,nrow=1000,ncol=3)
# or upload the true data alpha<-0.05
rank.gs(data,alpha,type=2)
rank.btmm Rank responses based on the Bradley-Terry model with the MM method
Description
Adopt the Bradley-Terry model to rank responses in a single response question or in a multiple response question with the MM method. This method associates each response with a value γ, and use the γ value to rank responses.
Usage
rank.btmm(data) Argument
data A m × n matrix (dij), where dij = 0 or 1. If the i th respondent selects the j th response, then dij = 1, otherwise dij = 0.
Value
The rank.btmm returns the associated γ values in the first line and the ranks ofthe responses in the second line.
References
Hunter DR (2004). MM algorithms for generalized Bradley-Terry models. The Annals of Statistics, 32, 384-406.
Examples
## This is an example to rank three responses in a multiple response
## question when the number of respondents is 1000. In this example,
## we do not use a real data, but generate data in the first three lines.
A <-sample(c(0,1),1000,p=c(0.37,0.63),replace=T)
B <-sample(c(0,1),1000,p=c(0.71,0.29),replace=T) C <-sample(c(0,1),1000,p=c(0.22,0.78),replace=T) D <-cbind(A,B,C)
data <-matrix(D,nrow=1000,ncol=3)
# or upload the true data rank.btmm(data)
rank.btqn Rank responses based on the Bradley-Terry model with the quasi-Newton accelerated MM method
Description
Adopt the Bradley-Terry model to rank responses in a single response question orin a multiple response question with quasi-Newton and the MM method. This method associates each response with a value γ, and use the γ value to rank responses.
Usage
rank.btqn(data) Argument
data A m × n matrix (dij), where dij = 0 or 1. If the i th respondent selects the j th response, then dij = 1, otherwise dij = 0.
Value
The rank.btqn returns the associated γ values in the first line and the ranks of the re-sponses in the second line.
References
Hunter DR (2004). MM algorithms for generalized Bradley-Terry models. The Annals of Statistics, 32, 384-406.
Examples
## This is an example to rank three responses in a multiple response
## question when the number of respondents is 1000. In this example,
## we do not use a real data, but generate data in the first three lines.
A <-sample(c(0,1),1000,p=c(0.37,0.63),replace=T) B <-sample(c(0,1),1000,p=c(0.71,0.29),replace=T) C <-sample(c(0,1),1000,p=c(0.22,0.78),replace=T) D <-cbind(A,B,C)
data <-matrix(D,nrow=1000,ncol=3)
# or upload the true data rank.btqn(data)
rank.btnr Rank responses based on the Bradley-Terry model with New-ton Raphson method
Description
Adopt the Bradley-Terry model to rank responses in a single response question or in a multiple response question with Newton-Raphson method. This method associates each response with a value γ, and use the γ value to rank responses.
Usage
rank.btnr(data)
Argument
data A m × n matrix (dij), where dij = 0 or 1. If the i th respondent selects the j th response, then dij = 1, otherwise dij = 0.
Value
The rank.btnr returns the associated γ values in the first line and the ranks of the re-sponses in the second line.
References
Hunter DR (2004). MM algorithms for generalized Bradley-Terry models. The Annals of Statistics, 32, 384-406.
Examples
## This is an example to rank three responses in a multiple response
## question when the number of respondents is 1000. In this example,
## we do not use a real data, but generate data in the first three lines.
A <-sample(c(0,1),1000,p=c(0.37,0.63),replace=T) B <-sample(c(0,1),1000,p=c(0.71,0.29),replace=T) C <-sample(c(0,1),1000,p=c(0.22,0.78),replace=T) D <-cbind(A,B,C)
data <-matrix(D,nrow=1000,ncol=3)
# or upload the true data rank.btnr(data)
rank.LN Rank responses under the Bayesian framework according to the loss function LN(d, n) = cF D + F N
Description
Rank responses of a single response question or a multiple response question under the Bayesian framework according to the loss function LN(d, n) = cF D + F N .
Usage
LN(data,response.number,prior.parameter,c) Argument
data A m × n matrix (dij), where dij = 0 or 1. If the i th res-pondent selects the j th response, then dij = 1, otherwise dij = 0.
response.number The number of the responses
prior.parameter The parameter vector of the Dirichlet prior distribution, where the vector dimension is 2response.number.
c The value of c in the loss function Value
The rank.LN returns the estimated probabilities of the responses being selected in the first line and the ranks of the responses in the second line.
References
Wang, H. and Huang, W. H. (2014). Bayesian Ranking Responses in Multiple Response Questions. Journal of the Royal Statistical Society: Series A (Statistics in Society), 177, 191-208.
Examples
##This is an example to rank three responses in a multiple response
##question when the number of respondents is 1000 and the value c is
##1. In this example, we do not use a real data, but generate data in
##the first three lines.
A <-sample(c(0,1),1000,p=c(0.37,0.63),replace=T) B <-sample(c(0,1),1000,p=c(0.71,0.29),replace=T) C <-sample(c(0,1),1000,p=c(0.22,0.78),replace=T) D <-cbind(A,B,C)
data <-matrix(D,nrow=1000,ncol=3)
# or upload the true data response.number <-3
prior.parameter <- c(5,98,63,7,42,7,7,7) c <-1
rank.LN(data,response.number,prior.parameter,c)
rank.LR Rank responses under the Bayesian framework according to the loss LR(d, n) = cF DR + F N R
Description
Rank responses of a single response question or a multiple response question under the Bayesian framework according to the loss function LR(d, n) = cF DR + F N R.
Usage
rank.LR(data,response.number,prior.parameter,c) Argument
data A m × n matrix (dij), where dij = 0 or 1. If the i th res-pondent selects the j th response, then dij = 1, otherwise dij = 0.
response.number The number of the responses
prior.parameter The parameter vector of the Dirichlet prior distribution, where the vector dimension is 2response.number.
c The value of c in the loss function Value
The rank.LR returns the estimated probabilities of the responses being selected in the first line and the ranks of the responses in the second line.
References
Wang, H. and Huang, W. H. (2014). Bayesian Ranking Responses in Multiple Response Questions. Journal of the Royal Statistical Society: Series A (Statistics in Society), 177, 191-208.
Examples
##This is an example to rank three responses in a multiple response
##question when the number of respondents is 1000 and the value c is
##0.33. In this example, we do not use a real data, but generate data
##in the first three lines.
A <-sample(c(0,1),1000,p=c(0.37,0.63),replace=T) B <-sample(c(0,1),1000,p=c(0.71,0.29),replace=T) C <-sample(c(0,1),1000,p=c(0.22,0.78),replace=T) D <-cbind(A,B,C)
data <-matrix(D,nrow=1000,ncol=3)
# or upload the true data response.number <-3
prior.parameter <- c(5,98,63,7,42,7,7,7) c <-0.33
rank.LR(data,response.number,prior.parameter,c)
rank.L2R Rank responses under the Bayesian framework according to the loss L2R(d, n) = (F DR, F N R)
Description
Rank responses of a single response question or a multiple response question under the Bayesian framework according to the loss function L2R(d, n) = (F DR, F N R).
Usage
rank.L2R(data,response.number,prior.parameter,e) Argument
data A m × n matrix (dij), where dij = 0 or 1. If the i th res-pondent selects the j th response, then dij = 1, otherwise dij = 0.
response.number The number of the responses
prior.parameter The parameter vector of the Dirichlet prior distribution, where the vector dimension is 2response.number.
e A cut point used in the loss function which depends on the economic costs.
Value
The rank.L2R returns the estimated probabilities of the responses being selected in the first line and the ranks of the responses in the second line.
References
Wang, H. and Huang, W. H. (2014). Bayesian Ranking Responses in Multiple Response Questions. Journal of the Royal Statistical Society: Series A (Statistics in Society), 177, 191-208.
Examples
##This is an example to rank three responses in a multiple response
##question when the number of respondents is 1000 and the value e is
##0.15. In this example, we do not use a real data, but generate data
##in the first three lines.
A <-sample(c(0,1),1000,p=c(0.37,0.63),replace=T) B <-sample(c(0,1),1000,p=c(0.71,0.29),replace=T) C <-sample(c(0,1),1000,p=c(0.22,0.78),replace=T) D <-cbind(A,B,C)
data <-matrix(D,nrow=1000,ncol=3)
# or upload the true data response.number <-3
prior.parameter <- c(5,98,63,7,42,7,7,7) e <-0.15
rank.L2R(data,response.number,prior.parameter,e)
6 Conclusion
In this thesis, our goal is ranking responses of a single response question or a multiple response question. We proposed some methods to solve this problem. For ranking responses, the simulation results in Section 4 show that the proposed methods have good performance.
Although these methods are not consistent for small sample size, their performances are very good for large sample size. In real applications, the sample size of responses is usually not very small, these methods are feasible in apply to real applications. According to the simu-lation results, the consistent rates of the Wald test and the generalized score test are larger than the method with Bradley-Terry model. Hence, we conclude that the Wald test and the generalized score test are more powerful than the method with Bradley-Terry model on ranking responses. We do not discuss Bayesian ranking responses method in the simulation study, because the other methods are under frequentist setup, but the Bayesian method is under the Bayesian framework. The codes of these methods has been written as an R package such that it is more convenient for readers to use them.
7 References
[1] Agresti, A. and Liu, I.-M. (1999) Modeling a categorical variable allowing arbitrarily many category choices. Biometrics 55, 936-943.
[2] Benjamini, Y. and Hochberg, Y. (1995) Controlling the false discovery rates: a practical and powerful approach to multiple testing. J. R. Statist. Soc. B, 57, 289–300.
[3] Berger, J. O. (1985) Statistical Decision Theory and Bayesian Analysis. 2nd edn. New York: Springer.
[4] Bilder, C. R., Loughin, T. M.and Nettleton, D. (2000) Multiple marginal independence testing for pick any/c variables. Communications in Statistics-Simulation and Compu-tation, 29(4), 1285-1316.
[5] Bradley, R. A. and Terry, M. E. (1952). Rank analysis of incomplete block designs. I.
Themethod of paired comparisons. Biometrika 39, 324–345.
[6] Decady, Y. J. and Thomas, D. H. (2000). A simple test of association for contingency tables with multiple column responses. Biometrics 56, 893-896.
[7] Hunter DR (2004), MM algorithms for generalized Bradley-Terry models, Annals of Statistics , 32, 386-408
[8] J.-L. Ke (2011) Algorithms for Ranking Responses in Multiple Response Questions. Na-tional Chiao Tung University
[9] Lange, K., Hunter, D.R., Yang I. (2000) Optimization transfer using surrogate objective functions (with discussion). J. Comput. Graphical Stat. 9: 1-59.
[10] Loughin, T. M. and Scherer, P. N. (1998). Testing for association in contingency tables with multiple column responses. Biometrics 54, 630-637.
[11] Umesh, U. N. (1995). Predicting nominal variable relationships with multiple responses.
Journal of Forecasting 14, 585-596.
[12] Wang, H. (2008). Ranking responses in multiple-choice questions. Journal of Applied Statistics, 35, 465-474
[13] Wang, H. and Huang, W. H. (2014). Bayesian Ranking Responses in Multiple Response Questions. Journal of the Royal Statistical Society, Series A (Statistics in Society), 177, 191-208