El denominador de la varianza

Llevo una semana impartiendo el curso de Introducción a la Estadística a los participantes en el MBA de la ESTE. El programa de este curso consiste en una presentación de la Estadística Descriptiva univariante con una breve incursión en la Estadística Bivariante de la mano de las tablas de contingencia y del análisis de correlación y de regresión.

Hoy le ha tocado el turno a la Estadística Descriptiva para una variable cuantitativa. En concreto, me tocaba presentar el concepto de varianza y después de media hora de teoría se me ha ocurrido que sería una buena cosa utilizar el ordenador para hacer unos cuantos números.

Como ya he comentado en algún post anterior he decidido utilizar R para explicar Estadística. La cosa no parecía demasiado complicada: “Necesitáis tener los datos referidos a vuestra variable en un vector numérico y ahora, debéis emplear la función var() para obtener la varianza.” les he dicho.

Entonces un alumno me dice: “Oye…este resultado no coincide con el que nos ha dado cuando hemos hecho los calculos a mano”. “¿Cómo que no?”, le he dicho, “eso es imposible.”… Y de repente me he acordado del dichoso denominador de la varianza.

En efecto, la varianza, que se define como el momento de segundo orden de una variable aleatoria respecto a su valor esperado, tiene como denominador el número de observaciones de la variable. Por tanto, si disponemos de un vector numérico en R de nombre dist -que recoge las distancias de frenado de 50 vehículos-, los siguientes comandos de R deberían ser equivalentes:

> sum((dist-mean(dist))^2)/(length(dist))

> var(dist)

Los resultados, sin embargo, son distintos:

> sum((dist-mean(dist))^2)/(length(dist))
[1] 650.7796
> var(dist)
[1] 664.0608

¿Por qué?

La razón es muy simple: la función var() de R -al igual que la función VAR() de Excel- considera que el colectivo que se está analizando es una muestra, y lo que es más importante, que lo que se desea es utilizar el valor del estadístico muestral para estimar el parámetro varianza poblacional. Como resulta que la varianza muestral -empleando como denominador el número de observaciones- no es un estimador insesgado de la varianza poblacional pero la cuasivarianza muestral -empleando como denominador el número de observaciones menos uno, es decir, los grados de libertad- sí lo es, estos paquetes estadísticos optan por emplear como denominador de la varianza el valor n-1 sin avisarnos de ello de manera explícita. Así, los dos siguientes comandos sí que son equivalentes:

> sum((dist-mean(dist))^2)/(length(dist)-1)
> var(dist)

Pero cuidado, lo que nos proporcionan no es la varianza del conjunto de datos sino la cuasivarianza muestral que deberá ser empleada en una situación de inferencia para estimar -de manera insesgada- la varianza poblacional.

En resumen:

• Si lo que deseamos es calcular la varianza de un conjunto de observaciones -sea este conjunto una población, una muestra o cualquier otra cosa- sin ningún ánimo inferencial debemos emplear como denominador el número de observaciones.
• Si, por el contrario, contamos con una muestra aleatoria seleccionada de una población y deseamos estimar de manera insesgada el verdadero valor de la varianza poblacional a partir de la información de la muestra entonces debemos emplear la cuasivarianza muestral, es decir, debemos utilizar el denominador n-1 en el cálculo de este famoso índice de dispersión.

La Estadística es cosa de locos

Hoy he asistido a una defensa de tesis doctoral en mi facultad. El doctorando -ahora ya doctor- es un alto cargo de una entidad financiera y todos los miembros del tribunal han alabado -en mi opinión muy merecidamente- su capacidad para conciliar sus elevadas responsabilidades laborales con sus actividades de investigación.

Lo que me ha gustado mucho menos ha sido el comentario de dos de los miembros del tribunal al señalar como muy loable el esfuerzo que le debía de haber supuesto al doctorando dejar de lado los modos y maneras propios de la empresa y zambullirse en las costumbres propias de los investigadores y, en particular, en esas complicadas técnicas estadísticas tan alejadas de la práctica habitual de la alta dirección de las entidades financieras.

En mi opinión -y estoy seguro de que es el caso del nuevo doctor- los altos directivos de las empresas -en general- y de las entidades financieras -en particular- deberían ser personas dotadas de elevadas habilidades y competencias en el análisis de datos. Me parece además muy grave que estas opiniones sean vertidas en un tribunal de tesis doctoral por académicos que deberían estar de acuerdo con esta postura. Quizás todo esto sea tan solo un paso más en el arrinconamiento y desprecio que los métodos cuantitativos vienen sufriendo en los últimos años -y si no, a los nuevos planes de estudio de LADE ADE me remito-.

¡Cuidado!

Cuando mi hija -hoy de tres años- tenía dos se lanzaba como loca a por cualquier cuchillo que encontrara encima de la mesa. Cuando eso ocurría alguien le decía “ez duzu labana hartu behar oso arriskutsua da eta” y poco a poco hemos ido consiguiendo que abandone esa peligrosa costumbre. Ahora, cuando se sienta a la mesa se dedica a recoger los cuchillos de todo el mundo y dejarlos a un lado diciendo que son muy peligrosos -tendremos que ver cómo recuperamos esta situación-.

El párrafo anterior puede resultar chocante en un blog dedicado al análisis de datos cuantitativos pero es fácil anticipar por dónde van los tiros: no cabe duda de que el acceso a los paquetes de análisis de datos es muy deseable ya que pone al alcance de cualquiera la posibilidad de analizar datos y obtener de ellos información sobre la que tomar decisiones. Operaciones que hace diez años eran impensables se hacen hoy en día de forma rutinaria y esto, sin duda, es muy bueno. Sin embargo, all that glitters is not gold, y esta democratización del software tiene también un efecto perverso: semejante potencia de cálculo en manos inexpertas es como un arma en manos de un niño: ¡peligro!.

Esto es lo que ocurrió la semana pasada: un compañero de facultad me pidió que echara un vistazo a ciertos resultados que le había enviado por correo un doctorando cuya tesis dirige. El doctorando había elaborado un cuestionario y lo había administrado a un conjunto de 188 encuestados. Su interés radicaba en evaluar si el cuestionario era fiable y válido.

Es bien conocido que uno de los índices más comúnmente empleados para evaluar la fiabilidad -mejor dicho, la consistencia interna- de un instrumento diseñado para medir un concepto es el alfa de Cronbach. Este índice se calcula de manera muy simple y muchos paquetes estadísticos -por ejemplo, SPSS- lo incluyen entre sus opciones.

Pues bien, el doctorando, empeñado en evaluar la fiabilidad de su cuestionario y orgulloso poseedor de una licencia de SPSS seleccionó Analizar / Fiabilidad / Alpha de Cronbach en su flamante SPSS v16 y realizó un análisis de fiabilidad de todo el cuestionario. El error es de dimensiones galácticas:

• En primer lugar, un análisis de consistencia interna debe realizarse de modo individual para cada uno de los diferentes conceptos o constructos que se trata de medir con el cuestionario. No tiene sentido un análisis global de la fiabilidad cuando se miden conceptos tan distintos como notoriedad de una marca de vehículos, actitud hacia la publicidad o valoración de ciertos atributos del vehículo.
• En segundo lugar, y aún mucho más grave, es el hecho de incluir entre las variables sujetas al análisis de fiabilidad algunas de naturaleza cualitativa tales como el sexo del encuestado, su nivel de estudios o la marca de vehículo adquirida. Resultaba muy curioso ver en los resultados como ¡la marca de vehículo promedio era 13,27 con una desviación estándar de 3,18!

¡Si es que hay que tener cuidado con el SPSS, que lo carga el diablo!

Estadísticas del blog

¿Será que mañana es el examen?

Selección de variables explicativas en la regresión

Uno de las cuestiones más importantes a la hora de encontrar el modelo de ajuste más adecuado para explicar la variabilidad de una característica cuantitativa es la correcta especificación del llamado modelo teórico. En otras palabras, debemos seleccionar de entre todas las variables candidatas a ser explicativas de la variable dependiente un subconjunto que resulte suficientemente explicativo -lo que podemos medir, por ejemplo, mediante el coeficiente de determinación- y también no demasiado complejo -es decir, con muchas variables explicativas-.

En la práctica, no obstante, la selección del subconjunto de variables explicativas de los modelos de regresión se deja en manos de procedimientos más o menos automáticos. Los procedimientos más usuales son los siguientes:

• Método backward: se comienza por considerar incluidas en el modelo teórico a todas las variables disponibles y se van eliminando del modelo de una en una según su capacidad explicativa. En concreto, la primera variable que se elimina es aquella que presenta un menor coeficiente de correlación parcial con la variable dependiente-o lo que es equivalente, un menor valor del estadístico t- y así sucesivamente hasta llegar a una situación en la que la eliminación de una variable más suponga un descenso demasiado acusado en el coeficiente de determinación.
• Método forward: se comienza por un modelo que no contiene ninguna variable explicativa y se añade como primera de ellas a la que presente un mayor coeficiente de correlación -en valor absoluto- con la variable dependiente. En los pasos sucesivos se va incorporando al modelo aquella variable que presenta un mayor coeficiente de correlación parcial con la variable dependiente dadas las independientes ya incluidas en el modelo. El procedimiento se detiene cuando el incremento en el coeficiente de determinación debido a la inclusión de una nueva variable explicativa en el modelo ya no es importante.
• Método stepwise: es uno de los más empleados y consiste en una combinación de los dos anteriores. En el primer paso se procede como en el método forward pero a diferencia de éste en el que cuando una variable entra en el modelo ya no vuelve a salir, en el procedimiento stepwise es posible que la inclusión de una nueva variable haga que otra que ya estaba en el modelo resulte redundante y sea “expulsada” de él.

El modelo de ajuste al que se llega partiendo del mismo conjunto de variables explicativas es distinto según cuál sea el método de selección de variables elegido. La consecuencia de este hecho resulta obvia: ninguno de los llamados métodos automáticos garantiza encontrar el modelo óptimo -en el sentido, por ejemplo de maximizar el coeficiente de determinación o cualquier otro criterio que nos parezca relevante-.

En efecto, consideremos los siguientes datos -lógicamente ficticios y preparados para mostrar el efecto que se desea- . Se trata de un conjunto de 40 observaciones de tres variables a las que llamaremos Y (la dependiente) y X1 y X2 (las explicativas). Si sobre estos datos aplicamos un procedimiento de tipo forward o de tipo stepwise el resultado será el siguiente:

a Variable dependiente: Y

Es decir, el modelo óptimo según los procedimientos forward o stepwise es el que no contiene ninguna variable explicativa. Esto es debido a que ninguna de las dos variables -cuando son consideradas de manera independiente- supera los criterios mínimos para ser incluida en el modelo -en concreto, que su coeficiente t lleve asociada una probabilidad crítica inferior a 0,05-. En efecto, los resultados para los modelos independientes son los siguientes:

Modelo		Coeficientes no estandarizados		Coeficientes estandarizados	t	Sig.
		B	Error típ.	Beta
1	(Constante)	15,504	1,337		11,596	,000
	X1	,026	,018	,221	1,399	,170
a Variable dependiente: Y

Modelo		Coeficientes no estandarizados		Coeficientes estandarizados	t	Sig.
		B	Error típ.	Beta
1	(Constante)	16,675	1,091		15,284	,000
	X2	,746	1,543	,078	,483	,632
a Variable dependiente: Y

Si por el contrario consideramos un procedimiento backward el resultado será el siguiente:

Modelo		Coeficientes no estandarizados		Coeficientes estandarizados	t	Sig.
		B	Error típ.	Beta
1	(Constante)	1,300	,317		4,097	,000
	X2	-63,555	1,230	-6,662	-51,688	,000
	X1	,792	,015	6,812	52,856	,000
a Variable dependiente: Y

El coeficiente de determinación para este modelo con dos variables explicativas es 0,987 y al coeficiente F asociado le corresponde una probabilidad crítica inferior a 0,001. Adicionalmente, a los estadísticos t asociados a cada una de las dos variables explicativas les corresponden probabilidades críticas muy reducidas. Hemos encontrado, por tanto, un buen modelo lineal para explicar el comportamiento de Y a partir del comportamiento de X1 y X2. El problema radica en que si hubieramos elegido de forma acrítica utilizar un procedimiento forward o stepwise, jamás lo habríamos encontrado.

Existen otras opciones de selección de subconjuntos de variables -como por ejemplo el llamado método fuwil o de selección de todos los posibles subconjuntos-. El inconveniente de este último es la explosión combinatoria que se produce cuando el número de variables candidatas a ser explicativas crece.

La paradoja de Simpson (II)

Las frecuencias observadas en una tabla de contingencia -como la que se mostraba en el post anterior relativa a la campaña de la DGT- pueden ser modelizadas adecuadamente mediante el llamado modelo log-lineal.

Recuérdese que en el post anterior se hablaba de tres variables distintas:

• MOMENTO: Antes o después de la campaña
• SEXO: Hombre o mujer
• RESULTADO: Positivo o negativo en el test de alcoholemia

El modelo log-lineal trata de explicar la frecuencia observada de cada una de las casillas de la tabla tridimensional a partir del efecto de los siguientes factores:

• Un termino constante (c)
• El efecto de la variable SEXO (S)
• El efecto de la variable MOMENTO (F)
• El efecto de la variable RESULTADO (R)
• El efecto de interacción entre SEXO y MOMENTO (SM)
• El efecto de interacción entre SEXO y RESULTADO (SR)
• El efecto de interacción entre MOMENTO y RESULTADO (MR)
• El efecto de interacción entre SEXO, MOMENTO y RESULTADO (SMR)

Es decir, el modelo que se propone es el siguiente:

Freq = c + S + M + R + SM + SR + MR +SMR

El modelo que se acaba de presentar es el llamado modelo saturado. Éste recoge todos los posibles factores y todas sus interacciones (de grado 2, 3, etcétera). Obviamente, el modelo saturado será el que mejor se ajuste a los datos pero, al mismo tiempo será el más complejo.

A la hora de ajustar un modelo log-lineal a los datos observados en el cruce de variables cualitativas debe seguirse un criterio que tome en cuenta no sólo la calidad del ajuste sino también la complejidad del modelo. El procedimiento que suele seguirse es de tipo backwards. Se comienza por el modelo saturado y se van eliminando las interacciones de orden más elevado hasta llegar a un modelo en el que la eliminación de otro efecto -lo que lo haría más simple- perjudique mucho la calidad de ajuste.

En particular, si tomamos los datos del post anterior e intentamos ajustar sobre ellos un modelo log-lineal, encontramos que el efecto de la interacción de tercer orden no es significativo (p-valor= 0,801) y que tampoco es significativa la interacción de orden 2 entre MOMENTO y RESULTADO. De este modo, el modelo log-lineal resultante es el siguiente:

Freq = c + S + M + R + SM + SR

¿Qué es lo que significa este modelo? ¿Explica la aparente contradicción en los resultados a la que se aludía en el post anterior?

Vayamos por partes. Veamos, en primer lugar, que quiere decir la significatividad de los coeficientes de las interacciones: el hecho de que el coeficiente de la interacción de SEXO y MOMENTO resulte significativo implica que las variables SEXO y MOMENTO no son condicionalmente independientes dado el RESULTADO; análogamente, el hecho de que el coeficiente de la interacción entre SEXO y RESULTADO resulte significativo quiere decir que SEXO y RESULTADO no son condicionalmente independientes dado el MOMENTO; por el contrario, el hecho de que el coeficiente de la interacción entre MOMENTO y RESULTADO no sea significativo implica que las variables MOMENTO y RESULTADO son condicionalmente independientes dado el SEXO.

En resumen, las variables MOMENTO y RESULTADO son condicionalmente independientes dado el SEXO. Eso quiere decir que podemos colapsar la tabla de contingencia tridimensional por cualquiera de las variables MOMENTO o RESULTADO y las conclusiones que obtendremos no serán equivocadas:

• Por ejemplo, si colapsamos la tabla tridimensional por la variable MOMENTO obtendremos una tabla bidimensional que cruzará las variables SEXO y RESULTADO. Pues bien, la conclusión respecto a la relación entre SEXO y RESULTADO será la misma si consideramos la tabla global (con los resultados de antes y de después de la campaña juntos) que si consideramos por un lado la tabla SEXO x RESULTADO antes y por otro lado la tabla SEXO x RESULTADO después.
• De modo análogo, si colapsamos la tabla tridimensional por la variable RESULTADO obtendríamos una tabla bidimensional que cruzaría SEXO con MOMENTO. En este caso, aunque con un menor interés, también encontraríamos que los resultados de la tabla colapsada y los resultados de las tablas SEXO x MOMENTO para los que han dado positivo y SEXO x MOMENTO para los que han dado negativo nos llevarían a las mismas conclusiones.

Sin embargo SEXO no es condicionalmente independiente con MOMENTO (dado el RESULTADO). En efecto, consideremos la tabla que cruza SEXO x MOMENTO para los que han dado positivo: antes de la campaña el 23,1% de los interceptados han sido hombres; después de la campaña el 89,3% de los encuestados interceptados han sido hombres. Análogamente, si consideramos la tabla que cruza SEXO x MOMENTO para los que han dado negativo veremos que el 14,9% de los interceptados antes de la campaña han sido hombres y que después de la campaña el 81,5% de los interceptados han sido hombres. En resumen, se ha tendido a parar muchos más hombres después de la campaña y muchos menos hombres antes de la campaña.

Además, SEXO tampoco es condicionalmente independiente de RESULTADO (dado el MOMENTO). En efecto, si consideramos la tabla antes de la campaña observamos que el 30% de los hombres ha dado positivo y sólo el 20% de las mujeres. En la tabla de resultados después de la campaña el 25% de los hombres ha dado positivo mientras que sólo da positivo el 15% de las mujeres. En resumen, las mujeres dan positivo en menor proporción que los hombres.

Si ahora juntamos los resultados de los dos párrafos anteriores la paradoja está servida. La mayor tendencia de los hombres a dar positivo (en comparación con las mujeres) unida a la mayor presencia de los hombres tras la campaña hace que la cifra global (resultado de colapsar la tabla tridimensional por la variable SEXO) nos lleve a conclusiones equivocadas.

En efecto, el hecho de que SEXO no sea condicionalmente independiente de RESULTADO (dado MOMENTO) ni de MOMENTO (dado RESULTADO) hace que no se pueda colapsar la tabla tridimensional por la variable SEXO. Si lo hacemos -y eso es exactamente lo que se propuso en el post anterior- podríamos llegar a conclusiones erróneas. Eso es justamente lo que ocurrió: decidimos colapsar la tabla tridimensional empleando para ello la variable menos adecuada y eso nos llevó a una paradoja: los datos para los hombres y para las mujeres tomados por separado hablaban de la eficacia de la campaña; al colapsarlos, la campaña parecía ineficaz.

¿Sería una buena política intentar ajustar un modelo log-lineal a los datos de las tablas de contingencia antes de aventurarnos a colapsarlas y arriesgarnos a llegar a conclusiones erróneas?

La paradoja de Simpson

Supongamos que la DGT ha lanzado una campaña para concienciar a la población de los peligros que supone conducir ebrio y que desea medir la eficacia de dicha campaña. Con ese fin la policía realiza controles de alcoholemia a 600 conductores antes de la emisión de los anuncios y a otros 600 tras dicha emisión. Los resultados son los siguientes:

• De los 600 conductores a los que se realiza el control de alcoholemia antes de la campaña, 130 (21,66%) presentan una tasa de alcohol superior al límite permitido; de los 600 conductores a los que se realiza el control de alcoholemia tras la campaña, 140 (23,33%) presentan una tasa de alcoholemia superior al límite. Conclusión: la campaña ha sido un fracaso.

Alguien en la DGT decide investigar los datos con mayor profundidad y comprueba que si se considera de manera separada a los hombres y a las mujeres las cifras podrían indicar otra cosa:

• En efecto: antes de la emisión se interceptó a 500 conductoras, de las que 100 (un 20%) dieron positivo. Tras la campaña se realizó el control de alcoholemia a 100 conductoras, de las que 15 (un 15%) dieron positivo. Entre las conductoras, la proporción de positivos ha disminuido por lo que la campaña puede considerarse un éxito.
• En cuanto a los hombres, antes de la campaña se realizó un control de alcoholemia a 100 conductores, de los que 30 (un 30%) dieron positivo. Tras la emisión de los anuncios se realizaron controles de alcoholemia a 500 conductores, de los que 125 (un 25%) dieron positivo. Entre los hombres también ha disminuido la proporción de positivos, por lo que la campaña puede considerarse un éxito entre los hombres.

(Nótese que las 100 mujeres unidas a los 30 hombres que dan positivo antes de la emisión constituyen los 130 conductores que dan positivo a los que se alude en el primer párrafo de los resultados. Análogamente ocurre con las 15 mujeres y los 125 hombres que dan positivo tras la campaña).

¿Cómo es posibe que si realizamos un análisis separando los datos por sexo la campaña sea un éxito pero cuando los agrupamos la campaña sea un fracaso? En realidad, ¿qué podemos decir acerca de la eficacia de la campaña? ¿Qué está pasando?