Spanish (pdf)
Article in xml format
Article references
Automatic translation
Send this article by e-mail
Cited by SciELO 
Cited by Google 
Similars in
SciELO 
Similars in Google 
Permalink
1 Introducción
Los diseños óptimos permiten utilizar de manera eficiente los recursos requeridos para realizar un experimento en las áreas de investigación donde se desea obtener alguna información de interés mediante el uso de la estadística. De este modo, buscan las mejores condiciones donde se debe realizar un experimento con el fin de que se alcancen algunas propiedades estadísticas, por ejemplo: minimizar la varianza de los estimadores de los parámetros, minimizar el volumen del elipsoide de confianza asociado a los parámetros del modelo, entre otros. La búsqueda de estos diseños requiere del cumplimiento de los siguientes supuestos: independencia, normalidad y homogeneidad de la varianza del término de error del modelo. En la literatura existen trabajos donde se proponen métodos para hallar diseños óptimos cuando los supuestos anteriores no se cumplen. Por ejemplo, en el caso donde la varianza del error no es constante hay varios aportes para la búsqueda de diseños D-óptimos tanto para modelos lineales como no lineales (Dette y Müller, 2013 [1], Dette y Wong, 1999 [2], Downing, Fedorov y Leonov, 2001 [3], Gaviria y López-Ríos, 2014 [4], Burridge y Sebastiani, 1994 [5], Atkinson y Cook, 1995 [6], Dette y Wong, 1996 [7], Fang y Wiens, 2000 [8]). En todos los trabajos consultados el criterio de optimalidad propuesto depende de la elección de un valor local para los parámetros del modelo. Una forma de evitar esta dependencia es considerar una distribución a priori para el vector de parámetros del modelo e incorporarla en el criterio de optimalidad que se va a optimizar. En el caso de modelos lineales, existen trabajos en esta dirección (Duarte y Wong, 2015 [9], Abebe, Tan, Van Breukelen, Serroyen y Berger, 2014 [10], Burghaus and Dette [11], Buzoianu y Kadane (2019) [12], Braess y Dette [13],Atkinson y Cook, 1995 [6], Woods, Overstall, Adamou y Waite, 2017 [14], Goudarzi, 2019 [15], Bodunwa y Fasoranbaku, 2020 [16]). Mosquera y López-Ríos, 2019 [17] investigan el efecto de diferentes distribuciones a priori en los diseños D-óptimos Bayesianos en el caso de un modelo no lineal donde el supuesto de incorrelación no se cumole. Adicionalmente, se han publicado recientemente trabajos relacionados con programación para la obtención de diseños D-óptimos Bayesianos (ver Duarte y Wong, 2015 [18] y Abebe et. al. 2017[19]).
En este artículo se investiga el comportamiento de los diseños D-óptimos en modelos no lineales heteroscedásticos cuando se incorpora una distribución a priori asociada a los parámetros del modelo con el fin de dar respuesta a las siguientes preguntas de investigación: ¿cómo se afecta el número de puntos de soporte del diseño D-óptimo con respecto a la elección de la distribución a priori? y ¿qué características deben tener las distribuciones a priori con el fin de que los diseños D-óptimos obtenidos sean robustos a la elección de esta distribución?.
En la sección 2 se dan algunos elementos básicos sobre los diseños óptimos considerando el modelo no lineal heteroscedástico donde se modela la varianza como una potencia de la media. Se presenta la matriz de información para este modelo, se define el criterio de optimalidad y se enuncia el teorema general de equivalencia. Luego, se define el criterio Dπ-optimalidad y se enuncia el teorema de equivalencia para este criterio. Finalmente, se define la Dπ-eficiencia, la cual permite comparar los diferentes diseños obtenidos. En la sección 3 se realiza la descripción de la metodología propuesta para el modelo considerado. En la sección 4 se halla el diseño D-óptimo local para el ejemplo mencionado. Luego, para evitar la dependencia de este vector, se construyen distribuciones a priori discretas y continuas, lo cual se explica en las subsección 4.4 en las secciones 4.4.1 y 4.4.2 se obtienen diseños Dπ-óptimos para el modelo inicial con tres y cuatro puntos de soporte y se comparan con el diseño D-óptimo local mediante el cálculo de la eficiencia. Finalmente, en la sección 4.4.3, para estudiar la dependencia del diseño a la elección de la apriori, se realiza un estudio vía simulación, donde se toman diferentes distribuciones que pertenecen a la familia de escala de modo que también pertenezca la distribución a priori continua construída a partir del vector de parámetros local y se puede concluir que el diseño es robusto a la elección de la distribución a priori.
2 Conceptos preliminares
Se considera el siguiente modelo donde la variable respuesta Y está relacionada con las variables explicativas
a través de la siguiente expresión:
donde se dispone de n corridas experimentales, x varía en un espacio compacto
dotado de una σ-álgebra
es una función no lineal en
y los
son independientes, tienen distribución normal con media cero y varianza
, donde
(Ritz, 2008 [20]).
Asociado al modelo (1), se define un diseño ξ como:
donde
es una medida de probabilidad definida en , tal que ξ tiene soporte finito, Sop(ξ) =
=
, d es el número de puntos de soporte de ξ y las observaciones
se hacen en
con pesos asociados aproximadamente proporcionales a
. Si 0 < wj ≤1, se dice que el diseño es continuo o aproximado (López-Ríos, Ramos, 2007) [21].
2.1 Matriz de información
Al considerar la estructura de la varianza de los errores del modelo definida por
, interesa hallar la expresión para la matriz de información de Fisher asociada al diseño ξ, de modo que ésta incorpore dicha estructura.
Ahora, para el modelo (1) donde
es una función positiva, se tiene que la matriz de información de Fisher para el diseño aproximado es (Downing, Fedorov y Leonov, 2001 [3]
donde
para θ = (α, β)T .)
De acuerdo a la expresión anterior se observa que la matriz de información para este modelo se puede escribir de la siguiente manera:
donde
donde.
para j = 1, . . . , d.
Gaviria y López-Ríos (2014) [4] muestran, como una consecuencia del resultado anterior, la forma de la matriz de información para el caso en el cual la estructura de la varianza del error es
=
:
donde
para j = 1, . . . , d.
Una vez encontrada la matriz de información de Fisher, interesa hallar el diseño que optimice algún funcional de valor real de esta matriz, denominado criterio de optimalidad.
2.2 Criterios de optimalidad
Como se puede observar en el modelo no lineal considerado en este artículo, donde los errores tienen una distribución normal y su varianza es una función de la media (ver ecuación (1)], la matriz de información depende del vector de parámetros local y por lo tanto el diseño depende de los valores que tomen las componentes de este vector. Una manera de evitar esta dependencia es considerar una distribución a priori para el vector de parámetros. De acuerdo a lo anterior se plantea el criterio Dπ-optimalidad de la siguiente manera:
Criterio Dπ-optimalidad
donde π(θ) es una distribución de probabilidad discreta. De manera análoga se plantea este criterio para el caso en el cual la distribución π es continua:
Se dice que un diseño ξ∗ es Dπ-óptimo si maximiza el criterio
dado en las expresiones (4) y (5).
2.3 Teorema de equivalencia
Uno de los objetivos principales de este artículo es extender el teorema de equivalencia de Kiefer y Wolfowitz [22], cuando se consideran distribuciones a priori para el vector de parámetros local en modelos no lineales heteroscedásticos. A continuación se da la definición de derivada direccional y se enuncian los teoremas de equivalencia para el criterio Dπ-optimalidad cuando la distribución a priori es discreta y cuando esta distribución es continua.
Deftnición 1. Dados dos diseños
, se define la derivada direccional de ξ en la dirección de η por:
Ahora, se puede probar que la derivada direccional puede expresarse como se indica en el siguiente lema (para detalles de la prueba ver [23])
Lema 1. Cuando la distribución a priori π es discreta, se cumple que
donde
. En el caso continuo, se tiene que
El siguiente teorema es válido tanto para distribuciones a priori discretas como continuas. Los detalles se pueden consultar en [24].
Teorema 1. Sean M (ξ, θ) la matriz de información del diseño ξ definida positiva,
el criterio Dπ-optimalidad asociado a la distribución a priori π, χ un conjunto compacto, p el número de parámetros del modelo y
el diseño que asigna toda la probabilidad al punto
. Si la derivada direccional
en la dirección de
satisface:
y además se cumple que
para todo punto de soporte de ξ∗,
entonces el diseño ξ∗ es Dπ-óptimo.
2.4 D-eftciencia
Dados un diseño ξ y un diseño ξ∗ Dπ-óptimo, es posible compararlos mediante la eficiencia:
La eficiencia muestra qué tan próximo está el diseño ξ de la información dada por el diseño Dπ-óptimo ξ∗. Se dice que el diseño ξ es comparable con el diseño Dπ-óptimo si su eficiencia es cercana a 1.
3 Metodología
A continuación se describe la metodología propuesta para hallar los diseños óptimos heteroscedásticos para el modelo propuesto en la ecuación (1):
1. Se determinan los valores locales de los parámetros a partir de la estimación de los mismos por máxima verosimilitud. Se encuentran los intervalos de confianza vía boostrap.
2. Con el valor local hallado en el paso anterior, se halla el diseño óptimo que maximice el determinante de la matriz de información dada en la ecuación (3).
3. Se construyen las diferentes distribuciones a priori descritas a continuación:
Distribución uniforme discreta para cada parámetro: Se halla al perturbar por una cantidad δ, tanto a la derecha como a la izquierda, cada componente del vector de valores locales, es decir, para cada parámetro θi, 1 ≤ i ≤ r, se obtienen dos perturbaciones θi(1 - δ) y θi(1+δ). Luego, se consideran todas las posibles combinaciones de las perturbaciones obtenidas, generando un conjunto de 2r vectores para cada δ. A este conjunto se le adiciona el vector de parámetros local y se le asigna a cada uno, una probabilidad de 1/(2r + 1). De este modo, la distribución a priori a considerar dará la misma probabilidad a cada uno de los puntos del espacio parametral.
Distribución uniforme continua: Para cada θi se considera una distribución a priori continua definida en el respectivo intervalo de confianza hallado en el paso 1. La distribución uniforme continua conjunta, π, es el producto de las respectivas distribuciones a priori individuales para cada componente θi.
Distribución a priori de la familia de escala: Se construye la familia de escala para cada uno de los θi, de tal forma que cada distribución uniforme continua considerada en (b) esté incluida en esta familia.
4. Para cada una de las distribuciones a priori dadas en el paso 3, se halla el respectivo diseño Dπ-óptimo y se calcula la respectiva eficiencia. Lo anterior se hace implementando un programa en el lenguage estadístico R, [25].
4 Estudio de un caso
En esta sección, se presenta un ejemplo donde se considera un modelo no lineal heteroscedástico y se halla el diseño D π -óptimo local. Luego, a partir del vector de parámetros local se construyen distribuciones a priori para evitar la dependencia de éste y se comparan los diseños obtenidos con estas distribuciones con el diseño D π -óptimo local mediante la D π -eficiencia.
4.1 Descripción de los datos
La metodología anterior se bosqueja a partir de un caso de estudio con datos tomados de Bates & Watts (1988). Los datos corresponden a las concentraciones de residuos de policlorobifenilos (PCB) en una clase de trucha del lago Cayuga, medidas en 28 peces y se tiene el registro exacto de la edad (en meses) debido a su almacenamiento anual. Todos los peces se cortaron, se molieron, se mezclaron y se tomaron muestras de cinco gramos. Se trataron las muestras y se estimaron los residuos de PCB en partes por millón (ppm) usando cromatografía en columna.
En el gráfico de dispersión, ver Figura 1, se puede observar que a medida que aumenta la edad de los peces también aumenta la variabilidad, por lo que la varianza del error no es constante.
4.2 Ajuste del modelo
En la Figura (1) se puede observar que no hay una relación lineal entre las variables; por lo cual, se propone ajustar el siguiente modelo no lineal al conjunto de datos:
donde 
y 
son los parámetros a estimar, los errores tienen distribución normal con media cero y varianza constante y son independientes. Ahora, uno de los supuestos que se requieren para hallar los diseños óptimos es que se cumpla la homogeneidad de la varianza de los errores. Al realizar los gráficos de: residuales contra la edad de los peces, residuales contra los valores ajustados y el test de Levene para la homogeneidad de varianza, se puede concluir que la varianza del error no es constante por lo que se propone considerar la varianza como una función de la media. El modelo propuesto es el siguiente
donde los errores se suponen independientes con distribución normal (NI), media cero y varianza
Ahora, interesa hallar las estimaciones de los parámetros α 0, α 1, τ y σ. Para ésto, se utiliza el método de mínimos cuadrados generalizados, mediante la función gnls de R [25]. Las estimaciones puntuales obtenidas se muestran en la Tabla (1):
Luego, se usará como vector de parámetros local:
Con la prueba de razón de verosimilitudes se comprueba que, en efecto, el modelo con función de varianza dada en (8) es mejor que el modelo (6), cuando se supone varianza constante, ver resultados en la Tabla 2:
4.3 Diseño D-óptimo local
Con el valor local θ0, se halla el diseño D-óptimo local a partir del uso de dos funciones, la función DEoptim de R, [25], para hallar los valores iniciales y luego, con estos valores se utiliza la función nlminb de R [25] y se obtiene el diseño D-óptimo local siguiente:
Con el respectivo teorema de equivalencia se verifica que en efecto este diseño es D-óptimo local.
4.4 Diseños Dπ-óptimos
En esta sección se obtienen diseños Dπ-óptimos para el modelo inicial con tres y cuatro puntos de soporte y se comparan con el diseño D-óptimo local mediante el cálculo de la eficiencia. Se estudia la dependencia del diseño a la elección de la a priori, por medio de un estudio vía simulación, donde se toman diferentes distribuciones que pertenecen a la familia de escala de modo que también pertenezca la distribución a priori continua construída a partir del vector de parámetros local y se puede concluir que el diseño es robusto a la elección de la distribución a priori.
4.4.1 Distribuciones a priori discretas. Para este estudio se van a considerar distribuciones uniformes discretas definidas en los diferentes intervalos generados a partir de las componentes del vector de parámetros local, los cuales se definen a continuación:
Con el vector de parámetros dado en (9) se generan las distribuciones a priori discretas uniforme de la siguiente manera:
Se toman diferentes valores de δ entre 0.40 y 0.95 con incrementos de 0.05, denotado por δ = 0.40 : 0.95 : 0.05, generando 11 valores de δ, con el fin de obtener diferentes perturbaciones de θ0.
Se multiplica cada componente del vector de parámetros local por (1 - δ) y (1 + δ) y se consideran todas las combinaciones posibles de estos valores, generando así un conjunto de 24 vectores para cada δ.
Al conjunto del numeral anterior se le adiciona el vector de parámetros local, obteniendo un conjunto, Θδ, con 17 vectores. Luego, se le asigna una probabilidad 1/17 a cada vector. De lo anterior, se obtiene una distribución a priori discreta uniforme para cada δ:
Con las diferentes distribuciones a priori definidas anteriormente se hallan los diseños Dπ-óptimos mostrados en las Tablas (3) y (4). Se puede observar que para perturbaciones de δ de hasta del 0.65 se obtuvieron diseños con tres puntos de soporte, mientras que para perturbaciones del vector de parámetros local mayores al 0.65 se consiguieron diseños con cuatro puntos de soporte.
Ahora, mediante la Dπ-eficiencia, se comparan algunos de los diseños Dπ-óptimos obtenidos, que tan competitivos son los diferentes diseños hallados con relación a un diseño de referencia. En la Figura (2) se muestran las Dπ-eficiencias considerando como diseño de referencia el diseño óptimo local.
Se observa que la mayoría de las eficiencias son cercanas a uno (ver Figura 2]. Luego, los diseños obtenidos con las distribuciones a priori uniformes discretas con las diferentes perturbaciones son casi tan eficientes como el diseño obtenido con el vector de parámetros local.
Figura 2 Dπ-eficiencia de los diseños mostrados en las Tablas (3) y (4) para diferentes valores de δ
Tabla 5 Comparación de los diseños mediante la Dπ-eficiencia para los diseños obtenidos con una perturbación δ tomando como diseño de referencia el diseño Dπ-óptimo con perturbación δ0
La Tabla 5 muestra todas las Dπ-eficiencias obtenidas al comparar todos los diseños Dπ-óptimos dos a dos y con el diseño D-óptimo local. Se puede observar que al comparar el diseño D-óptimo local con el diseño Dπ-óptimo obtenido al tomar una distribución uniforme construida a partir de una perturbación del 95 % se obtuvo una Dπ-eficiencia del 97 %, por lo que se puede concluir que el diseño óptimo local es tan eficiente como el diseño Dπ-óptimo con δ = 0.95.
4.4.2 Distribuciones a priori continuas. En general, cuando se tienen distribuciones a priori continuas, en el proceso de maximización se debe evaluar la expresión dada en (5). Para ello, se utiliza simulación Monte Carlo como se indica a continuación: Se generan NS números aleatorios de la distribución a priori π(θ); por la ley débil de los grandes números, para la distribución π(θ) y para la función log |M (ξ, θ)|, se tiene que
donde la notación 
indica que la serie {a
n
} converge en probabilidad a b y
Es decir, 
converge en probabilidad al criterio Dπ-óptimo 
. El procedimiento anterior se aplicó para el caso de distribuciones uniforme continuas tomando NS = 1000, descritas más adelante.
Mediante simulaciones Bootstrap (Rana, Midi e Imon, [26]), se obtienen los intervalos de confianza para el vector de parámetros local dados en la Tabla (6):
Los números aleatorios de las distribuciones a priori uniformes continuas a partir del vector de parámetros local se generan de la siguiente manera:
1 Se generan números aleatorios de la distribución uniforme continua en (0, 1), vía simulación, mediante la función runif de R, [ 25 ].
2 Luego, por medio de la transformación gi : [0, 1] → [ai, bi] definida por
para i = 1, 2, 3 y 4, se generan números aleatorios de las cuatro distribuciones a priori continuas, 
, para cada una de las cuatro componentes del vector de parámetros local en los intervalos de confianza 
y 
.
3 Los números aleatorios generados en el numeral 2, serán los números aleatorios de la distribución a priori uniforme continua final definida como la distribución de probabilidad conjunta de las distribuciones definidas en el numeral anterior:
donde 
Con los números aleatorios generados anteriormente, se evalúa el criterio Dπ-optimalidad y se obtiene el diseño Dπ-óptimo siguiente:
El diseño anterior, sugiere que el 49 % de las corridas experimentales se debe realizar en peces de un año, aproximadamente el 49 % en peces de doce años y el 1.5 %, en peces de ocho años y cuatro meses aproximadamente.
En la Figura (3) se muestra la derivada direccional para el diseño obtenido con 1000 puntos. De esta gráfica, se puede observar que el diseño obtenido, es Dπ-óptimo, ya que la derivada es cero en los puntos de soporte y es estrictamente menor que cero en los demás puntos.
4.4.3 Estudio de robustez a la elección de la distribución a priori En esta sección se estudia la robustez a la selección de la distribución a priori vía simulación, donde se considera como diseño base el diseño Dπ-óptimo obtenido con el vector de parámetros local y una familia de escala donde las distribuciones se construyen a partir de una distribución uniforme continua en (0, 1), donde se fija el extremo izquierdo de los intervalos de confianza de las componentes del vector de parámetros local obtenidos en la sección anterior, ai, para i = 1, 2, 3 y 4 y se toma como parámetro de escala el valor (bi − ai)k, .
De la sección 4.4.2, el diseño Dπ-óptimo obtenido estaba dado por:
Ahora, sean , i = 1, 2, 3 y 4, la distribución uniforme continua en [0, 1] obtenida mediante simulaciones en R, con función de densidad de probabilidad fi(x), para cada una de las cuatro componentes del vector de parámetros local. Se realiza un procedimiento análogo al descrito en la sección 4.4.2 para construir estas distribuciones.
La familia de escala para el parámetro θi considerada en este caso, es la familia que tiene como parámetro de escala (bi -ai)k, donde k es un número real positivo:
con fi la función de densidad de probabilidad de la distribución uniforme . Se considera esta familia con el fin de que cada una de las distribuciones uniforme continuas, construidas en el paso 2 de la sección anterior, pertenezca a la familia, , lo cual se obtiene al tomar k = 1 y así poder compararla con otras distribuciones, que también pertenezcan a la familia, , obtenidas con diferentes valores de k.
Para la construcción de la familia de localización y escala, se siguen los siguientes pasos:
1. En la sección anterior se definieron los intervalos de confianza para , τ y σ, los cuales se denotaron por: y .
2. Luego, se generan 1000 números aleatorios de la distribución uniforme con la función runif de R, [ 25 ], para cada una de las cuatro componentes del vector de parámetros local. Denote por U i los valores aleatorios generados.
3. Sean Z i = a i + k(b i - a i )U i , las variables aleatorias construidas a partir de las variables U i . Entonces Z i tiene f.d.p.
con f (·) la función de densidad de la distribución .
4. Finalmente, de manera análoga a como se hallaron números aleatorios de las distribuciones a priori continuas, se generan números aleatorios de la distribución a priori uniforme continua final definida como la distribución de probabilidad conjunta de las distribuciones definidas en el numeral anterior:
donde
Después de generar los números aleatorios de las distribuciones a priori uniformes continuas, y a partir de un programa elaborado en R, [25], se obtienen los diseños D π -óptimos, mostrados en la Tabla 7:
Tabla 8 D π -eficiencia de los diseños obtenidos de la familia de escala con respecto al diseño D π -óptimo.
En la Tabla (8) se muestran las D π -eficiencias que resultan de comparar los diseños obtenidos de la familia de escala, con k = 0.9, 1.1, 1.2, 1.3 con respecto al diseño D π -óptimo. Este último, se halla, con la respectiva a priori que pertenece a esta misma familia tomando k = 1. De los resultados arrojados, se puede notar que todas las D π -eficiencias están por encima del 98 %, permitiendo concluir que los diseños obtenidos con las diferentes distribuciones consideradas pertenecientes a la familia son casi tan eficientes como el diseño obtenido con la distribución de referencia, es decir, cuando k = 1. Lo anterior, permite afirmar que el diseño no depende de la elección de la distribución a priori, es decir, el diseño es robusto a la elección de la a priori, al menos si se selecciona la distribución a priori entre las distribuciones investigadas.
5 Conclusiones
En este artículo se extendió el teorema de equivalencia de Kiefer y Wolfowitz [22] al caso heteroscedástico cuando se considera en el criterio de optimalidad información a priori de los parámetros del modelo. Se propuso una metodología para la construcción de distribuciones a priori tanto continuas como discretas, las cuales, a partir de un ejemplo, permitieron incrementar el número de puntos experimentales del diseño de dos hasta cuatro puntos, lo cual es beneficioso para el usuario cuando se necesitan hacer pruebas de bondad de ajuste para el modelo que se está considerando. Se mostró que los nuevos diseños encontrados son tan competitivos como el diseño óptimo local, que tiene un reducido número de puntos, con la ventaja que los diseños construidos tienen un mayor número de puntos experimentales. Se exhibió una familia de distribuciones de escala en la cual los diseños óptimos obtenidos son robustos a la elección de la distribución a priori, con respecto al cálculo de la eficiencia.
