Generar un experimento de elección discreta con un diseño D-eficiente usando R

Una vez que hemos visto la base teórica de los diseños D-eficientes para los experimentos de elección discreta, podemos pasar a generar nuestro diseño D-eficiente en R. Si nunca has usado R, no te preocupes, voy a empezar desde el principio.

Vamos a necesitar R (la versión simple del programa estadístico) y RStudio (que incorpora una interfaz visual a la versión simple. Instala en orden los dos siguientes programas:

Una vez que tengamos instalados ambos programas, abrimos RStudio y creamos un nuevo «R Script» en File/New File/R Script. Los R Script son libretas en las que guardamos código para poder ejecutarlo en cualquier momento. Con Ctrl+S podemos guardar el script en nuestra carpeta de trabajo (así, una vez hemos seleccionado donde queremos guardarlo, cada vez que pulsemos Ctrl+S, se guardará una versión actualizada en el mismo lugar).

A partir de aquí, vamos a ir viendo código que tendremos que pegar en nuestro R Script para después ejecutar. Para ejecutar el código que hayamos pegado (o escrito) solo tenemos que seleccionarlo y pulsar Ctrl+Intro (o «Run» en la propia interfaz del programa). Al final de este post puedes descargar el R Script al completo, pero te recomiendo que vayas paso a paso pegando tu propio código.

Voy a comenzar con una línea de código opcional pero muy útil.

setwd(dirname(rstudioapi::getActiveDocumentContext()$path))

Al copiar y ejecutar ese código en R, el directorio de trabajo del programa cambia automáticamente al mismo lugar en el que está guardado el Script. De esta forma podemos trabajar con directorios relativos (no hace falta que entremos mucho en esto, pero facilitará nuestra vida).

Ahora vamos a utilizar la librería «idefix» de R (damos las gracias a Frits Traets, Daniel Gil Sánchez y Martina Vandebroek https://www.jstatsoft.org/article/view/v096i03). Para ello, primero la instalamos y después la activamos:

install.packages("idefix")
library(idefix)

Ahora vamos a comenzar a definir las características del diseño. La librería funciona con una sola fórmula, pero a esa fórmula hay que insertarle toda la información que necesita para saber cómo tiene que ser el diseño. Lo primero, lógicamente, son los atributos y niveles. Los atributos y niveles deben ser insertados en la fórmula como un vector. Por ejemplo, imaginemos que tenemos un diseño con 4 atributos con 2, 3, 4 y 5 niveles cada atributo. El vector que queremos insertar sería:

(2, 3, 4, 5)

Voy a crear este vector en R y lo voy a llamar «niveles»:

niveles <- c(2,3,4,5)

El siguiente paso es indicarle al programa cómo queremos que nos codifique la salida del diseño. En general, existen dos tipos de codificaciones («dummy» y «efectos»), ambas equivalentes aunque la dummy es sensiblemente más sencilla de entender. Para especificar que queremos que el diseño use la codificación dummy necesitamos especificar un vector del mismo tamaño que el anterior, pero en lugar de indicar cada nivel por atributo, incluiremos una «D» de dummy:

codif <- c("D","D","D","D")

Ahora pasaremos a un punto un poco más complejo. Cómo explicamos en la entrada anterior, ¿Qué es un diseño D-eficiente?, los algoritmos que crean los diseños D-eficientes requieren coeficientes previos. Si no tenemos, lo que estamos diseñando es un piloto y no tenemos coeficientes previos (aunque muchos investigadores usen este diseño piloto como encuesta final). Cuando no tenemos coeficientes previos, los coeficientes previos que le daremos al algoritmo serán ceros. Fácil, ¿no? Bien, pues ahora tenemos que especificar un vector con los coeficientes previos (en nuestro caso ceros, pero si ya se hubiera realizado un piloto serían los coeficientes de este piloto). Vale, pues ahora hay que responder a una segunda pregunta, ¿cuántos coeficientes tendrá mi modelo? Es fácil, los experimentos de elección discreta, al estimarse los modelos, siempre omiten un nivel por atributo (que será restringido a ser cero y servirá como base para el resto). Por tanto, podemos generalizar cuántos coeficientes tendremos de la siguiente forma

  • Para diseños sin opción nula: l-k
  • Para diseños con opción nula: l-k+1

donde k es el número de atributos y l es el número de niveles. Aquí hemos realizado la primera distinción, si queremos choice sets con o sin opción nula, el camino que deberemos tomar es diferente. Yo voy a comenzar explicando cómo hacerlo para choice sets sin opción nula y después explicaré el otro caso.

Choice sets sin opción nula

Siguiendo con la especificación de atributos y niveles que hemos hecho antes, tendríamos que crear un vector de coeficientes previos con 14-4=10 ceros:

previos <- c(0,0,0,0,0,0,0,0,0,0)

Después de esto, podemos automatizar el resto del código que necesitamos introducir como:

I <- diag(length(previos))
set.seed(1234)
sim <- MASS::mvrnorm(n=100, mu=previos, Sigma=I)

Lo que hemos realizado en las anteriores líneas es crear una matriz identidad con de tamaño 10x10 (por el número de coeficientes previos) que nos ha servido para calcular una distribución normal por cada coeficiente previo con media el coeficiente previo. El comando intermedio «set.seed(1234)» es una contraseña para que esta distribución normal siempre se simule de la misma forma y a todos (tanto a mi como a ti) nos salga el mismo diseño (no quiero entrar mucho en esto, pero los diseños siempre están sometidos a cierta aleatoriedad por lo que ejecutando el mismo código, sin establecer esa «contraseña», tendríamos diseños diferentes pero muy similares en eficiencia). Finalmente, ya podemos ejecutar el comando que nos va a generar el diseño eficiente:

dis <- CEA(lvls=niveles, coding=codif, n.alts=2, n.sets=16, par.draws=sim, no.choice=FALSE, best=TRUE)

Vamos paso a paso. Lo que estamos haciendo es guardar en el almacén «dis» los resultados de la función «CEA», que es la que nos genera el diseño D-eficiente. Dentro de la función, escribimos las opciones:

  • lvls: hace referencia al vector de atributos y niveles que hemos creado al principio.
  • coding: hace referencia al vector que especificaba qué codificación íbamos a usar.
  • n.alts: ¿qué número de alternativas por choice set queremos? Yo siempre uso 2.
  • n.sets: ¿qué número de choice sets queremos? Aquí ponemos el que queramos. Yo siempre uso 16 y luego divido la encuesta en dos bloques (porque 16 preguntas por encuestado pueden ser muchas).
  • par.draws: son las simulaciones de la distribución normal que hemos creado antes.
  • no.choice: FALSE porque no queremos usar alternativa nula.
  • best: nos pregunta si deseamos encontrar el mejor diseño (pues claro, ¿no?)

Una vez ejecutamos este código, se nos ha creado un almacen (un data frame) al que podremos acceder para ver el diseño.

Este almacen contiene cuatro elementos, nuestro diseño se encuentra dentro del elemento «design» al que podemos acceder escribiendo:

dis$design

La salida es la siguiente:

             Var12 Var22 Var23 Var32 Var33 Var34 Var42 Var43 Var44 Var45
set1.alt1      1     0     1     1     0     0     0     1     0     0
set1.alt2      0     0     0     0     0     0     0     0     0     0
set2.alt1      0     1     0     0     0     1     0     0     0     0
set2.alt2      0     0     0     0     0     0     0     0     0     1
set3.alt1      0     0     1     1     0     0     0     0     0     0
set3.alt2      1     0     0     0     0     0     0     1     0     0
set4.alt1      1     0     0     0     1     0     0     0     0     0
set4.alt2      1     0     1     0     0     1     0     0     1     0
set5.alt1      1     1     0     1     0     0     0     0     0     1
set5.alt2      1     0     0     0     0     1     0     0     0     0
set6.alt1      1     0     0     1     0     0     0     0     1     0
set6.alt2      0     0     1     0     1     0     0     0     0     1
set7.alt1      0     0     0     0     0     1     0     1     0     0
set7.alt2      1     1     0     0     0     0     0     0     1     0
set8.alt1      0     0     0     1     0     0     0     0     0     0
set8.alt2      1     0     1     0     1     0     0     0     1     0
set9.alt1      1     1     0     0     0     0     0     1     0     0
set9.alt2      1     0     0     0     0     1     0     0     0     1
set10.alt1     0     0     0     1     0     0     0     0     1     0
set10.alt2     1     0     1     0     0     0     0     0     0     0
set11.alt1     1     0     0     0     1     0     0     0     1     0
set11.alt2     0     0     1     0     0     0     1     0     0     0
set12.alt1     1     0     1     1     0     0     0     0     0     0
set12.alt2     0     1     0     0     0     1     1     0     0     0
set13.alt1     1     0     0     0     0     0     1     0     0     0
set13.alt2     0     1     0     0     1     0     0     1     0     0
set14.alt1     1     1     0     0     1     0     1     0     0     0
set14.alt2     0     0     0     0     0     0     0     0     1     0
set15.alt1     1     0     0     0     1     0     1     0     0     0
set15.alt2     0     1     0     0     0     0     0     0     0     1
set16.alt1     0     0     1     0     1     0     1     0     0     0
set16.alt2     1     1     0     0     0     1     0     0     0     1

Como vemos, las filas contienen un título que nos indican a qué choice set (set) se refieren y a qué alternativa (alt) describen. La interpretación de las columnas es un poco diferente. La columna Var12 hace referencia al nivel 2 del atributo 1. Como puede observarse, todas las columnas empiezan por VarX2 porque VarX1 es el nivel omitido (restringido a cero) del que antes hemos hablado. El resto de columnas toma el valor 1 cuando el nivel de su atributo está presente en la alternativa representada por la fila y 0 cuando no. Cuando en una fila todos los niveles de un atributo toman el valor cero (por ejemplo en el set1.alt2 el atributo 2 -Var2X), significa que el nivel que está presente en la alternativa es el omitido.

Con esta matriz, ya tendríamos listo nuestro diseño D-eficiente. Sin embargo, os voy a mostrar cómo exportar el diseño a un excel:

D <- as.data.frame(dis$design)
D <- cbind(" "=rownames(D), D)
#install.packages("writexl")
library(writexl)
write_xlsx(D, "diseño.xls")

Si hemos seguido todos los pasos adecuadamente, en la misma carpeta en la que está guardado nuestro Script se ha debido de generar un .xls con los resultados del diseño.

Choice sets con opción nula

Si antes hemos creado un vector de coeficientes previos con 10 ceros, para este caso necesitamos 11 ceros.

previos <- c(0,0,0,0,0,0,0,0,0,0,0)

Igualmente, necesitaremos crear las simulaciones de la misma forma en la que lo hemos hecho antes:

I <- diag(length(previos))
set.seed(1234)
sim <- MASS::mvrnorm(n=100, mu=previos, Sigma=I)

Y ahora aquí tenemos un punto importante. Tenemos que separar la lista de simulaciones que hemos creado en «sim» en dos partes, una para el coeficiente de la alternativa nula y otra para el resto de coeficientes. Sin entrar mucho en detalles, esto se puede hacer con el siguiente código:

sim <- list(sim[,1:1], sim[, 2:(length(previos))])

Y, por último, tenemos que crear un vector en el que le indiquemos al algoritmo la posición en la que se encontrará la alternativa nula:

constante <- c(0,0,1)

Por fin, podemos usar el comando para generar el diseño experimental:

dis <- CEA(lvls=niveles, coding=codif, n.alts=3, n.sets=16, alt.cte=constante, par.draws=sim, no.choice=TRUE, best=TRUE)

Nótense varios cambios con respecto al caso anterior:

  • Ahora tenemos tres alternativas en n.alts (las normales + nula)
  • En alt.cte especificamos el vector de posición de la alternativa nula.
  • no.choice tiene el valor TRUE esta vez.

La interpretación de los resultados y la exportación a excel es igual que en el caso anterior. Para exportar, ejecutamos el siguiente código:

D <- as.data.frame(dis$design)
D <- cbind(" "=rownames(D), D)
#install.packages("writexl")
library(writexl)
write_xlsx(D, "diseño.xls")

Descargar el código: enlace

One thought on “Generar un experimento de elección discreta con un diseño D-eficiente usando R

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *