# The following lines of code construct choice sets for Case 3 BWS
# using rotation.design() in the package support.CEs; create a dataset
# for the analysis using bws3.dataset(); and then conduct a maxdiff model
# analysis on the basis of a conditional logit model specification
# using clogit() in the package survival.
# Load packages
library(support.CEs)
library(survival)
# Reproduce old results
if(getRversion() >= "3.6.0") RNGkind(sample.kind = "Rounding")
# Create choice sets
dsgn <- rotation.design(
attribute.names = list(
A = c("A1", "A2", "A3"), B = c("B1", "B2", "B3"),
C = c("C1", "C2", "C3"), D = c("D1", "D2", "D3")),
nalternatives = 4, nblocks = 1,randomize = TRUE, seed = 9876)
# Prepare response variables
dat <- data.frame(
IDvar = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10),
BLOCK = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1),
q1b1 = c(2, 1, 1, 1, 1, 2, 1, 1, 1, 1),
q1w1 = c(4, 2, 4, 3, 4, 4, 4, 3, 3, 4),
q1b2 = c(1, 4, 2, 2, 2, 3, 3, 2, 2, 2),
q1w2 = c(3, 3, 3, 4, 3, 1, 2, 4, 4, 3),
q2b1 = c(1, 4, 2, 1, 1, 4, 4, 4, 1, 2),
q2w1 = c(3, 3, 3, 3, 3, 1, 3, 2, 2, 3),
q2b2 = c(4, 1, 4, 4, 4, 3, 1, 1, 4, 4),
q2w2 = c(2, 2, 1, 2, 2, 2, 2, 3, 3, 1),
q3b1 = c(4, 2, 1, 4, 2, 2, 4, 4, 3, 4),
q3w1 = c(1, 3, 3, 1, 1, 1, 1, 1, 4, 1),
q3b2 = c(2, 1, 2, 3, 4, 3, 2, 3, 1, 3),
q3w2 = c(3, 4, 4, 2, 3, 4, 3, 2, 2, 2),
q4b1 = c(2, 2, 2, 1, 2, 3, 4, 4, 2, 1),
q4w1 = c(3, 4, 3, 3, 4, 4, 3, 2, 4, 2),
q4b2 = c(1, 1, 4, 2, 1, 1, 1, 1, 1, 3),
q4w2 = c(4, 3, 1, 4, 3, 2, 2, 3, 3, 4),
q5b1 = c(1, 1, 1, 3, 3, 2, 3, 3, 3, 1),
q5w1 = c(2, 2, 2, 2, 4, 4, 2, 2, 2, 2),
q5b2 = c(3, 3, 3, 1, 1, 1, 1, 1, 4, 3),
q5w2 = c(4, 4, 4, 4, 2, 3, 4, 4, 1, 4),
q6b1 = c(1, 3, 3, 1, 2, 4, 3, 3, 3, 3),
q6w1 = c(4, 2, 1, 3, 1, 1, 1, 4, 1, 1),
q6b2 = c(3, 1, 2, 4, 4, 3, 2, 1, 2, 4),
q6w2 = c(2, 4, 4, 2, 3, 2, 4, 2, 4, 2),
q7b1 = c(4, 4, 2, 2, 3, 4, 4, 4, 2, 1),
q7w1 = c(3, 2, 4, 1, 1, 3, 3, 3, 1, 3),
q7b2 = c(2, 3, 3, 4, 4, 2, 2, 2, 4, 4),
q7w2 = c(1, 1, 1, 3, 2, 1, 1, 1, 3, 2),
q8b1 = c(3, 2, 2, 2, 1, 3, 2, 2, 3, 3),
q8w1 = c(4, 4, 1, 4, 4, 4, 3, 1, 2, 2),
q8b2 = c(2, 1, 3, 1, 2, 2, 4, 3, 1, 4),
q8w2 = c(1, 3, 4, 3, 3, 1, 1, 4, 4, 1),
q9b1 = c(4, 1, 3, 4, 3, 3, 3, 1, 4, 4),
q9w1 = c(2, 4, 1, 3, 4, 1, 1, 2, 3, 1),
q9b2 = c(1, 2, 2, 2, 2, 4, 2, 3, 2, 2),
q9w2 = c(3, 3, 4, 1, 1, 2, 4, 4, 1, 3))
# Store names of response variables
rsp.vars <- list(
q1 = c("q1b1", "q1w1", "q1b2", "q1w2"),
q2 = c("q2b1", "q2w1", "q2b2", "q2w2"),
q3 = c("q3b1", "q3w1", "q3b2", "q3w2"),
q4 = c("q4b1", "q4w1", "q4b2", "q4w2"),
q5 = c("q5b1", "q5w1", "q5b2", "q5w2"),
q6 = c("q6b1", "q6w1", "q6b2", "q6w2"),
q7 = c("q7b1", "q7w1", "q7b2", "q7w2"),
q8 = c("q8b1", "q8w1", "q8b2", "q8w2"),
q9 = c("q9b1", "q9w1", "q9b2", "q9w2"))
# Store names of attributes
attributes <- c("A", "B", "C", "D")
# Create a dataset
bws3dat <- bws3.dataset(
data = dat, id = "IDvar", response = rsp.vars,
choice.sets = dsgn, categorical = attributes, model = "maxdiff")
# Fit the model
clogit(RES ~ A2 + A3 + B2 + B3 + C2 + C3 + D2 + D3 + strata(STR), bws3dat)
Run the code above in your browser using DataLab