Utility in Time Description in Priority Best–Worst Discrete Choice Models: An Empirical Evaluation Using Flynn’s Data

Abstract

Discrete choice models (DCMs) are applied in many fields and in the statistical modelling of consumer behavior. This paper focuses on a form of choice experiment, best–worst scaling in discrete choice experiments (DCEs), and the transition probability of a choice of a consumer over time. The analysis was conducted by using simulated data (choice pairs) based on data from Flynn’s (2007) ‘Quality of Life Experiment’. Most of the traditional approaches assume the choice alternatives are mutually exclusive over time, which is a questionable assumption. We introduced a new copula-based model (CO-CUB) for the transition probability, which can handle the dependent structure of best–worst choices while applying a very practical constraint. We used a conditional logit model to calculate the utility at consecutive time points and spread it to future time points under dynamic programming. We suggest that the CO-CUB transition probability algorithm is a novel way to analyze and predict choices in future time points by expressing human choice behavior. The numerical results inform decision making, help formulate strategy and learning algorithms under dynamic utility in time for best–worst DCEs.

1. Introduction

Modelling choice behaviors has always been a very challenging part in statistical research. Best–worst scaling (BWs) choice pairs in model behaviors has been coined by Louviere and Woodworth (1990) [1]. With the addition of the attribute and attribute levels, the logistic distributions are convenient tools from analytic and practical aspects; they rely on assumptions of ordinality in the responses that can only be met under limited situations. Furthermore, the computational issues related to the large choices in the best–worst cases hamper patterns in the preferences. Discrete choice experiments (DCE) are one of popular models used to quantify the inspiration of attributes, which characterize the choice options. In the literature, Train (2009) [2] present multiple models based on different assumptions about the distribution of random components. In some of his suggested models, the error terms seem to be homogeneous and uncorrelated. Lancsar et al. (2013) [3] used a form of choice experiment, “Best–Worst Discrete Choice Experiments”, which is designed in a way by asking respondents to choose not only the best option in a choice set but also the worst option followed by the best from the remaining options and so on until the implied preference ordering of the option is obtained. Marley and Louviere (2005) [4] and Marley et al. (2008) [5] discuss some of the theoretical aspects of attribute-level best–worst experiments and Flynn et al. (2007) [6] provides detailed instructions for the analysis of such experiments. Following Marley et al. (2008) [5], Street and Knox (2012) [7] constructed a probability function called best–worst choice probability which can be used to find the probability of choice from one state to another. Working et al. (2018) [8] used a conditional logit approach to calculate the transition (of conditional choice) probabilities over time assuming that the consumer behavior is a time-dependent process. The assumption of “making two choices in two different time points is independent process for an individual” cannot be justified in general. As mentioned in Potoglou et al. (2011) [9], the interaction between attributes and choice pairs needs further piloting and investigation. The dependence in the modelling of the choices made over time is not fully included in many research venues. The dependence in the selection is a wide and fast field, but the theory is biased towards the estimation method. The premise of the traditional choice models assumes of mutually exclusive alternatives; that is, if a choice is made, alternates are not selected anymore. An alternative would be to build patterns of the preferences and their weighted measures. That option opens a new framework for choices, especially those that are made over time. Many factors shape the choices made. We include a priority in the selection over times and build it in the design. Specific association in the choices is highlighted to show the need to develop progressive models where temporal dependence in the choices is integrated, and utility is computed. In that sense, the innovation we propose is built under the copula theory, as it captures the dependence among the choices made over time periods irrespective of the marginal distributions. We utilized a completely new copula structure call ‘CO-CUB’ to calculate transition probabilities in the BWs choices over time.

The paper is organized as follows. In Section 2, DCMs are reviewed. The design for the BW choice pairs and their utilities are also included under a priority setting. The transition probabilities under CO-CUB are described in Section 3. We present the model’s theory overtime sequence to quantify and measure consumer behavior and derive utilities using the Markov decision process (MDPs) in Section 4. Section 5 is about the application. Empirical data from Flynn et al. (2007) [6] are simulated and our methods are tested, which have the ability to estimate selected significant parameters such as the “feeling” and “uncertainty” of human beings when making a choice. The literature on Copula for DCE is at its infancy by utilizing our newly constructed probability model to analyze the choices of individuals, comparisons are made, and we end up with a conclusion.

2. Preliminary Results for Best–Worst Discrete Choice Modelling

Best–worst scaling experiments are modified DCEs to elicit further information about the best and worst product, or best and worst attributes and attribute-levels of a product. Respondents are asked make choices between different attributes within a same product; that is the alternatives are the attributes. In the experiment, each of the alternatives is described by a set of $K$ attributes or characteristics of the product or scenario being modelled. Further, each attribute contains $l_k$ levels. That is, each product is represented by a profile $X=\left(A_1,A_2,\dots ,A_K\right),$ where $A_i$ represents an attribute level of a choice. Consider $K$ profiles such that $P_1=\left(A_{11},A_{12},\dots ,A_{1K}\right),P_2=\left(A_{21},A_{22},\dots ,A_{2K}\right),\dots ,P_G=\left(A_{G1},A_{G2},\dots ,A_{GK}\right)$. For every profile, the choice set is then, ${\mathrm{C}}_X=\left[\left(x_1,x_2\right),\dots ,\left(x_1,x_K\right),\left(x_2,x_3\right),\dots \left(x_{K-1},x_K\right),\left(x_1,x_2\right),\dots ,\left(x_K,x_{K-1}\right)\right],$ where the first attribute level is the best and the second is the worst, called best–worst pairs. This model is designed in such way that the choice pair $\left(x_a,x_b\right)$, $x_a$, and $x_b$ cannot be acquired from the same attribute. From the choice set ${\mathrm{C}}_X$, the respondent determines a BW choice pair from $\tau$ such pairs, where $\tau =K\left(K-1\right)$. If we consider $G$ products, that is, $G$ profiles and $G$ associated choice sets, then we can define a finite set $S$ that is $S=\left\{{\mathrm{C}}_1,{\mathrm{C}}_2,\dots ,{\mathrm{C}}_G\right\}$ such that ${\mathrm{C}}_i\ne {\mathrm{C}}_j$ for all $i\ne j$. The BW choice pairs are illustrated in Figure 1 and shows all the possible combinations.

Total number of attribute levels

$$L={\sum }_{k=1}^K{l}_k$$

and the total number of unique attribute-level pairs

$$Z={\sum }_{k=1}^K{l}_k\left(L-l_k\right),$$

(1)

where $l_k$ is the number of levels for $k^{th}$ attribute, where $1\le k\le K$. Equation (1) is also described in Street and Knox (2012) [7].

2.1. Design of Experiment

Assume that the profiles have $K$ attributes, and each attribute has $l_k$ levels, where $k=1,2,\dots ,K$. The simulation has performed in such a way that each person has to make choices in all sets. The debate concerning optimal designs for DCEs could be enhanced by a closer connection to the established literature on the optimal design of experiments as suggested in Sun et al. (2023) [10]. We used orthogonal design here for the simulation. There are three standard models in best–worst choice model design experiments which are: (1) paired, (2) marginal, and (3) marginal sequential models. In this manuscript, we used the paired model Aizaki and Fogarty (2019) [11] and Marley and Louviere (2005) [4], which assumes that the difference in utility between the two levels represents the greatest utility difference among all $\tau$ utility differences.

2.1.1. Paired Model with Attribute and Level Variables

Flynn et al. (2007) [6] used a paired model in his “Quality of life” analysis. In this approach, the attribute variables are created as dummy-coded variables (attribute-specific constants), with reversed signs when the attributes are treated as the worst and the level variables are created as effect coding with a base level for each attribute. The signs of the effect-coded level variables are also reversed when levels are treated as the worst. When estimating this model, an arbitrary attribute variable is omitted, and the coefficient of the omitted attribute variable is normalized to be zero. That is, the coefficients of the remaining attribute variables are estimated relative to the omitted attribute. Unlike dummy coding, effect coding allows the calculation of a coefficient of the base level in each attribute: it is the negative of the sum of the remaining coefficients in each attribute.

Let us consider the general scenario to the formulation for the paired model. Consider the $K$ attributes $A_1,A_2,\dots ,A_K$; each has $l_k$ levels. Assume that attribute $A_K$ is omitted, and that the last level $\left(l_k\right)$ in each attribute (i.e., levels $A_{1l_k},A_{2l_k},\dots$ and $A_{K{l}_k}$) is the base level. As in Aizaki and Fogarty (2019) [11], the systematic component of the utility function is then:

$$V={\beta }_{A_1}{M}_{A_1}+\cdots +{\beta }_{A_{K-1}}{M}_{A_{K-1}}+{\beta }_{A_{11}}{M}_{A_{11}}+\cdots +{\beta }_{A_{1\left(l_k-1\right)}}{M}_{A_{1\left(l_k-1\right)}}+\cdots +{\beta }_{A_{\mathrm{K}\left(l_k-1\right)}}{M}_{A_{\mathrm{K}\left(l_k-1\right),}}$$

(2)

where the $\beta$s are the coefficients to be estimated and $M$s are the design matrices associated with the attributes and attribute levels. The coefficient of the attribute variable corresponding to $A_K$ is zero, and the coefficients of the base level $A_{j{l}_k},j=1,2,\dots ,K$ can be calculated by $A_{j{l}_k}=-\sum _{i=1}^{k-1}{\beta }_{A_{ji}}$ for $j=1,2,\dots ,K$. Detailed calculations are shown in the application section. $M$ is composed of indicators for the best and worst attributes and attribute levels. Let $M$ be the $Z\times P$ design matrix, where $Z={\sum }_{k=1}^K{l}_k\left(L-l_k\right)$ and $P=K{l}_k-1$. Columns correspond to attribute $A_K$ and the base levels ($A_{1l_k},A_{2l_k},\dots$ and $A_{K{l}_k}$) are omitted from the design matrix. Consider the choice pair $\left(x_{ij},x_{i{j}^{\prime }}\right)$ from the choice set ${\mathrm{C}}_i,$ for $i=1,2,\dots ,G,j\ne j^{\prime }=1,2,\dots ,K,$ and $1\le x_{ij}\le l_k$. For

$$M_{A_h}=\{\begin{array}{c}1 , \mathrm{if} x_{ij} ϵ A_h \mathrm{for} h=1,2,\dots ,K-1 \\ -1,\mathrm{if} x_{i{j}^{\prime }} ϵ A_h \mathrm{for} h=1,2,\dots ,K-1 \\ 0, \mathrm{o}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{w}\mathrm{i}\mathrm{s}\mathrm{e}\end{array}$$

(3)

Let $M_{A_h{x}_{ih}}$ be the data for the attribute level $1\le x_{ih}\le l_h, h=1,2,\dots ,K-1$ (non-base levels) within attribute $A_J$, $\forall j=1,2,\dots ,K$. Referring to the choice pair $\left(x_{ij},x_{i{j}^{\prime }}\right),$ the corresponding data for the attribute-levels are given by,

$$M_{A_h{x}_{ih}}=\{\begin{array}{c}1 , \mathrm{if} x_{ij} ϵ A_h \mathrm{for} h=1,2,\dots ,K-1 \\ -1,\mathrm{if} x_{i{j}^{\prime }} ϵ A_h \mathrm{for} h=1,2,\dots ,K-1 \\ 0, \mathrm{o}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{w}\mathrm{i}\mathrm{s}\mathrm{e}\end{array}.$$

(4)

Within attribute $A_J$, $\forall j=1,2,\dots ,K$, the base levels ($A_{1l_k},A_{2l_k},\dots$ and $A_{K{l}_k}$)

$$M_{A_1{l}_k}=M_{A_2{l}_k}=\cdots =M_{A_K{l}_k}=\{\begin{array}{c}-1 , \mathrm{if} x_{ij} ϵ A_K \\ 1 ,\mathrm{if} x_{i{j}^{\prime }} ϵ A_K \\ 0, \mathrm{o}\mathrm{t}\mathrm{h}\mathrm{e}\mathrm{r}\mathrm{w}\mathrm{i}\mathrm{s}\mathrm{e}\end{array}.$$

(5)

Table 1. Design Matrix based on Equations (3)–(5).

B	W	${\mathit{M}}_{\mathit{A}}$	${\mathit{M}}_{\mathit{B}}$	${\mathit{M}}_{\mathit{C}}$	${\mathit{M}}_{\mathit{D}}$	${\mathit{M}}_{\mathit{A}1}$	${\mathit{M}}_{\mathit{A}2}$	${\mathit{M}}_{\mathit{A}3}$	${\mathit{M}}_{\mathit{B}1}$	${\mathit{M}}_{\mathit{B}2}$	${\mathit{M}}_{\mathit{B}3}$	${\mathit{M}}_{\mathit{C}1}$	${\mathit{M}}_{\mathit{C}2}$	${\mathit{M}}_{\mathit{C}3}$	${\mathit{M}}_{\mathit{D}1}$	${\mathit{M}}_{\mathit{D}2}$	${\mathit{M}}_{\mathit{D}3}$	${\mathit{M}}_{\mathit{E}1}$	${\mathit{M}}_{\mathit{E}2}$	${\mathit{M}}_{\mathit{E}3}$
${\mathit{A}}_{\mathbf{3}}$	${\mathit{D}}_{\mathbf{2}}$	$1$	$0$	$0$	$-1$	$0$	$0$	$1$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$-1$	$0$	$0$	$0$	$0$
${\mathit{A}}_{\mathbf{1}}$	${\mathit{B}}_{\mathbf{3}}$	$1$	$-1$	$0$	$0$	$1$	$0$	$0$	$0$	$0$	$-1$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$
${\mathit{E}}_{\mathbf{3}}$	${\mathit{C}}_{\mathbf{3}}$	$0$	$0$	$-1$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$-1$	$0$	$0$	$0$	$0$	$0$	$1$
${\mathit{A}}_{\mathbf{2}}$	${\mathit{C}}_{\mathbf{4}}$	$1$	$0$	$-1$	$0$	$0$	$1$	$0$	$0$	$0$	$0$	$1$	$1$	$1$	$0$	$0$	$0$	$0$	$0$	$0$
${\mathit{D}}_{\mathbf{2}}$	${\mathit{E}}_{\mathbf{3}}$	$0$	$0$	$0$	$1$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$1$	$0$	$0$	$0$	$-1$
${\mathit{A}}_{\mathbf{1}}$	${\mathit{B}}_{\mathbf{4}}$	$1$	$-1$	$0$	$0$	$1$	$0$	$0$	$1$	$1$	$1$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$
${\mathit{A}}_{\mathbf{4}}$	${\mathit{B}}_{\mathbf{2}}$	$1$	$-1$	$0$	$0$	$-1$	$-1$	$-1$	$0$	$-1$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$
${\mathit{A}}_{\mathbf{1}}$	${\mathit{E}}_{\mathbf{4}}$	$1$	$0$	$0$	$0$	$1$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$1$	$1$	$1$
$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$	$⋮$
${\mathit{E}}_{\mathbf{4}}$	${\mathit{A}}_{\mathbf{4}}$	$-1$	$0$	$0$	$0$	$1$	$1$	$1$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$0$	$-1$	$-1$	$-1$

shows the design matrix for the based on data from Flynn et al. (2007) [6] ‘Quality of Life experiment.’ The estimated parameters of the choices are given in

Table 2. Estimated parameters.

Parameters	${\widehat{\mathit{\beta }}}_{\mathit{F}\mathit{l}\mathit{y}\mathit{n}\mathit{n}}$	$\mathit{S}{\mathit{E}}_{\mathit{F}\mathit{l}\mathit{y}\mathit{n}\mathit{n}}$	${\widehat{\mathit{\beta }}}_{\mathit{S}\mathit{i}\mathit{m}}$	$\mathit{S}{\mathit{E}}_{\mathit{S}\mathit{i}\mathit{m}}$		Parameters	${\widehat{\mathit{\beta }}}_{\mathit{F}\mathit{l}\mathit{y}\mathit{n}\mathit{n}}$	$\mathit{S}{\mathit{E}}_{\mathit{F}\mathit{l}\mathit{y}\mathit{n}\mathit{n}}$	${\widehat{\mathit{\beta }}}_{\mathit{S}\mathit{i}\mathit{m}}$	$\mathit{S}{\mathit{E}}_{\mathit{S}\mathit{i}\mathit{m}}$
Constant	−0.3067	0.075	*	*		Security All	0.9511	*	0.8199	*
Attachment	0.8105	0.0803	0.8655	0.0575		Enjoyment None	−0.8888	0.1286	−0.8327	0.0759
Security	*	*	*	*		Enjoyment Little	−0.3367	0.1632	−0.3535	0.0773
Enjoyment	0.2632	0.101	0.2684	0.0563		Enjoyment Lot	0.6561	0.1493	0.5862	0.0781
Role	0.1908	0.0974	0.2294	0.0565		Enjoyment All	0.5695	*	0.60007	*
Control	0.1076	0.0971	0.1217	0.0564		Role None	−0.8956	0.1239	−0.7697	0.0761
Attachment_None	−1.9678	0.1129	−1.9633	0.0775		Role Few	−0.0277	0.1532	0.1086	0.0776
Attachment_Little	0.1694	0.1012	0.1811	0.0756		Role Many	0.4435	0.1363	0.2616	0.0784
Attachment_Lot	0.9053	0.0905	1.0493	0.0811		Role All	0.4798	*	0.3995
Attachment_All	0.8932	*	0.7329	*		Control None	−0.8085	0.1122	−0.7149	0.0763
Security_None	−0.6123	0.118	−0.6953	0.0764		Control Few	0.0835	0.1596	0.1505	0.0774
Security Little	−0.3761	0.1302	−0.2449	0.077		Control Many	0.278	0.1376	0.2854	0.0785
Security Lot	0.0373	0.1153	0.1203	0.0779		Control All	0.4471	*	0.2791	*

${\widehat{\beta }}_{Flynn}$—Estimated parameters from Flynn’s experiment. ${\widehat{\beta }}_{Sim}$—Estimated parameters from the simulation. The * represent the base levels in our simulated data.

. The BW choice pairs are then arranged according to their popularity (with highest number of picks) in

Table 3. Arrangement of BW choice pairs in choice set 1, based on their popularity.

Best Attribute	Worst Attribute	Best Level	Worst Level	No. of Picks	Rank
Attachment	Security	Attach_all	Security_none	2629	1
Attachment	Role	Attach_all	Role_few	1204	2
Enjoyment	Security	Enjoyment_all	Security_none	1075	3
Attachment	Control	Attach_all	Control_many	994	4
Control	Security	Control_many	Security_none	700	5
Attachment	Enjoyment	Attach_all	Enjoyment_all	629	6
Role	Security	Role_few	Security_none	509	7
Enjoyment	Role	Enjoyment_all	Role_few	498	8
Enjoyment	Control	Enjoyment_all	Control_many	400	9
Control	Role	Control_many	Role_few	341	10
Role	Control	Role_few	Control_many	206	11
Control	Enjoyment	Control_many	Enjoyment_all	162	12
Security	Role	Security_none	Role_few	130	13
Role	Enjoyment	Role_few	Enjoyment_all	122	14
Enjoyment	Attachment	Enjoyment_all	Attach_all	117	15
Security	Control	Security_none	Control_many	89	16
Control	Attachment	Control_many	Attach_all	62	17
Security	Enjoyment	Security_none	Enjoyment_all	58	18
Role	Attachment	Role_few	Attach_all	51	19
Security	Attachment	Security_none	Attach_all	24	20

.

2.2. Utility Function–Attribute–Level Best–Worst Design

The utility function as described in McFadden (1974) [12] for the $i^{th}$ consumer/individual selecting the $j^{th}$ choice is given as:

$$U_{ij}=V_{ij}+{\epsilon }_{ij}$$

(6)

where $U_{ij}$ is the utility for $i^{th}$ consumer selecting $j^{th}$ choice. $V_{ij}$ captures the systematic component and ${\epsilon }_{ij}$ captures the error component. The distribution of the error terms is studied in McFadden (1974) [12]. He proposed the type I extreme value distribution (Gumbel distribution) for the error terms, which modelled data from the conditional logit model. However, this utility function can only accommodate single choice cases. Under the random utility theory, the probability of an alternative is based on the utility as defined in Equation (6); Flynn et al. (2007, 2008) [6,13] and Louviere et al. (2008) [14] provided the utility function for best–worst choice models. Here, we consider the choice set $\mathrm{C}$, $\left(x_{ij},x_{i{j}^{\prime }}\right)$ being the chosen pair, and the utility for choosing this pair within set $\mathrm{C}$ is then given by:

$$U_{ij{j}^{\prime }}=V_{ij{j}^{\prime }}+{\epsilon }_{ij{j}^{\prime }},$$

(7)

where $U_{ij{j}^{\prime }}$ is the utility value for $i^{th}$ consumer selecting $j^{th}$ choice as the best and $j^{\prime th}$ as the worst. $V_{ij{j}^{\prime }}$ is the corresponding systematic component and ${\epsilon }_{ij{j}^{\prime }}$ is the error term. The systematic component in Equation (7) can be calculated under the above-described paired model (Equation (2)).

3. Transition Probability

In the context of choice models, we say that the choices are not independent of an irrelevant alternative (IIA) and they are also time dependent. Therefore, we may not achieve an accurate overall performance by not accounting for the relationship between the present and future decisions and their outcomes. Our interest is with discrete time finite horizon Markov Decision Processes (MDP); that is, $t=1,2,\dots T$ where $T$ is a fixed number of time periods. The rewards are maximized by the best sequential decisions over time, making MDPs a dynamic optimization tool as used in Blanchet (2016) [15] to identify the right choices of substitution behaviors of consumers. Let $s_t\in S$ be the states occupied at time $t,r_t\left(s_t\right)$ which is the reward associated with $s_t$, and $d_t \left(r_t, s_t\right)$ is the decision based on the possible rewards and states at time $t$. The decision process maps the movement from one state to another over time $t$ based on rewards received and an optimal decision set. As the decision process is Markovian, the transition probability to the next state, $S^{\prime }=S_{t+1}$ based solely on the decision made at the current state, $S=S_t$ is $P_{S{S}^{\prime }}=P\left[S_{t+1}=S^{\prime }|S_t=S\right]$, where $t=1,\dots ,T$. The conditional probability of state transition from $S$ to $S^{\prime }$ is defined as:

$$P\left(s_{t+1}|s_t\right)=P\left(s_{t+1}=s^{\prime }|s_t=s\right)=P_{s{s}^{\prime }}$$

(8)

The transition probability of state $s^{\prime }$ is not homogeneous. As described in Equation (8), it depends on the action taken in the state $S$. Factors such as ‘feeling’ and ‘uncertainty’ of an individual can be highly affected to this probability. Further, the human behavior has a tendency to select options that have immediate gratification or those with the nearest range rather than selecting alternatives further away in time. Those are some of challenging constraints in terms of building a transition probability function. Further, the dependent structure of adjacent time probabilities needs to be handled carefully. Based on McFadden (1974) [12], the conditional logit model on the error term, Working et al. (2018) [7] establish a time-dependent probability formula while avoiding the IIA condition. It is reasonable as an initiation, but the dependent structure of the choice at the next time point or the above-mentioned ‘close at hand’ constraint cannot be seen in that transition matrix. In this manuscript, we introduced copula-based distribution with a CUB marginal, which is a combination of a discreet uniform and shifted binomial distributions introduced by Piccolo (2003) [16]. This model is a totally new approach for discrete choice modelling.

3.1. Copula Methods

Many stochastic models define the relation between $X$ and $Y$ through expectations, $E\left(Y|X=x\right)=\alpha +\beta x+\epsilon$. Copulas are functions that connect multivariate distributions to their one-dimensional margins. If $F$ is an $m$-dimensional cumulative distribution function with one-dimensional margins, $F_1,\dots ,F_d,$ then there exists an $m$-dimensional copula $C$ such that $F\left(y_1,\dots ,y_d\right)=C\left(F_1\left(y_1\right),\dots ,F_d\left(y_d\right)\right)$. The case $d=2$ has attracted special attention; Trivedi and Zimmer (2005) [17]. Copula modelling demonstrates that practical implementation and estimation is relatively straightforward despite the complexity of its theoretical foundations. The copula approach for constructing joint distributions has gained popularity in recent years in applied research in finance, civil engineering, medicine, climate, and weather research and many more. Copula models have been widely used for modelling and characterizing dependence structure in multivariate data. The extension to handle discrete and continuous forms of data offers contingency to model the probability of best–worst pair selection probability.

3.2. CUB Probability Model

Piccolo (2003) [16], Piccolo and Simone (2019) [18], and D’Elia and Piccolo (2005) [19] introduce a radical model which is a Combination of a discrete Uniform and a shifted Binomial distribution (CUB). CUB jointly considers two latent components, called ‘feeling’ and ‘uncertainty’. General references for the statistical background for CUB are included in Iannario (2012) [20] and Iannario and Piccolo (2010) [21]. The stingy parameterization and the ease of estimation and interpretation make CUB models a very useful tool for choice data analyses. Piccolo and Simone (2019) [18] asserts that there is a possibility to fit statistical models to ordinal/choice data without the need of explanatory covariates. Our interest is to use CUB marginals for copula distribution to construct a transition probability matrix which can be used to identify human choice behavior in future time points, given the present state of choice. An interesting attempt at defining a bivariate CUB distribution using a multivariate model with fixed margins is made by Andreis and Ferrari (2013) [22]. The CUB is a class of mixture models, possibly involving covariates, developed as a new approach for modelling discrete choice processes. CUB models consider two components, which are the ‘uncertainty’ and ‘feeling’ of humans when making a choice. The intrinsic ‘uncertainty’ in choosing an item is modelled through a discrete uniform variable. The latent process leading to the choice is directed by the subjective ‘feeling’ and modelled using a shifted binomial distribution. The probability of observing a particular response $r=1,2, \dots m$, with $m$ known. Assuming $m$ > 3 to ensure identifiability, the mixture probability model is described in Piccolo (2003) [16] and D’Elia and Piccolo (2005) [19]. It is given as:

$$P\left(R=r|\pi ,\xi \right)=\pi \left(\begin{array}{c}m-1\\ r-1\end{array}\right){\xi }^{m-r}{\left(1-\xi \right)}^{r-1}+\left(1-\pi \right)\frac{1}{m }; r=1,2,\dots ,m,$$

(9)

with $\pi ϵ (0; 1]$ and $\xi ϵ (0; 1]$.

Here, $\pi$ defines the weights of discrete uniform shifted binomial and it is inversely related to the amount of uncertainty in the choice decisions. The $\xi$ is related to personal preferences and measures the strength of feeling. In CUB models, there is non-uniqueness of the copula representation and the non-uniqueness stems from the fact that marginal distribution functions, which are not strictly monotonically increasing, rather monotonically non-decreasing, do not possess an inverse in the usual sense, rather a pseudo-inverse; Nelsen (2007) [23], Genest et al. (2007) [24]. Because of that reason, we must be precisely cautioned when handling copulas with non-continuous margins, which is in the CUB model.

3.3. The CO-CUB Model (Copula-Based CUB)

Andreis and Ferrari (2013) [22] defines a multidimensional extension of CUB in Equation (9), called CO-CUB model, as a multivariate copula with discrete margins, each following a CUB distribution. A $k$-dimensional $\left(k \ge 2\right)$ CO-CUB model with copula $C$ is a multivariate discrete variable with margins $R_i ~CUB\left({\pi }_i,{\xi }_i\right), i=1,\dots ,k,$ each with support $\left\{1, \dots ,m_i\right\}, m_i>3$. Joint distribution function is given by:

$$\mathsf{\Psi }\left(r_1,\dots ,r_k,\underset{\_}{\pi },\underset{\_}{\xi },\underset{\_}{\alpha }\right)=P\left(R_1\le r_1,\dots ,R_k\le r_k,\underset{\_}{\pi },\underset{\_}{\xi },\underset{\_}{\alpha }\right)=C_{\underset{\_}{\alpha }}\left[F_1\left(r_1;{\pi }_1,{\xi }_1\right),\dots ,F_k\left(r_k;{\pi }_k,{\xi }_k\right)\right],$$

(10)

where $\underset{\_}{\pi }={\left({\pi }_1,\dots ,{\pi }_k\right)}^{\prime }$ and $\underset{\_}{\xi }={\left({\xi }_1,\dots ,{\xi }_k\right)}^{\prime }$ as in Equation (6), for a particular choice of copula $C$, characterized by a parameter $\underset{\_}{\alpha }={\left({\alpha }_1,\dots ,{\alpha }_d\right)}^{\prime }$ taking values in some real $d$-dimensional space $\mathsf{\alpha }$ defining the dependence structure of its components. $F_i\left(r_i;{\pi }_i,{\xi }_i\right)$ stands for the distribution function of the $i^{th}$ margin; that is, $F_i\left(r_i\right)=P\left(R_i\le r_i\right)$ and the support of the CO-CUB variable is the grid $\left\{1,\dots ,m_1\right\}\times \dots \times \left\{1,\dots ,m_k\right\}$. The whole parameter set for a $k$-dimensional CO-CUB is, then, the ordered triplet $\left(\underset{\_}{\pi },\underset{\_}{\xi },\underset{\_}{\alpha }\right)ϵ$ ${\left[0,1\right]}^k\times {\left[0,1\right]}^k$ $\times \mathsf{\alpha }$, having the following interpretation: by the definition of copula, being margins of CUB, parameters $\left(\underset{\_}{\pi },\underset{\_}{\xi }\right)$ hold the same interpretation as in the unidimensional case, while for what concerns the copula parameter $\underset{\_}{\alpha }$, its interpretation as a dependence measure relates to the specific copula $C$ adopted in Equation (10). Andreis and Ferrari (2013) [22] justified that the Placket copula distribution is one of the better candidates for the CO-CUB model.

3.4. The CO-CUB Model with Plackett Copula without Covariates

The Plackett copula is defined as:

$$\begin{gathered} C_{\alpha }\left(u_{{\pi }_1,{\xi }_1},v_{{\pi }_2,{\xi }_2}\right)=\frac{A_{\alpha }\left(u_{{\pi }_1,{\xi }_1},v_{{\pi }_2,{\xi }_2}\right)-\sqrt{A_{\alpha }^2\left(u_{{\pi }_1,{\xi }_1},v_{{\pi }_2,{\xi }_2}\right)-4\alpha u_{{\pi }_1,{\xi }_1}{v}_{{\pi }_2,{\xi }_2}\left(\alpha -1\right)}}{2\left(\alpha -1\right)}, \\ \mathrm{for} \alpha \in \left(0,1\right){\displaystyle \cup }\left(1,+\infty \right), \end{gathered}$$

(11)

Here, $u_{{\pi }_1,{\xi }_1}=F_1\left(r_1\right)=\sum _{i=1}^{r_1}{\pi }_1\left\{\left(\begin{array}{c}{m}_1-1\\ i-1\end{array}\right){\xi }^{m_1-i}{\left(1-{\xi }_1\right)}^{i-1}+\left(1-{\pi }_1\right)\frac{1}{m_1 }\right\}$, where $m_1$ is the (fixed) number of categories of the first CUB model. This is similar for $v_{{\pi }_2,{\xi }_2}$ and we will from now on assume $m_1=m_2=m>3$. The copula probability mass function can be obtained as $c_{\alpha }\left(r_1,r_2\right)=C_{\alpha }\left(b_1,b_2\right)-C_{\alpha }\left(a_1,b_2\right)-C_{\alpha }\left(b_1,a_2\right)+C_{\alpha },$ where $a_i=F_i\left(r_i\right)$ and $b_i=F_i\left(r_i-1\right)$ for i = 1, 2. Model parameters are then ${\pi }_1,{\xi }_1,{\pi }_2,{\xi }_2$ and $\alpha$. Estimation can be performed by inference for the margin method (IFM), as shown in Joe and Xu (1996) [25].

3.5. Setting up the Transition Probability Matrix

Let $F\left(x,y\right)=C\left(x,y\right)$ be the cumulative function and $\frac{f_{XY}\left(x,y\right) }{f_X\left(x\right)f_Y\left(y\right)}=c\left(u,v\right)$, where $f_{Y=j|X=i}=\frac{f\left(y=j,x=i\right)}{f\left(x=i\right)}$. Then, $f_{Y|X=i}=c\left(u,v\right)f_Y\left(y=j\right)$. For the notation convenient, let the choice state at time point $t+1,$ $s_{i\left(t+1\right)}:=s_i^{\prime }$ and the state at $t, s_{i\left(t\right)}:=s_i$, where $s_i^{\prime },s_i$ are ordered labels (ranks) of the BW choice pairs and $s_i^{\prime },s_i=1,2,\dots ,\tau$; $i=1,2,\dots ,G,$ and time $t=1,2,\dots ,T.$ The transition probability is denoted as $P_{is{s}^{\prime }}^t=P^t\left(s_i^{\prime }|s_i\right)$. Since we do not have publicly available data from Flynn et al. (2007) [6] experiment, we have to simulate data under the paired model and distributed them among different choice sets. We simulated a number of different datasets under the same set of estimated regression parameters to prioritize choice pairs based on their popularity. That is, for a given choice set $\mathrm{C},$ we counted the number of times an individual picked the choice pair $\left(x_{ij},x_{i{j}^{\prime }}\right)$. Finally, all the choices are arranged according to their popularity; that is, from the highest count to lowest count. This proceeding gave elements in a choice set, an ordinal power and it is a necessary requirement to incorporate the CUB model for choice data. Further, we adjusted the function for the repetition of selections in order to keep the transition probability formula as a probability distribution function. Then, the transition probability formula can be written as

$$c^{*}\left(u=g,v=h\right)=\{\begin{array}{c}\frac{c\left(u=g,v=h\right) n}{C_p} ;C_p\ne 0\\ 0 ;C_p=0\end{array}$$

(12)

where $n$ is the number of repetitions for a selection of the pair $\left(x_{ij},x_{i{j}^{\prime }}\right)$ in the selected choice set and $C_p=\sum _{r=1}^{N}c\left(u=g,v\right).I\left(c\left(u=g,v\right)\right)$, $c$ is the density function of the Placket copula as shown in Equation (11). Let $N$ be the size of the random sample, $g, h=1,2,\dots ,\tau$, and $\tau$ be the number of BW choice pairs considered in each choice set. Then, the probability transition matrix for the $i^{th}$ choice set is in the form of

$$\begin{gathered} c_i^{*}\left(u=g,v=h\right)=\left(\begin{array}{ccc}{P}_{i11}& P_{i12}& \begin{array}{cc}\dots & P_{i1\tau }\end{array}\\ P_{i21}& P_{i22}& \begin{array}{cc}\dots & P_{i2\tau }\end{array}\\ \begin{array}{c}⋮\\ P_{i\tau 1}\end{array}& \begin{array}{c}⋮\\ P_{i\tau 2}\end{array}& \begin{array}{c}\begin{array}{cc}\ddots & ⋮\end{array}\\ \begin{array}{cc}\dots & P_{i\tau \tau }\end{array}\end{array}\end{array}\right)=\left(P_{is{s}^{\prime }}\right) \\ i=1,2,\dots ,G;, s, s^{\prime }=1,2,\dots ,\tau =K\left(K-1\right). \end{gathered}$$

(13)

$P_{is{s}^{\prime }}$ is denoted as the transition probability when a consumer changing his choice selection state $s$ at time $t$ to $s^{\prime }$ at time $t+1$. For all the transitions, $\sum _{s^{\prime }=1}^{\tau }{P}_{is{s}^{\prime }}=1$ (e.g., the sum of the row values in Equation (13) is 1) and if there are no possible transitions $\sum _{s^{\prime }=1}^{\tau }{P}_{is{s}^{\prime }}=0$. The ordering scheme is such that for $\left|s_i^{\prime }-s_i\right|\le \mathsf{\omega }$, $\sum P_{s{s}^{\prime }}>\sum P_{s{s}^{″}}$ for any other $s^{″}\ne s^{\prime }$ for the $i^{th}$ choice set, where $\mathsf{\omega }$ is some small-time epoch number. That is, in simple words the row sum of the probability of the region $\left|s_i^{\prime }-s_i\right|\le \mathsf{\omega }$ is the highest portion than any probability sum in the same row of the probability transition matrix in the $i^{th}$ choice set. We introduce a constraint $\left|s_i^{\prime }-s_i\right|\le \mathsf{\omega }$ for our copula model because the nature of human behavior is that they have a greater tendency to grab options that are close at hand or those within the nearest range rather than picking up alternatives further away in time. The data were simulated based on that choice behavior assumption; that is, $\left|s_i^{\prime }-s_i\right|\le \mathsf{\omega }$. The transition matrix may be either stationary or dynamic in nature. In this manuscript, we considered only the stationary case. It is important to state that the transition probability matrix is calculated by considering all $N$ individuals in the experiment. The row sum of the transition matrix is kept as 0 if there is no transition from $s$ to $s^{\prime }$. These proprieties can be seen in

Table 4. Probability transition matrix for choice set 1; $\left(n=100\right)$.

$( {\mathit{P}}_{\mathit{i}\mathit{s}{\mathit{s}}^{\prime }})=$
	1	2	3	4	5	6	7	8	9	10	11	12	13	14	15	16	17	18	19	20
1	0.189	0.226	0.252	0.332	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
2	0.058	0.760	0.000	0.115	0.067	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
3	0.034	0.085	0.783	0.000	0.062	0.036	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
4	0.016	0.160	0.369	0.397	0.058	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
5	0.000	0.000	0.000	0.000	0.882	0.000	0.000	0.118	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
6	0.000	0.000	0.371	0.000	0.177	0.311	0.000	0.142	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
7	0.000	0.000	0.000	0.269	0.159	0.093	0.415	0.064	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
8	0.000	0.000	0.000	0.000	0.000	0.327	0.000	0.673	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
9	0.000	0.000	0.000	0.000	0.000	0.421	0.000	0.000	0.579	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
10	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.661	0.339	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
11	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
12	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
13	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
14	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000	0.000	0.000	0.000
15	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
16	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000	0.000
17	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000
18	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000
19	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
20	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000

and

Table 5. Probability transition matrix for choice set 1; $\left(n=500\right)$.

$({\mathit{P}}_{\mathit{i}\mathit{s}{\mathit{s}}^{\prime }})=$
	1	2	3	4	5	6	7	8	9	10	11	12	13	14	15	16	17	18	19	20
1	0.441	0.266	0.178	0.116	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
2	0.101	0.557	0.197	0.112	0.033	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
3	0.113	0.171	0.631	0.046	0.012	0.026	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
4	0.189	0.191	0.219	0.355	0.022	0.016	0.007	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
5	0.000	0.076	0.434	0.083	0.377	0.016	0.000	0.015	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
6	0.000	0.000	0.263	0.075	0.121	0.459	0.000	0.053	0.027	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
7	0.000	0.000	0.000	0.204	0.146	0.052	0.549	0.048	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
8	0.000	0.000	0.000	0.000	0.053	0.037	0.034	0.730	0.071	0.037	0.038	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
9	0.000	0.000	0.000	0.000	0.000	0.091	0.000	0.000	0.865	0.044	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
10	0.000	0.000	0.000	0.000	0.000	0.000	0.078	0.000	0.000	0.750	0.171	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
11	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
12	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.116	0.000	0.000	0.884	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
13	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.830	0.170	0.000	0.000	0.000	0.000	0.000	0.000
14	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000	0.000	0.000	0.000
15	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000	0.000	0.000
16	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000	0.000
17	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.136	0.000	0.717	0.147	0.000	0.000
18	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000
19	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
20	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000

.

4. Utility in Time with Dynamic Programming

In this part, we incorporate MDPs and Dynamic Programming to calculate utility in time with the help of transition probabilities. Bellman (1954, 1956) [26,27] introduced a dynamic programming technique to evaluate the utility value function, also known as Bellman’s equation. At each time step, dynamic programming uses numerical methods to evaluate the value function moving backwards in time. Rust (1994, 2008) [28,29] and Ellickson (2011) [30] presented a numerical method with steps for evaluating DCEs as MDPs.

The value function for the DCEs is defined by at each time $t$ as:

$$V_t\left(X_t,{\epsilon }_t\right)=\underset{d_t\in D}{\max }{{\displaystyle \sum }}_{t^{\prime }=t}^T{P}_{s{s}^{\prime }}^{t^{\prime }}\left({\gamma }^{t^{\prime }-t}U\left(X_{t^{\prime }},d_{t^{\prime }}\right)+\epsilon \left(d_{t^{\prime }}\right)|X_t\right),$$

(14)

with $P_{s{s}^{\prime }}=P\left(s_{t+1}|s_t\right)=P\left(s_{t+1}=s^{\prime }|s_t=s\right)$. This may also be defined in terms of the transition in decisions as $P\left(d_{t+1}|d_t,x_t\right)$,where $d_t\in D$ is the decision at time $t, U\left(x_t ,d_t\right)$ is the derived expected utility, and $\epsilon$ is the associated error term at time $t$, where $t=1, 2,\dots , T$. The discount utility rate in Equation (14) is given by γ ∈ (0, 1); see Feinberg and Schwartz (1994) [31]. There exist $\tau$ value functions for each of the $\tau$ alternatives in the experiments evaluated at each time point $t=1, 2,\dots , T$. The sum goes from $t$ up to $T$ because it is evaluated using a backwards recursive method; that is, we start at the last time point and work our way backwards to earlier time points.

$$\begin{gathered} V_T\left(X_T\right)=U\left(X_T,d_T\right) \\ {V}_{T-1}\left(X_{T-1},d_{T-1}\right)=U\left(X_{T-1},d_{T-1}\right)+{\displaystyle \sum }_{d_t\in D }\gamma V_T\left(X_T\right)P\left(d_T|d_{T-1}\right), \\ {V}_{T-2}\left(X_{T-2},d_{T-2}\right)=U\left(X_{T-2},d_{T-2}\right)+{\displaystyle \sum }_{d_t\in D }\gamma V_{T-1}\left(X_{T-1}\right)P\left(d_{T-1}|d_{T-2}\right), \\ ⋮ \\ {V}_t\left(X_t,d_t\right)=U\left(X_t,d_t\right)+{\sum }_{t^{\prime }=t+1}^T{\sum }_{d_t\in D }\gamma V_{t^{\prime }}\left(X_{t^{\prime }}\right)P\left(d_{t^{\prime }}|d_{t^{\prime }}\right), \end{gathered}$$

(15)

where $U_t$ is utility in time (UiT) and $V_t$ is the value function. Based on Working et al. (2018) [8], the initial utility for the selected choice pair can be calculated using:

$$V_{\mathrm{T}}=\left(\left({\widehat{\beta }}_{\mathrm{A}}+{\widehat{\beta }}_{{\mathrm{A}}_k}\right)-\left({\widehat{\beta }}_{{\mathrm{A}}^{\prime }}+{\widehat{\beta }}_{A_k^{\prime }}\right)\right)\equiv U\left(X_T,d_T\right),$$

(16)

where ${\widehat{\beta }}_{\mathrm{A}},{\widehat{\beta }}_{{\mathrm{A}}_k}$ are estimated parameters for best attribute and attribute level and ${\widehat{\beta }}_{{\mathrm{A}}^{\prime }},{\widehat{\beta }}_{A_k^{\prime }}$ are the estimated parameters for worst attribute and attribute level, respectively. Instead of this simple model, in this manuscript we used the paired model which is one of the more appropriate models for the paired data. In fact, to distribute the utility over time for the future time points, further adjustment is needed for the value function. We adjusted the value function in Equation (12) for the time in such a way that

$$V_{\mathrm{T}}=\left({\mathrm{M}}_{+}\widehat{\beta }\right){\theta }_b^t+\left({\mathrm{M}}_{-}\widehat{\beta }\right){\theta }_w^t\equiv U\left(X_T,d_T\right), \mathrm{where} M=M_{+}+M_{-}$$

(17)

$M$ is the design matrix under the paired model which contains +1, −1, and 0 described in Equations (3)–(5). We broke $M$ into $M_{+}$ and $M_{-}$, which contains $+1$s and $-1$s. Further, we introduce arbitrary time-sensitive constants ${\theta }_b^t$, ${\theta }_w^t$ to make the initial systematic component vary with time.

5. Application and the Design of the Experiment

A pilot best–worst study was conducted in summer 2005 among people aged 65 and over with the aim of informing a larger quality of life valuation exercise by Flynn et al. (2007) [6]. In this survey, there are $K=5$ attributes, which are Attachment, Security, Role, Enjoyment, and Control. Each attribute with $l_k=4$ levels. That is

• Attachment: Attach_none, Attach_little, Attach_lot, Attach_all
• Security: Security_none, Security_little, Security_lot, Security_all
• Role: Enjoyment_none, Enjoyment_little, Enjoyment_lot, Enjoyment_all
• Enjoyment: Role_none, Role_few, Role_many, Role_all
• Control: Control_none, Control_few, Control_many, Control_all

We simulated a new data set for $100$ individuals and estimated regression parameters by using a conditional logit model with the help of R package ‘support.BWS2’; see Aizaki and Fogarty (2019) [11]. The application is from a balance design where all five attributes have four- levels each; there are 320 unique BW choice pairs. We simulated data in such a way that 16 choice sets were used to cage all 320 unique BW choice pairs (see Equation (1)), by keeping each choice set with 20 choice pairs. Each person of the 100 sample should make a choice in each choice set. Figure 2 shows the orthogonal design we used to simulate the data. The shape of the questionnaire design based on the orthogonal array is shown in Figure 3. Such a questionnaire could be used in the case of a real-life data collection process. Table 1 shows the design matrix for the paired model. The attribute ‘Control’ is omitted and the fourth level in each attribute (‘Attach_all’, ‘Security_all’, ‘Enjoyment_all’, ‘Role_all’ and ‘Control_all’) is considered to be the base level.

The parameter estimated by Flynn et al. (2007) [6] and conditional logit regression estimates for the simulation are presented in Table 2. Based on Equation (17), systematic components for each choice and then the corresponding utility can be calculated. Estimations for base levels are calculated by $A_{j{l}_k}=-{\sum }_{i=1}^{k-1}{\beta }_{A_{ji}}$ as mentioned in Section 2.1.1. As per the iterative dynamic programming method introduced in Equation (15), with initial utility as described in Equation (16), utility values for a consecutive 10 time periods are calculated recursively by utilizing the transition probabilities and a discount factor γ = 0.95. As a demonstration, the results for choice set 1 is presented in this manuscript.

The BW choice pairs are arranged in order, as mentioned in Section 3.4. We ran the choice selections algorithm 10,000 times within each choice set to get a proper ordering of elements. Table 3 shows the arrangement for the BW pairs for choice set 1. This arrangement is organized according to the popularity (highest number of picks) of the BW choice pairs. That information is then used to calculate the transition probabilities at the next time point. The probability transition matrix can be calculated utilizing our new model as derived in Equation (12). It is a $20\times 20$ matrix since each choice sets contain 20 choice pairs. Figure 4 shows the highest probability corresponding to each BW choice pair in choice set 1 for three different samples sizes, which $n=100$, $n=500$, and n = 1000; that is, the highest probable transition slot at time t to t + 1. Table 4, Table 5 and

Table 6. Probability transition matrix for choice set 1; $\left(n=1000\right)$.

$({\mathit{P}}_{\mathit{i}\mathit{s}{\mathit{s}}^{\prime }})=$
	1	2	3	4	5	6	7	8	9	10	11	12	13	14	15	16	17	18	19	20
1	0.574	0.217	0.123	0.086	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
2	0.131	0.583	0.161	0.085	0.040	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
3	0.197	0.217	0.498	0.052	0.017	0.019	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
4	0.236	0.217	0.127	0.329	0.050	0.029	0.012	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
5	0.000	0.204	0.327	0.140	0.288	0.017	0.017	0.008	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
6	0.000	0.000	0.148	0.158	0.093	0.497	0.012	0.038	0.053	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
7	0.000	0.000	0.000	0.187	0.111	0.045	0.548	0.046	0.063	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
8	0.000	0.000	0.000	0.000	0.021	0.017	0.034	0.780	0.054	0.037	0.057	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
9	0.000	0.000	0.000	0.000	0.000	0.044	0.021	0.044	0.820	0.047	0.024	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
10	0.000	0.000	0.000	0.000	0.000	0.000	0.031	0.000	0.000	0.860	0.071	0.000	0.038	0.000	0.000	0.000	0.000	0.000	0.000	0.000
11	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.145	0.000	0.000	0.800	0.055	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000
12	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.122	0.000	0.000	0.808	0.069	0.000	0.000	0.000	0.000	0.000	0.000	0.000
13	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.927	0.073	0.000	0.000	0.000	0.000	0.000	0.000
14	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000	0.000	0.000	0.000
15	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000	0.000	0.000	0.000
16	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.893	0.000	0.000	0.107	0.000
17	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.086	0.000	0.820	0.094	0.000	0.000
18	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000	0.000
19	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000	0.000
20	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	0.000	1.000

show the corresponding transition probability matrices for choice set 1 for sample sizes 100, 500, and 1000, respectively. The constraint region, |s_i’ − s_i| ≤ 3, is highlighted in the probability transition matrix and the highest probability is squared. Nobody picked choice pairs 11, 12, 15, 19, and 20 at time t or t + 1 when n = 100. That has happened only to choice pair 19 in the n = 500 scenario. However, for n = 1000 all the choice pairs have been selected.

Table 7. Expected utilities in time (UiT) for choice set 1.

Choice	t = 1	t = 2	t = 3	t = 4	t = 5	t = 6	t = 7	t = 8	t = 9	t = 10
1	1.2346	1.2613	1.3110	1.3983	1.5419	1.7660	2.1003	2.5814	3.2520	4.1576
2	0.9365	1.0047	1.1107	1.2680	1.4947	1.8143	2.2575	2.8627	3.6758	4.7452
3	0.8028	0.7827	0.7642	0.7572	0.7739	0.8287	0.9381	1.1201	1.3917	1.7638
4	1.2050	1.2500	1.3260	1.4446	1.6215	1.8769	2.2366	2.7325	3.4027	4.2880
5	0.5933	0.6260	0.6770	0.7530	0.8628	1.0178	1.2323	1.5241	1.9134	2.4200
6	0.9295	0.9248	0.9264	0.9418	0.9804	1.0536	1.1745	1.3574	1.6166	1.9616
7	0.9481	0.9799	1.0310	1.1087	1.2231	1.3867	1.6152	1.9276	2.3447	2.8856
8	1.0941	1.1007	1.1168	1.1484	1.2032	1.2902	1.4205	1.6063	1.8600	2.1914
9	1.1719	1.1818	1.2027	1.2404	1.3022	1.3971	1.5357	1.7300	1.9920	2.3310
10	0.7217	0.7911	0.8881	1.0144	1.1734	1.3711	1.6166	1.9242	2.3144	2.8145
11	0.3888	0.4186	0.4630	0.5224	0.5973	0.6890	0.8003	0.9358	1.1030	1.3135
12	0.1282	0.1306	0.1330	0.1335	0.1301	0.1200	0.0997	0.0648	0.0096	−0.0736
13	0.4320	0.4390	0.4430	0.4389	0.4205	0.3802	0.3089	0.1959	0.0295	−0.2022
14	0.7528	0.7431	0.7243	0.6911	0.6365	0.5516	0.4256	0.2452	−0.0044	−0.3383
15	0.3652	0.3564	0.3409	0.3160	0.2779	0.2214	0.1397	0.0234	−0.1402	−0.3686
16	1.4454	1.4596	1.4725	1.4786	1.4711	1.4417	1.3807	1.2765	1.1157	0.8841
17	0.8189	0.8620	0.9079	0.9462	0.9642	0.9469	0.8772	0.7351	0.4989	0.1469
18	1.5412	1.4548	1.3269	1.1446	0.8907	0.5427	0.0715	−0.5588	−1.3894	−2.4609
19	1.0774	1.0947	1.1169	1.1398	1.1576	1.1636	1.1492	1.1046	1.0170	0.8699
20	1.2778	1.2731	1.2608	1.2341	1.1843	1.0997	0.9652	0.7611	0.4613	0.0295

shows the utilities (UiT) calculated for time points t = 1, 2,…, 10. Figure 5 shows the utility–time distribution of selected choices pairs (2, 7, 14, 16, and 18) in choice set 1. It highlights the dynamic utility in time.

6. Conclusions

In general, the choices are made at different time points and are not IIA for any individual. They depend on not only qualitative factors such as priority setting, feeling, and uncertainty, but also quantitative covariates such as measure of sequential selection probability in the process. By considering a specific popular model, we capture and monitor the dependent structure and calculate the transition probability between consecutive choices at adjacent time points. We were able to include parameters like feeling and uncertainty to the model. This new approach to calculating transition probabilities extend the ideas proposed by Piccolo (2003) [16] and Piccolo et al. (2019) [32]. Further, we could justify our initial assumption of the human choice behavior of $\left|s_i^{\prime }-s_i\right|\le \mathsf{\omega }$, that they have a greater tendency to grab options that are close at hand or those within the nearest range rather than picking up alternatives further away in time.

Data to quantify the BW in the DCE models are mainly found outside of academia. In this paper, we have simulated the data. Such an issue is one limitation of the paper, but the approach has the advantage of bringing and completing the ranking. The paper could also gain if more covariates were given. Based on Flynn’s ‘Quality of life experiment’, we choose value $\mathsf{\omega }$ as 3. The accuracy of the CO-CUB model can be further improved relaxing that constraint and by adding more latent covariates and updating knowledge of the ranking labels of the choice pairs within the choice set. We leave such ideas for future research.