Un estudio
de Monte Carlo sobre la
estabilidad
y el sesgo de los
coeficientes
en el análisis discriminante
Juan Camacho Rosales
Introducción
En la
interpretación del análisis discriminante y el análisis de varianza
multivariado se pueden utilizarlos coeficientes típicos discriminantes y los
coeficientes estructura (Bray y Maxwell,1982; Camacho, 1990; Huberty,1984;
Share,1984). Los coeficientes típicos indican la contribución relativa de cada
variable a la función discriminante y los coeficientes estructura indican la
proporción de varianza que comparte cada variable con la función discriminante.
Si bien ambos coeficientes tienen una interpretación diferente, existen
criterios estadísticos para juzgar la bondad de estos coeficientes. Estos
criterios estadísticos se basan en el conocimiento de las propiedades
muestrales de los coeficientes.
El
propósito del presente estudio es estudiar algunas características muestrales
de estos dos coeficientes y de un tercer tipo de coeficiente propuesto por
Huberty, (1984) como otro indicador de la "importancia" de las
variables. Este coeficiente es análogo al coeficiente de correlación
semiparcial al cuadrado en la Regresión Múltiple. Darlington (1968), al
discutir indicadores de la "importancia" de las variables en la
Regresión Múltiple, lo denominó utilidad
de la variable e indica la contribución única de cada variable al
coeficiente de correlación múltiple, es decir, indica lo que disminuye el
coeficiente de correlación múltiple si se elimina esa variable del análisis.
Thomdikey
Weiss (1973) estudiaron la estabilidad de los coeficientes típicos y los coeficientes
estructura en el marco de la correlación canónica. Utilizando validación
cruzada de datos reales concluyeron que los coeficientes estructura eran más
estables y útiles que los coeficientes típicos. Sin embargo, Barcikowski y
Stevens (1975) en otro estudio con la misma finalidad que el anterior no
hallaron una clara superioridad de ningún tipo de coeficiente. Estos autores
hicieron un estudio de Monte Carlo utilizando las matrices de varianzas y
covarianzas (intra grupo) de ocho estudios obtenidos de la literatura
psicológica, con un número de variables que iba de 7 a 41. Utilizaron como
criterio de la estabilidad de los coeficientes el valor del coeficiente de concordancia
W de Kendall y concluyeron que se necesitaba una relación de 42/1 a 68/1 entre
número de sujetos y número de variables para encontrar resultados fiables.
Thorndike (1976) critica el estudio anterior en base al uso de una medida (W de
Kendall) que se basa en rangos y que, por lo tanto, puede omitir importantes
diferencias. Argumenta que la información interesante para la comparación entre
los coeficientes se ha de basar en las características muestrales de los coeficientes.
Además Thorndike critica el estudio anterior en base a que no se hace ninguna
manipulación sistemática de las condiciones del estudio y propone nuevas
dimensiones para la comparación de los coeficientes: la razón entre el tamaño
de la muestra y el número de variables y la correlación media (intragrupo)
entre variables. Barcikowski y Stevens (1976), finalmente, señalan que la
importancia relativa de cada una de estas dimensiones para causar inestabilidad
podría ser importante para determinar situaciones en las que una razón menor
entre el número de sujetos y el número de variables podría ser suficiente para
una interpretación exacta.
Varios
autores han realizado estudios comparativos de los coeficientes dentro del marco
del análisis discriminante o del Manova (Borgen y Seling, 1978; Camacho, 1990;
Spector, 1977; Wilkinson,1977) pero estos estudios se han limitado a analizar
el efecto de una o dos replicaciones. Sin embargo, ponen de relieve una
dimensión de estudio propia del análisis discriminante: el número de grupos.
Para la simulación del análisis discriminante, donde se resuelve una ecuación
del tipo A v = k B v, B se manipula mediante el nivel de correlación media
intragrupo y A mediante el perfil de medias de los grupos.
Debido al
gran número de dimensiones implicadas, el estudio presente se limitará a estudiar
algunas situaciones básicas. La dimensión tamaño de la muestra partido por el
número de variables se estudiará mediante la manipulación del número de
variables y del número de sujetos por grupo. La correlación media entre
variables se mantendrá constante. Finalmente, el número de grupos se
manipulará, pero de forma que sólo exista una dimensión subyacente, es decir,
todos los grupos tendrán igual perfil. Como criterios estadísticos se
utilizarán el sesgo (diferencia entre el parámetro y la media de las
estimaciones) y la estabilidad (desviación típica de las estimaciones). Se
estudiará el sesgo y la estabilidad de los coeficientes típicos, los
coeficientes estructura y los coeficientes de correlación semiparcial al
cuadrado. Esto permitirá comparar las características muestrales de los
coeficientes en diferentes situaciones experimentales.
Método (Diseño del estudio de Monte Carlo)
Modelos
Los modelos unidimensionales se
construyeron a partir de las siguientes variables experimentales:
1)el número de variables: 3, 6 ó 9;
2) el número de grupos: 3 ó 6;
3) el tamaño de la muestra en cada grupo: 25, 50, 100 ó
200.
Estas variables experimentales
se agruparon en seis modelos de cuatro tamaños cada uno:
|
Modelo |
A |
B |
C |
D |
E |
F |
|
Número de variables |
3 |
6 |
9 |
3 |
6 |
9 |
|
Número de grupos |
3 |
3 |
3 |
6 |
6 |
6 |
La
proporción de la varianza común (p) entre la variable i y la variable j se
mantuvo constante e igual a 0.50 (correlación media igual a 0.70). Por lo
tanto, la covarianza entre la variable i y la variable j es igual a 283
√(S1 * S2 * p) = √
(400 * 400 * 0,50)
ya que en la matriz B la
desviación típica se mantuvo constante para todas las variables, y fue igual a
20. El perfil de diferencia entre grupos para el caso de tres grupos era X1i =
0, X2i =10 y X3i, = 20, donde Xli es
la media del grupo 1 en la variable i. Esto producía una matriz A de rango 1.
En la matriz de sumas de cuadrados y productos cruzados (SSCP) debida a la
diferencia entre los grupos, todos los elementos
eran
igual a 5000 (para n=25). En el caso de 6 grupos se eligió el perfil de manera
que los elementos de la matriz de sumas de cuadrados y productos cruzados fuese
equivalente, es decir, igual a 10000 (para n=25). El perfil para cada variable
fue
Xl1 = 0, X2i =
4, X3i = 9, X4i = 13, X5i = 19 y X6i = 23. (En el presente estudio
fue exactamente de 9737.5 para n=25).
Organización del Estudio de Monte Carlo
A continuación se describe paso
a paso como se llevó a cabo el estudio. Estos pasos son una descripción de dos
programas FORTRAN ideados para el estudio (los programas se pueden obtener a
través del autor).
1. En primer lugar, se hallaron las matrices SSCP A y B, para el modelo
A y n=25. En general, el elemento b8ij, de B8 es igual a
bij * k * (n -1) donde k es el número de grupo, n el número de
sujetos por grupo y bij es el elemento correspondiente de la matriz B (de
varianzas y covarianzas intragrupo). Y el elemento a8ij de A8
es igual a aij * n. El elemento bij se mantuvo igual a
283, el elemento bij iguala 400 y el elemento aij igual a
200. Por ejemplo, para k=3 y n=25, b8ij, = 283 * (3*(25 ‑ 1))
= 20942 y a8ij= 200 * (25) = 5000.
2. Se
resolvió el sistema de ecuaciones As v = k Be v, mediante
la subrutina EIGZS de IMSL y se hallaron los parámetros mediante los siguientes
algoritmos: coeficientes típicos Dt = √(Bii) V,
donde Bii es la diagonal de la matriz B y V la matriz de valores
propios (por columnas); los coeficientes estructura Sw = √(Bii))-1
V; y los coeficientes semiparciales P2ri = d2i
/ rii, donde d2i es el coeficiente típico de
la variable i elevado al cuadrado y rii es el elemento i de la
diagonal de R-1w, donde Rw = (√(Bii))
-1 B (√(Bii)) -1 es la matriz de
correlaciones intragrupo.
3. Cada modelo tenia una matriz B que definía la matriz de varianzas y
covarianzas de la distribución normal multivariada, de la que se obtenían 300
muestras. Estas muestras se obtuvieron también para cada tamaño de grupo. Se
utilizó la subrutina GGNSM de IMSL para el muestreo. Esta subrutina se basa en
el cálculo de una cierta matriz C tal que B = C C' y la matriz de puntuaciones
se ha¡ la mediante X = Z C, donde X es la matriz de puntuaciones con matriz de
var‑covar B y Z es la matriz de puntuaciones normales (0,I) (Sheuer,
1962).
4. Utilizando la subrutina GGNSM
de IMSL se obtuvieron k (3 6 6) matrices X. N(mu,B), (en este paso el subíndice
i indica el grupo), de las que se hallaron B.i y mui, donde B8i es
la SSCP intragrupo y mui es el vector de medias. A continuación se halló Bs
= ΣB8i , y As (la SSCP intergrupos) se obtenía a
partir de la suma de la matriz [mul, mu2, ...,
mukj a la matriz del perfil de medias característico del modelo. Por
ejemplo, para k=3 y 3 variables la matriz del perfil era:
|
Variable |
1 |
2 |
3 |
|
Grupo 1 |
0 |
0 |
0 |
|
Grupo 2 |
10 |
10 |
10 |
|
Grupo 3 |
20 |
20 |
20 |
5. Mediante los algoritmos del
paso 2 se resolvía el sistema de ecuaciones As v = k Bs v y se hallaban los
coeficientes típicos, los coeficientes estructura y los coeficientes
semiparciales. (Coeficientes típicos y coeficientes estructura en valor
absoluto).
6. Se repitieron 300 veces los
pasos 4 y 5. La subrutina GGNSM está diseñada para asegurar resultados
independientes.
7. Se repitieron los pasos 1 a 6
para el resto de los modelos. En total se hallaron 7200 (300 x 3 x 2 x 4)
análisis discriminantes.
Para probar la estabilidad del procedimiento de
generación aleatoria se obtuvieron los datos del modelo A con diversas semillas
aleatorias y con 1000 repeticiones en vez de 300 y en ambos casos los
resultados fueron similares. Una vez obtenidas las trescientas estimaciones de
cada coeficiente en cada modelo y tamaño de muestra se calculó la media, la desviación
típica y el sesgo de las distribuciones de los coeficientes. El sesgo era la
diferencia entre la media de las 300 estimaciones y el parámetro.
Resultados
Para conseguir una
representación concisa de los resultados se halló la media de los sesgos en
valor absoluto para cada grupo de coeficientes, se restó el valor del parámetro
y se multiplicó por cien el resultado. Los datos (Ver Tabla 1) se presentan sin
decimales, de manera que una diferencia en valor absoluto menor que 0.005 se
iguala a cero.
Tabla
1.100 veces la diferencia entre la media de las estimarionesy el
correspondiente parámetro para los coeficientes típicos (Dd, los coeficientes
estructura (SW) y los
coeficientes semiparciales (P,) para todos los modelos y tamaños de muestra.
|
Tamaño de los Grupos |
Modelo
A Dt Sw
Pr |
Modelo
B Dt Sw
Pr |
Modelo
C Dt Sw
Pr |
Modelo
D Dt Sw
Pr |
Modelo
E Dt Sw
Pr |
Modelo
F Dt Sw
Pr |
|
25 |
9 8 7 |
20 16 6 |
26 22 5 |
4 4 4 |
12 9 4 |
17 14 3 |
|
50 |
3 4 3 |
11 8 3 |
1512 2 |
2 2 2 |
5 4 2 |
9 6 1 |
|
100 |
0 2 2 |
5 4 2 |
8 6 1 |
1 1 1 |
2 2 1 |
4 3 0 |
|
200 |
1 1 1 |
2 2 1 |
4 3 0 |
0 0 0 |
8 1 0 |
2 2 8 |
Se
produce un sesgo creciente a medida que aumenta el número de variables, tanto
en los coeficientes típicos como en los coeficientes estructura, aunque el
sesgo de estos últimos es ligeramente menor. Sin embargo, los coeficientes
semiparciales muestran la tendencia contraria, el sesgo se hace menor a medida
que intervienen más variables.
Comparando
ahora el sesgo dependiendo del número de grupos, se observa una clara tendencia
de todos los coeficientes a disminuir el sesgo en los modelos con seis grupos.
El
aumento del número de sujetos por grupo tiene también un claro efecto reductor
sobre el sesgo, por ejemplo, en el modelo C, el sesgo de los coeficientes
típicos disminuye de 26 con n=25 a 4 con n=200.
Barcikowski
y Stevens (1975) propusieron una razón entre el número de sujetos y el número
de variables de 42/1 a 68/1 para encontrar resultados fiables. Si en el
presente estudio consideramos como buenos los resultados con sesgo menor o
igual a 5 se puede estudiar la bondad de la anteriorrazón. La linea insertada
en la Tabla 1(qua tambien se puede insertar en la Tabla 2) limita las
situaciones en las que la razón es menor de 50/1 (zona superior, exactamente
33/1 en el caso mayor) o igual o mayor a 50/1 (zona inferior). Se puede
observar que para los coeficientes típicos y los coeficientes estructura se
obtuvieron sesgos menores a 5 cuando la razón es igual o mayor a 50/1. Los
coeficientes semiparciales, sin embargo, muestran buenos resultados incluso en
el caso de una razón igual a 8/1, en el caso del modelo C. Estos resultados
obtenidos por los coeficientes semiparciales se pueden deber a lo extremo de
los parámetros que se tienen que estimar. Los parámetros son 0.057, 0.013 y
0.005 (para los modelos A, B y C, respectivamente), mientras que en el caso de
los coeficientes típicos los parámetros son 0.372, 0.191 y 0.129, y en el caso
de los coeficientes estructura son 0.897, 0.870 y 0.860. Los modelos D, E y F
tenían los mismos parámetros que los modelos A, B y C, respectivamente. Por
otro lado, al ser los coeficientes semiparciales cantidades elevadas al
cuadrado se reduce la magnitud de los valores y, por lo tanto, de sus sesgos.
De manera que se puede concluir que todos los coeficientes muestran resultados
poco sesgados a partí r de una razón de 50/1 entre el número de sujetos y el
número de variables.
Tabla 2. Desviaciones típicas medias de los
coeficientes típicos (Dt), los coeficientes estructura (Sw)
y coeficientes semilos (Pr) para todos los modelos y tamaños de
muestra. (Todos los valores son decimales: 33 es 0’33, 25 es 0’25, etc.)
|
Tamaño de los Grupos |
Modelo
A Dt Sw
Pr |
Modelo
B Dt Sw
Pr |
Modelo
C Dt Sw
Pr |
Modelo
D Dt Sw
Pr |
Modelo
E Dt Sw
Pr |
Modelo
F Dt Sw
Pr |
|
25 |
33 15
15 |
29 1610 |
28 16
08 |
26 1010 |
23 11 07 |
23 12
06 |
|
50 |
25 09
10 |
22 10
06 |
21 10
05 |
20 06
06 |
17 07
04 |
16 07
03 |
|
100 |
10 06
07 |
17 07
04 |
16 07
03 |
16 05
05 |
14 05
03 |
12 05
02 |
|
200 |
15 05
05 |
13 05
02 |
12 05
02 |
11 03
04 |
11 03
02 |
10 04
01 |
En la Tabla 2 se pueden ver las
desviaciones típicas medias de los diversos tipos de coeficientes. Se halló la
desviación típica media para cada grupo de coeficientes. La variabilidad de las
estim aciones disminuye cuando aumenta el número de variables, el número de
grupos o el tamaño de los grupos. Sin embargo, la variabilidad no muestra una
tendencia similar a la del sesgo en cuanto a la razón entre el número de
sujetos y el número de variables. Ello se debe en parte a la clara diferencia
entre las desviaciones típicas de los diversos coeficientes. Cuando la razón es
50/1 o superic la desviación típica de los coeficientes estructura y de los
coeficientes semiparcial se mantiene por debajo de 0.10. Sin embargo, en los
coeficientes típicos la desviacic típica en esta zona llega a 0.26 para el
modelo D. Los coeficientes típicos tienen w variabilidad mayor que los otros
dos tipos de coeficientes. Los coeficientes semiparcial tienen menor
variabilidad que los coeficientes estructura aunque en los modelos A D son
iguales. Si se considera el error cuadrático medio (desviación típica
asintótica) como indice de la eficiencia de las estimaciones, los coeficientes
estructura y los coeficientes semiparciales son más eficientes que los
coeficientes típicos. El en cuadrático medio es igual a la varianza de las
estimaciones más el sesgo elevado cuadrado.
De los
tres tipos de coeficientes, los coeficientes seiuiparciales son los que
presentan menor sesgo. El sesgo de los coeficientes típicos y de los
coeficiente estructura es similar, ligeramente menor el de estos últimos. También
son k coeficientes semiparciales los que menor variabilidad presentan, los
coeficienU típicos son los que presentan mayor dispersión y con los modelos de
tres variables presenta una igualación entre las desviaciones típicas de los
coeficientes semiparciales y los coeficientes estructura. Aunque estos
resultados admiten una generalizacit limitada, debido a la limitación en el
número de situaciones estudiadas, se pues concluir, en función del sesgo y la
variabilidad de los coeficientes, que el orden de bondad estadística de mejor a
peor es: coeficientes semiparciales, coeficienn estructura y coeficientes
típicos.
La razón
de Barcikowski y Stevens usada en correlación canónica se puede utilizar en
análisis discriminante de forma que la razón entre el número de variables y
número de sujetos sea de 50/1 ó más. En la presente investigación el modelo
utilizaba 6 grupos de 25 sujetos, de manera que la razón era 6 x 25 / 3 = 50,
igual razón que para el modelo A con n=50.
Para
futuras investigaciones, en cuanto a la comparación de las caracterlsticas
estadísticas de los coeficientes, hay una serie de condiciones de estudio que
necesitarían ser investigadas. En el modelo de análisis discriminante
representado por la ecuación A v =k B v, intervienen las matrices A (intergrupo)
y B (intragrupo). L. condiciones de estudio futuras para la simulación de la
matriz B serian: 1) correlacit media entre variables (en el presente estudio se
mantuvo constante); 2) homocedasáicida 3) igualdad de las matrices que producen
B; 4) dimensionalidad de la matriz B. L; condiciones de estudio para la matriz
A: 1) dimensionalidad de la matriz (en el presente estudio el rango era uno,
pero manipulando los perfiles se pueden obtener rangos mayores); 2) tamaño del
efecto producido por cada variable ‑valor F individual(con esta condición
se podría controlar explícitamente el poder separador de cada variable,
haciendo, por ejemplo, que variables importantes sean las que tienen una F
mayor).
Existe otra posibilidad, para
estudios posteriores, que consistiría en obtener las matrices A y B de estudios
reales y emplear los coeficientes obtenidos de estas matrices como parámetros.
La simulación se basaría en la generación aleatoria de muestras normales
utilizando la matriz B como modelo.
(1) El coeficiente es u2ji
=b2ji / rii, donde rii es el
elemento i de la diagonal de la inversa de la matriz de correlaciones
intragrupo, b2ji es el coeficiente típico de la variable
i en la función j, elevado al cuadrado. El coeficiente u2ji
tiene en cuenta el efecto de la colinearidad ya que rii = [1‑R2ii],
y si R2ii es el coeficiente de correlación múltiple entre
la variable i y el resto de las variables, u2ji es una
función inversa de la colinearidad.
BARCIKOWSKI, R. S. & STEVENS,
J. P (1975): A Monte‑Carlo Study of the Stability oí Canonical
Correlations, Canonical Weights and Canonical Varíate‑Variable
Correlations. Multivariate Behavioral
Research,10, 353‑365.
BARCIKOWSKI, R. S. & STEVENS,
J. P (1976): Studying Canonical Analysis: A Reply to Thbomdike's. Comment.
Multivariate Behavioral Research,11,
255‑258.
BORGEN, F. H. & SELING, M. J.
(1978): Uses oí Discriminant Analysis Following MANOVA: MultivariateStatistics
for Multivariate Purposes. lournal ofApplied
Psychology, 63, 689‑697. BRAY, J. H. & MAXWELL, S. E. (1982):
Analyzing and Interpreting Significant MANOVA. Review ofEducational Research,
52, 340‑367.
CAMACHO, J. (1990): Interpretación del Manova: Análisis de la
Importancia de las Variables. Dependientes.Qurriculum,1,
107‑120.
HUBERTY, C. J. (1984): Issues in
the Use and Interpretation of Discriminant Analysis. Psychological. Bulletin, 95, 156‑171.
IMSL Subroutines, (1981): IMSL,
Inc, Houston, Texas.
SHARE, D. L. (1984): Interpreting
the Output of Multivariate Analysis: A Discussion oí Current. Approaches.British
lournal of Psychology, 75, 349‑362.
SHEUER, E. M. (1962): On the
Generation of Normal Random Vectors. Technometrics,
4, 278‑281.
SPECTOR, P E. (1977): What to Do
With Significant Multivariate Effects in MANOVA. lournal ofApplied Psychology,62,158‑163.
THORNDIKE, R. B. M. (1976)
Studying Canonical Analysis: Comments on Barcikowski and Stevens. Multivariate Behavioral Research,11, 249‑253. '
THORNDIKE, R. B. M. & WEISS,
D. J. (1973): A Study‑of the Stability of Canonical Correlations and
Canonical Components. Educational and Psychological Measurement, 33, 123‑134.
WILKINSON, L. (1975)
Response Variable Hypothesis in the Multivariate Analysis of Variance. Psychological
Bulletin, 82, 408‑412.