16.1: ANOVA Factorial 1- Diseños Balanceados, Sin Interacciones
- Page ID
- 151322
\( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)
\( \newcommand{\id}{\mathrm{id}}\) \( \newcommand{\Span}{\mathrm{span}}\)
( \newcommand{\kernel}{\mathrm{null}\,}\) \( \newcommand{\range}{\mathrm{range}\,}\)
\( \newcommand{\RealPart}{\mathrm{Re}}\) \( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)
\( \newcommand{\Argument}{\mathrm{Arg}}\) \( \newcommand{\norm}[1]{\| #1 \|}\)
\( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)
\( \newcommand{\Span}{\mathrm{span}}\)
\( \newcommand{\id}{\mathrm{id}}\)
\( \newcommand{\Span}{\mathrm{span}}\)
\( \newcommand{\kernel}{\mathrm{null}\,}\)
\( \newcommand{\range}{\mathrm{range}\,}\)
\( \newcommand{\RealPart}{\mathrm{Re}}\)
\( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)
\( \newcommand{\Argument}{\mathrm{Arg}}\)
\( \newcommand{\norm}[1]{\| #1 \|}\)
\( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)
\( \newcommand{\Span}{\mathrm{span}}\) \( \newcommand{\AA}{\unicode[.8,0]{x212B}}\)
\( \newcommand{\vectorA}[1]{\vec{#1}} % arrow\)
\( \newcommand{\vectorAt}[1]{\vec{\text{#1}}} % arrow\)
\( \newcommand{\vectorB}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vectorC}[1]{\textbf{#1}} \)
\( \newcommand{\vectorD}[1]{\overrightarrow{#1}} \)
\( \newcommand{\vectorDt}[1]{\overrightarrow{\text{#1}}} \)
\( \newcommand{\vectE}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash{\mathbf {#1}}}} \)
\( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)
\(\newcommand{\avec}{\mathbf a}\) \(\newcommand{\bvec}{\mathbf b}\) \(\newcommand{\cvec}{\mathbf c}\) \(\newcommand{\dvec}{\mathbf d}\) \(\newcommand{\dtil}{\widetilde{\mathbf d}}\) \(\newcommand{\evec}{\mathbf e}\) \(\newcommand{\fvec}{\mathbf f}\) \(\newcommand{\nvec}{\mathbf n}\) \(\newcommand{\pvec}{\mathbf p}\) \(\newcommand{\qvec}{\mathbf q}\) \(\newcommand{\svec}{\mathbf s}\) \(\newcommand{\tvec}{\mathbf t}\) \(\newcommand{\uvec}{\mathbf u}\) \(\newcommand{\vvec}{\mathbf v}\) \(\newcommand{\wvec}{\mathbf w}\) \(\newcommand{\xvec}{\mathbf x}\) \(\newcommand{\yvec}{\mathbf y}\) \(\newcommand{\zvec}{\mathbf z}\) \(\newcommand{\rvec}{\mathbf r}\) \(\newcommand{\mvec}{\mathbf m}\) \(\newcommand{\zerovec}{\mathbf 0}\) \(\newcommand{\onevec}{\mathbf 1}\) \(\newcommand{\real}{\mathbb R}\) \(\newcommand{\twovec}[2]{\left[\begin{array}{r}#1 \\ #2 \end{array}\right]}\) \(\newcommand{\ctwovec}[2]{\left[\begin{array}{c}#1 \\ #2 \end{array}\right]}\) \(\newcommand{\threevec}[3]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \end{array}\right]}\) \(\newcommand{\cthreevec}[3]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \end{array}\right]}\) \(\newcommand{\fourvec}[4]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}\) \(\newcommand{\cfourvec}[4]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}\) \(\newcommand{\fivevec}[5]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}\) \(\newcommand{\cfivevec}[5]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}\) \(\newcommand{\mattwo}[4]{\left[\begin{array}{rr}#1 \amp #2 \\ #3 \amp #4 \\ \end{array}\right]}\) \(\newcommand{\laspan}[1]{\text{Span}\{#1\}}\) \(\newcommand{\bcal}{\cal B}\) \(\newcommand{\ccal}{\cal C}\) \(\newcommand{\scal}{\cal S}\) \(\newcommand{\wcal}{\cal W}\) \(\newcommand{\ecal}{\cal E}\) \(\newcommand{\coords}[2]{\left\{#1\right\}_{#2}}\) \(\newcommand{\gray}[1]{\color{gray}{#1}}\) \(\newcommand{\lgray}[1]{\color{lightgray}{#1}}\) \(\newcommand{\rank}{\operatorname{rank}}\) \(\newcommand{\row}{\text{Row}}\) \(\newcommand{\col}{\text{Col}}\) \(\renewcommand{\row}{\text{Row}}\) \(\newcommand{\nul}{\text{Nul}}\) \(\newcommand{\var}{\text{Var}}\) \(\newcommand{\corr}{\text{corr}}\) \(\newcommand{\len}[1]{\left|#1\right|}\) \(\newcommand{\bbar}{\overline{\bvec}}\) \(\newcommand{\bhat}{\widehat{\bvec}}\) \(\newcommand{\bperp}{\bvec^\perp}\) \(\newcommand{\xhat}{\widehat{\xvec}}\) \(\newcommand{\vhat}{\widehat{\vvec}}\) \(\newcommand{\uhat}{\widehat{\uvec}}\) \(\newcommand{\what}{\widehat{\wvec}}\) \(\newcommand{\Sighat}{\widehat{\Sigma}}\) \(\newcommand{\lt}{<}\) \(\newcommand{\gt}{>}\) \(\newcommand{\amp}{&}\) \(\definecolor{fillinmathshade}{gray}{0.9}\)Cuando discutimos el análisis de varianza en el Capítulo 14, asumimos un diseño experimental bastante simple: cada persona cae en uno de varios grupos, y queremos saber si estos grupos tienen medias diferentes en alguna variable de resultado. En esta sección, voy a discutir una clase más amplia de diseños experimentales, conocidos como diseños factoriales, en tenemos más de una variable de agrupación. Dé un ejemplo de cómo este tipo de diseño podría surgir arriba. Otro ejemplo aparece en el Capítulo 14, en el que estábamos viendo el efecto de diferentes drogas en el estado de ánimo.ganancia
experimentada por cada persona. En ese capítulo sí encontramos un efecto significativo del fármaco, pero al final del capítulo también hicimos un análisis para ver si había un efecto de la terapia. No encontramos uno, pero hay algo un poco preocupante en tratar de ejecutar dos análisis separados tratando de predecir el mismo resultado. Tal vez en realidad hay un efecto de la terapia en el aumento del estado de ánimo, pero no lo pudimos encontrar porque estaba siendo “oculto” por el efecto de la droga? En otras palabras, vamos a querer realizar un único análisis que incluya tanto el fármaco
como la terapia
como predictores. Para este análisis cada persona se clasifica de manera cruzada por el medicamento que se le administró (un factor con 3 niveles) y qué terapia recibió (un factor con 2 niveles). Nos referimos a esto como un diseño factorial 3×2. Si tabulamos de forma cruzada fármaco
por terapia
, usando la función xtabs ()
(ver Sección 7.1), obtenemos la siguiente tabla: 230
load("./rbook-master/data/clinicaltrial.rdata")
xtabs( ~ drug + therapy, clin.trial )
## therapy
## drug no.therapy CBT
## placebo 3 3
## anxifree 3 3
## joyzepam 3 3
Como puede ver, no sólo tenemos participantes correspondientes a todas las combinaciones posibles de los dos factores, lo que indica que nuestro diseño está completamente cruzado, resulta que hay un número igual de personas en cada grupo. En otras palabras, tenemos un diseño equilibrado. En esta sección hablaré sobre cómo analizar datos a partir de diseños equilibrados, ya que este es el caso más sencillo. La historia de diseños desequilibrados es bastante tediosa, así que la pondremos a un lado por el momento.
¿Qué hipótesis estamos probando?
Al igual que el ANOVA unidireccional, el ANOVA factorial es una herramienta para probar ciertos tipos de hipótesis sobre medias poblacionales. Entonces, un lugar sensato para comenzar sería ser explícito sobre cuáles son realmente nuestras hipótesis. Sin embargo, antes incluso de que podamos llegar a ese punto, es realmente útil tener alguna notación limpia y simple para describir las medias poblacionales. Debido a que las observaciones se clasifican de manera cruzada en términos de dos factores diferentes, hay bastantes medios diferentes que uno podría interesarle. Para ver esto, comencemos pensando en todos los diferentes medios de muestra que podemos calcular para este tipo de diseño. En primer lugar, está la idea obvia de que podríamos estar interesados en esta tabla de medios grupales:
aggregate( mood.gain ~ drug + therapy, clin.trial, mean )
## drug therapy mood.gain
## 1 placebo no.therapy 0.300000
## 2 anxifree no.therapy 0.400000
## 3 joyzepam no.therapy 1.466667
## 4 placebo CBT 0.600000
## 5 anxifree CBT 1.033333
## 6 joyzepam CBT 1.500000
Ahora bien, este resultado muestra una tabulación cruzada de las medias grupales para todas las combinaciones posibles de los dos factores (por ejemplo, personas que recibieron el placebo y ninguna terapia, personas que recibieron el placebo mientras recibían TCC, etc.). Sin embargo, también podemos construir tablas que ignoran uno de los dos factores. Eso nos da salida que se ve así:
aggregate( mood.gain ~ drug, clin.trial, mean )
## drug mood.gain
## 1 placebo 0.4500000
## 2 anxifree 0.7166667
## 3 joyzepam 1.4833333
aggregate( mood.gain ~ therapy, clin.trial, mean )
## therapy mood.gain
## 1 no.therapy 0.7222222
## 2 CBT 1.0444444
Pero claro, si podemos ignorar un factor ciertamente podemos ignorar ambos. Es decir, también podríamos estar interesados en calcular el aumento promedio del estado de ánimo en los 18 participantes, independientemente de qué medicamento o terapia psicológica se les administró:
mean( clin.trial$mood.gain )
## [1] 0.8833333
¡En este punto tenemos 12 diferentes medios de muestra para hacer un seguimiento! Es útil organizar todos estos números en una sola tabla, que se vería así:
knitr::kable(tibble::tribble(
~V1, ~V2, ~V3, ~V4,
NA, "no therapy", "CBT", "total",
"placebo", "0.30", "0.60", "0.45",
"anxifree", "0.40", "1.03", "0.72",
"joyzepam", "1.47", "1.50", "1.48",
"total", "0.72", "1.04", "0.88"
), col.names = c("", "no therapy", "CBT", "total"))
sin terapia | CBT | total | |
---|---|---|---|
NA | sin terapia | CBT | total |
placebo | 0.30 | 0.60 | 0.45 |
ansioso | 0.40 | 1.03 | 0.72 |
joyzepam | 1.47 | 1.50 | 1.48 |
total | 0.72 | 1.04 | 0.88 |
Ahora, cada una de estas diferentes medias es, por supuesto, una estadística de muestra: es una cantidad que pertenece a las observaciones específicas que hemos hecho durante nuestro estudio. De lo que queremos hacer inferencias son los parámetros poblacionales correspondientes: es decir, los verdaderos medios tal y como existen dentro de alguna población más amplia. Esas medias poblacionales también se pueden organizar en una tabla similar, pero necesitaremos un poco de notación matemática para hacerlo. Como de costumbre, usaré el símbolo μ para denotar una media poblacional. Sin embargo, debido a que hay muchos medios diferentes, voy a necesitar usar subíndices para distinguir entre ellos.
Así es como funciona la notación. Nuestra tabla se define en términos de dos factores: cada fila corresponde a un nivel diferente de Factor A (en este caso fármaco
), y cada columna corresponde a un nivel diferente de Factor B (en este caso terapia
). Si dejamos que R denote el número de filas en la tabla, y C denota el número de columnas, podemos referirnos a esto como un ANOVA factorial R×C. En este caso R=3 y C=2. Usaremos letras minúsculas para hacer referencia a filas y columnas específicas, por lo que μ rc se refiere a la media de la población asociada con el résimo nivel del Factor A (es decir, el número de fila r) y el nivel césimo del Factor B (número de columna c). 231 Entonces los medios poblacionales ahora están escritos así:
knitr::kable(tibble::tribble(
~V1, ~V2, ~V3, ~V4,
"placebo", "$\\mu_{11}$", "$\\mu_{12}$", "",
"anxifree", "$\\mu_{21}$", "$\\mu_{22}$", "",
"joyzepam", "$\\mu_{31}$", "$\\mu_{32}$", "",
"total", "", "", ""
), col.names = c( "", "no therapy", "CBT", "total"))
sin terapia | CBT | total | |
---|---|---|---|
placebo | μ 11 | μ 12 | |
ansioso | μ 21 | μ 22 | |
joyzepam | μ 31 | μ 32 | |
total |
Bien, ¿qué pasa con las entradas restantes? Por ejemplo, ¿cómo debemos describir la ganancia promedio del estado de ánimo en toda la población (hipotética) de personas a las que se les podría dar Joyzepam en un experimento como este, independientemente de si estaban en TCC? Usamos la notación “punto” para expresar esto. En el caso de Joyzepam, fíjese que estamos hablando de la media asociada a la tercera fila de la tabla. Es decir, estamos promediando a través de dos medias celulares (es decir, μ 31 y μ 32). El resultado de este promedio se denomina media marginal, y se denotaría μ 3. en este caso. La media marginal para TCC corresponde a la media poblacional asociada a la segunda columna de la tabla, por lo que usamos la notación μ .2 para describirla. La gran media se denota μ.. porque es la media obtenida promediando (marginando 232) sobre ambos. Entonces nuestra tabla completa de medios poblacionales se puede anotar así:
knitr::kable(tibble::tribble(
~V1, ~V2, ~V3, ~V4,
"placebo", "$\\mu_{11}$", "$\\mu_{12}$", "$\\mu_{1.}$",
"anxifree", "$\\mu_{21}$", "$\\mu_{22}$", "$\\mu_{2.}$",
"joyzepam", "$\\mu_{31}$", "$\\mu_{32}$", "$\\mu_{3.}$",
"total", "$\\mu_{.1}$", "$\\mu_{.2}$", "$\\mu_{..}$"
), col.names=c( NA, "no therapy", "CBT", "total"))
NA | sin terapia | CBT | total |
---|---|---|---|
placebo | μ 11 | μ 12 | μ 1. |
ansioso | μ 21 | μ 22 | μ 2. |
joyzepam | μ 31 | μ 32 | μ 3. |
total | μ .1 | μ .2 | μ.. |
Ahora que tenemos esta notación, es sencillo formular y expresar algunas hipótesis. Supongamos que el objetivo es averiguar dos cosas: en primer lugar, ¿la elección del medicamento tiene algún efecto en el estado de ánimo y, en segundo lugar, la TCC tiene algún efecto en el estado de ánimo? Estas no son las únicas hipótesis que podríamos formular por supuesto, y veremos un ejemplo realmente importante de un tipo diferente de hipótesis en la Sección 16.2, pero estas son las dos hipótesis más simples de probar, y así comenzaremos ahí. Considera la primera prueba. Si la droga no tiene efecto, entonces esperaríamos que todos los medios de la fila fueran idénticos, ¿verdad? Entonces esa es nuestra hipótesis nula. Por otro lado, si el medicamento sí importa entonces debemos esperar que estos medios de fila sean diferentes. Formalmente, escribimos nuestras hipótesis nulas y alternativas en términos de igualdad de medios marginales:
knitr::kable(tibble::tribble(
~V1, ~V2,
"Null hypothesis $H_0$:", "row means are the same i.e. $\\mu_{1.} = \\mu_{2.} = \\mu_{3.}$",
"Alternative hypothesis $H_1$:", "at least one row mean is different."
), col.names = c("", ""))
Hipótesis nula H 0: | las medias de fila son las mismas es decir μ 1. =μ 2. =μ 3. |
Hipótesis alternativa H 1: | al menos una media de fila es diferente. |
Vale la pena señalar que estas son exactamente las mismas hipótesis estadísticas que formamos cuando ejecutamos un ANOVA unidireccional sobre estos datos en el Capítulo 14. En ese entonces usé la notación μ P para referirme a la ganancia media del estado de ánimo para el grupo placebo, con μ A y μ J correspondientes a las medias del grupo para los dos fármacos, y la hipótesis nula fue μ P =μ A =μ J. Entonces en realidad estamos hablando de la misma hipótesis: es solo que el ANOVA más complicado requiere una notación más cuidadosa debido a la presencia de múltiples variables de agrupación, así que ahora nos estamos refiriendo a esta hipótesis como μ 1. =μ 2. =μ 3. . Sin embargo, como veremos en breve, aunque la hipótesis es idéntica, la prueba de esa hipótesis es sutilmente diferente debido a que ahora estamos reconociendo la existencia de la segunda variable de agrupación.
Hablando de la otra variable de agrupación, no te sorprenderá descubrir que nuestra segunda prueba de hipótesis se formula de la misma manera. Sin embargo, como estamos hablando de la terapia psicológica más que de drogas, nuestra hipótesis nula ahora corresponde a la igualdad de los medios de columna:
knitr::kable(tibble::tribble(
~V1, ~V2,
"Null hypothesis $H_0$:", "column means are the same, i.e., $\\mu_{.1} = \\mu_{.2}$",
"Alternative hypothesis $H_1$:", "column means are different, i.e., $\\mu_{.1} \\neq \\mu_{.2}"
), col.names = c("", ""))
Hipótesis nula H 0: | las medias de columna son las mismas, es decir, μ .1 =μ .2 |
Hipótesis alternativa H 1: | las medias de columna son diferentes, es decir, $ {.1} {.2} |
Ejecutando el análisis en R
Las hipótesis nulas y alternativas que describí en la última sección deberían parecer muy familiares: básicamente son las mismas que las hipótesis que estábamos probando en nuestros ANOVA unidireccionales más simples en el Capítulo 14. Así que probablemente esperas que las pruebas de hipótesis que se utilizan en el ANOVA factorial sean esencialmente las mismas que las pruebas F del Capítulo 14. Estás esperando ver referencias a sumas de cuadrados (SS), cuadrados medios (MS), grados de libertad (df), y finalmente una estadística F que podemos convertir en un valor p, ¿verdad? Bueno, tienes toda y completamente razón. Tanto es así que voy a apartarme de mi enfoque habitual. A lo largo de este libro, generalmente he adoptado el enfoque de describir primero la lógica (y hasta cierto punto las matemáticas) que sustenta un análisis particular; y solo luego introducir los comandos R que usarías para producir el análisis. Esta vez voy a hacerlo al revés, y mostrarte primero los comandos R. La razón para hacer esto es que quiero destacar las similitudes entre la sencilla herramienta ANOVA unidireccional que discutimos en el Capítulo 14, y las herramientas más complicadas que vamos a utilizar en este capítulo.
Si los datos que intentas analizar corresponden a un diseño factorial equilibrado, entonces ejecutar tu análisis de varianza es fácil. Para ver lo fácil que es, empecemos por reproducir el análisis original del Capítulo 14. En caso de que lo hayas olvidado, para ese análisis estábamos usando solo un factor (es decir, fármaco
) para predecir nuestra variable de resultado (es decir, mood.gain
), y así esto fue lo que hicimos:
model.1 <- aov( mood.gain ~ drug, clin.trial )
summary( model.1 )
## Df Sum Sq Mean Sq F value Pr(>F)
## drug 2 3.453 1.7267 18.61 8.65e-05 ***
## Residuals 15 1.392 0.0928
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Tenga en cuenta que en esta ocasión he usado el nombre model.1
como etiqueta para mi objeto aov
, ya que estoy planeando crear bastantes otros modelos también. Para empezar, supongamos que también tengo curiosidad por saber si la terapia
tiene una relación con el estado de ánimo.ganancia
. A la luz de lo que hemos visto en nuestra discusión sobre la regresión múltiple en el Capítulo 15, probablemente no te sorprenderá que todo lo que tenemos que hacer es extender la fórmula: en otras palabras, si especificamos mood.gain ~ drug + therapy
como nuestro modelo, probablemente obtendremos lo que buscamos:
model.2 <- aov( mood.gain ~ drug + therapy, clin.trial )
Esta salida también es bastante simple de leer: la primera fila de la tabla reporta un valor de suma de cuadrados (SS) entre grupos asociado con el factor fármaco
, junto con un valor df correspondiente entre grupos. También calcula un valor cuadrático medio (MS), un estadístico F y un valor p. También hay una fila correspondiente al factor de terapia
, y una fila correspondiente a los residuos (es decir, la variación dentro de los grupos).
No solo todas las cantidades individuales son bastante familiares, las relaciones entre estas diferentes cantidades se han mantenido sin cambios: al igual que vimos con el ANOVA unidireccional original, tenga en cuenta que el valor cuadrático medio se calcula dividiendo SS por el df correspondiente. Es decir, sigue siendo cierto que
\(\ MS = {SS \over df}\)
independientemente de si estamos hablando de drogas
, terapia
o los residuos. Para ver esto, no nos preocupemos por cómo se calculan las sumas de los valores cuadrados: en cambio, tomemos con fe que R ha calculado correctamente los valores SS, e intentemos verificar que todo el resto de los números tengan sentido. Primero, tenga en cuenta que para el factor fármaco
, dividimos 3.45 por 2, y terminamos con un valor cuadrático medio de 1.73. Para el factor terapia
, solo hay 1 grado de libertad, por lo que nuestros cálculos son aún más simples: dividir 0.47 (el valor SS) por 1 nos da una respuesta de 0.47 (el valor MS).
Pasando a la estadística F y a los valores p, observe que tenemos dos de cada uno: uno correspondiente al factor fármaco
y el otro correspondiente al factor de terapia
. Independientemente de cuál estemos hablando, el estadístico F se calcula dividiendo el valor cuadrático medio asociado al factor por el valor cuadrático medio asociado a los residuales. Si usamos “A” como notación taquigráfica para referirnos al primer factor (factor A; en este caso fármaco
) y “R” como notación taquigráfica para referirnos a los residuales, entonces el estadístico F asociado al factor A se denota F A, y se calcula de la siguiente manera:
\(\ F_A={MS_A \over MS_R}\)
y existe una fórmula equivalente para el factor B (es decir, terapia
). Tenga en cuenta que este uso de “R” para referirse a los residuos es un poco incómodo, ya que también usamos la letra R para hacer referencia al número de filas en la tabla, pero solo voy a usar “R” para significar residuos en el contexto de SSR y MSR, así que ojalá esto no sea confuso. De todos modos, para aplicar esta fórmula al factor fármacos
, tomamos el cuadrado medio de 1.73 y lo dividimos por el valor cuadrático medio residual de 0.07, lo que nos da un estadístico F de 26.15. El cálculo correspondiente para la variable terapéutica
sería dividir 0.47 por 0.07 lo que da 7.08 como el estadístico F. No es sorprendente, por supuesto, estos son los mismos valores que R ha reportado en la tabla ANOVA anterior.
La última parte de la tabla ANOVA es el cálculo de los valores p. Una vez más, aquí no hay nada nuevo: para cada uno de nuestros dos factores, lo que estamos tratando de hacer es probar la hipótesis nula de que no hay relación entre el factor y la variable de resultado (voy a ser un poco más preciso sobre esto más adelante). Para ello, hemos seguido (aparentemente) una estrategia similar a la que hicimos de una manera ANOVA, y hemos calculado un estadístico F para cada una de estas hipótesis. Para convertir estos en valores p, todo lo que necesitamos hacer es señalar que el que la distribución de muestreo para el estadístico F bajo la hipótesis nula (que el factor en cuestión es irrelevante) es una distribución F: y que los valores de dos grados de libertad son aquellos correspondientes al factor, y los correspondientes a los residuales. Para el factor droga
estamos hablando de una distribución F con 2 y 14 grados de libertad (voy a discutir grados de libertad con más detalle más adelante). En contraste, para el factor de terapia
la distribución muestral es F con 1 y 14 grados de libertad. Si realmente quisiéramos, podríamos calcular el valor p nosotros mismos usando la función pf ()
(ver Sección 9.6). Solo para demostrar que no está pasando nada gracioso, así es como se vería para la variable droga
:
pf( q=26.15, df1=2, df2=14, lower.tail=FALSE )
## [1] 1.871981e-05
Y como puede ver, este es efectivamente el valor p reportado en la tabla ANOVA anterior.
En este punto, espero que puedan ver que la tabla ANOVA para este análisis más complicado correspondiente al modelo.2
debe leerse de la misma manera que la tabla ANOVA para el análisis más simple para modelo.1
. En definitiva, nos está diciendo que el ANOVA factorial para nuestro diseño 3×2 encontró un efecto significativo del fármaco (F 2,14 =26.15, p<.001) así como un efecto significativo de la terapia (F 1,14=7.08, p=.02). O bien, para usar la terminología más técnicamente correcta, diríamos que hay dos efectos principales de la droga y la terapia. Por el momento, probablemente parezca un poco redundante referirse a estos como efectos “principales”: pero en realidad tiene sentido. Más adelante, vamos a querer hablar sobre la posibilidad de “interacciones” entre los dos factores, por lo que generalmente hacemos una distinción entre los efectos principales y los efectos de interacción.
la suma de los cuadrados calculados?
En la sección anterior tenía dos objetivos: primero, mostrarles que los comandos R necesarios para hacer ANOVA factorial son más o menos los mismos que usamos para un ANOVA unidireccional. La única diferencia es el argumento de fórmula
para la función aov ()
. En segundo lugar, quería mostrarles cómo se ve la tabla ANOVA en este caso, para que se pueda ver desde el principio que la lógica y estructura básica detrás del ANOVA factorial es la misma que la que sustenta el ANOVA unidireccional. Intenta aferrarte a ese sentimiento. Es realmente cierto, en la medida en que el ANOVA factorial se construye más o menos de la misma manera que el modelo ANOVA unidireccional más simple. Es solo que este sentimiento de familiaridad empieza a evaporarse una vez que empiezas a indagar en los detalles. Tradicionalmente, esta sensación reconfortante es reemplazada por un impulso de asesinar a los autores de los libros de texto de estadística.
Bien, comencemos a mirar algunos de esos detalles. La explicación que di en la última sección ilustra el hecho de que las pruebas de hipótesis para los efectos principales (de fármaco y terapia en este caso) son pruebas F, pero lo que no hace es mostrarte cómo se calculan los valores de la suma de cuadrados (SS). Tampoco te dice explícitamente cómo calcular grados de libertad (valores df) aunque eso es algo sencillo por comparación. Supongamos por ahora que solo tenemos dos variables predictoras, Factor A y Factor B. Si usamos Y para referirnos a la variable de resultado, entonces usaríamos Yrci para referirnos al resultado asociado con el i-ésimo miembro del grupo rc (es decir, nivel/fila r para Factor A y nivel/columna c para Factor B). Así, si usamos\(\ \bar{Y}\) para referirnos a una media muestral, podemos usar la misma notación que antes para referirnos a medias grupales, medias marginales y grandes medias:\(\ \bar{Y_{rc}}\) es decir, es la media muestral asociada al résimo nivel del Factor A y al nivel césimo del Factor B,\(\ \bar{Y_r}\). sería la media marginal para el résimo nivel del Factor A,\(\ \bar{Y_{.c}}\) sería la media marginal para el nivel césimo del Factor B, y\(\ \bar{Y_{..}}\) es la gran media. Es decir, nuestras medias muestrales pueden organizarse en la misma tabla que las medias poblacionales. Para nuestros datos de ensayos clínicos, esa tabla se ve así:
knitr::kable(tibble::tribble(
~V1, ~V2, ~V3, ~V4,
"placebo", "$\\bar{Y}_{11}$", "$\\bar{Y}_{12}$", "$\\bar{Y}_{1.}$",
"anxifree", "$\\bar{Y}_{21}$", "$\\bar{Y}_{22}$", "$\\bar{Y}_{2.}$",
"joyzepam", "$\\bar{Y}_{31}$", "$\\bar{Y}_{32}$", "$\\bar{Y}_{3.}$",
"total", "$\\bar{Y}_{.1}$", "$\\bar{Y}_{.2}$", "$\\bar{Y}_{..}$"
), col.names = c("", "no therapy", "CBT", "total"))
sin terapia | CBT | total | |
---|---|---|---|
placebo | \(\ \bar{Y_{11}}\) | \(\ \bar{Y_{12}}\) | \(\ \bar{Y_{1.}}\) |
ansioso | \(\ \bar{Y_{21}}\) | \(\ \bar{Y_{22}}\) | \(\ \bar{Y_{2.}}\) |
joyzepam | \(\ \bar{Y_{31}}\) | \(\ \bar{Y_{32}}\) | \(\ \bar{Y_{3.}}\) |
total | \(\ \bar{Y_{.1}}\) | \(\ \bar{Y_{.2}}\) | \(\ \bar{Y_{..}}\) |
Y si miramos los medios de muestra que mostré antes, tenemos\(\ \bar{Y_{11}}\) =0.30,\(\ \bar{Y_{12}}\) =0.60 etc. En nuestro ejemplo de ensayo clínico, el factor fármacos
tiene 3 niveles y el factor de terapia
tiene 2 niveles, y entonces lo que estamos tratando de ejecutar es un ANOVA factorial 3×2. Sin embargo, seremos un poco más generales y diremos que el Factor A (el factor de fila) tiene niveles R y el Factor B (el factor de columna) tiene niveles C, y entonces lo que estamos corriendo aquí es un ANOVA factorial R×C.
Ahora que tenemos nuestra notación recta, podemos calcular la suma de los valores cuadrados para cada uno de los dos factores de una manera relativamente familiar. Para el Factor A, nuestra suma de cuadrados entre grupos se calcula evaluando la medida en que las medias marginales (filas)\(\ \bar{Y_{1.}}\),\(\ \bar{Y_{2.}}\) etc, son diferentes de la gran media\(\ \bar{Y_{..}}\). Hacemos esto de la misma manera que lo hicimos para ANOVA unidireccional: calcular la suma de la diferencia cuadrada entre los\(\ \bar{Y_{i.}}\) valores y los\(\ \bar{Y_{..}}\) valores. Específicamente, si hay N personas en cada grupo, entonces calculamos esto:
\(\ SS_A=(N \times C) \sum_{r=1}^R(\bar{Y_{r.}} - \bar{Y_{..}})^2\)
Al igual que con el ANOVA unidireccional, la parte 233 más interesante de esta fórmula es el\(\ (\bar{Y_{r.}} - \bar{Y_{..}})^2\) bit, que corresponde a la desviación cuadrada asociada al nivel r. Todo lo que hace esta fórmula es calcular esta desviación cuadrada para todos los niveles R del factor, sumarlos y luego multiplicar el resultado por N×C La razón de esta última parte es que hay múltiples celdas en nuestro diseño que tienen el nivel r en el Factor A: de hecho, hay C de ellas, ¡una correspondiente a cada nivel posible del Factor B! Por ejemplo, en nuestro ejemplo de juguete, hay dos células diferentes en el diseño correspondientes al medicamento ansifree
: una para personas con no.terapia
y otra para el grupo CBT
. No sólo eso, dentro de cada una de estas celdas hay N observaciones. Entonces, si queremos convertir nuestro valor SS en una cantidad que calcule la suma de cuadrados entre grupos sobre una base “por observación”, tenemos que multiplicar por N×C. La fórmula para el factor B es, por supuesto, lo mismo, solo con algunos subíndices barajados:
\(\ SS_B = (N \times R)\sum_{c=1}^C(\bar{Y_{.c}} - \bar{Y_{..}})^2\)
Ahora que tenemos estas fórmulas, podemos verificarlas contra la salida R de la sección anterior. Primero, observe que calculamos todas las medias marginales (es decir, medias marginales de fila\(\ \bar{Y_{r.}}\) y medias marginales de columna\(\ \bar{Y_{.c}}\)) antes usando agregado ()
, y también calculamos la media grande. Repetimos esos cálculos, pero esta vez guardaremos los resultados en variables para que podamos usarlos en cálculos posteriores:
drug.means <- aggregate( mood.gain ~ drug, clin.trial, mean )[,2]
therapy.means <- aggregate( mood.gain ~ therapy, clin.trial, mean )[,2]
grand.mean <- mean( clin.trial$mood.gain )
Bien, ahora calculemos la suma de cuadrados asociados al efecto principal de la droga
. Hay un total de N=3 personas en cada grupo, y C=2 tipos diferentes de terapia. O, para decirlo de otra manera, hay 3×2=6 personas que recibieron alguna droga en particular. Entonces nuestros cálculos son:
SS.drug <- (3*2) * sum( (drug.means - grand.mean)^2 )
SS.drug
## [1] 3.453333
No en vano, este es el mismo número que obtienes cuando buscas el valor SS para el factor drogas en la tabla ANOVA que presenté anteriormente. Podemos repetir el mismo tipo de cálculo para el efecto de la terapia. Nuevamente hay N=3 personas en cada grupo, pero como hay R=3 medicamentos diferentes, esta vez observamos que hay 3×3=9 personas que recibieron TCC, y otras 9 personas que recibieron el placebo. Entonces nuestro cálculo es ahora:
SS.therapy <- (3*3) * sum( (therapy.means - grand.mean)^2 )
SS.therapy
## [1] 0.4672222
y, una vez más, no nos sorprende ver que nuestros cálculos son idénticos a la salida del ANOVA.
Entonces así es como se calculan los valores SS para los dos efectos principales. Estos valores SS son análogos a los valores de suma de cuadrados entre grupos que calculamos al hacer ANOVA unidireccional en el Capítulo 14. Sin embargo, ya no es buena idea pensar en ellos como valores SS entre grupos, solo porque tenemos dos variables de agrupación diferentes y es fácil confundirse. Para construir una prueba F, sin embargo, también necesitamos calcular la suma de cuadrados dentro de grupos. De acuerdo con la terminología que usamos en el capítulo de regresión (Capítulo 15) y la terminología que R utiliza al imprimir la tabla ANOVA, comenzaré a referirme al valor SS dentro de grupos como la suma residual de cuadrados SS R.
La forma más fácil de pensar sobre los valores residuales de SS en este contexto, creo, es pensarlo como la variación sobrante en la variable de resultado después de tomar en cuenta las diferencias en las medias marginales (es decir, después de eliminar SSA y SSB). Lo que quiero decir con eso es que podemos comenzar calculando la suma total de cuadrados, que voy a etiquetar SST. La fórmula para esto es prácticamente la misma que para el ANOVA unidireccional: tomamos la diferencia entre cada observación Y rci y la gran media\(\ \bar{Y_{..}}\), cuadramos las diferencias y las sumamos todas
\(\ SS_T=\sum_{r=1}^R\sum_{c=1}^C\sum_{i=1}^N(Y_{rci} - \bar{Y_{..}})^2\)
La “triple suma” aquí parece más complicada de lo que es. En las dos primeras sumas, estamos sumando en todos los niveles del Factor A (es decir, sobre todas las filas posibles r en nuestra tabla), a través de todos los niveles del Factor B (es decir, todas las columnas posibles c). Cada combinación rc corresponde a un solo grupo, y cada grupo contiene N personas: así que también tenemos que sumar a todas esas personas (es decir, todos los valores i). En otras palabras, todo lo que estamos haciendo aquí es sumar todas las observaciones en el conjunto de datos (es decir, todas las combinaciones posibles de rci).
En este punto, conocemos la variabilidad total de la variable de resultado SS T, y sabemos cuánto de esa variabilidad se puede atribuir al Factor A (SS A) y cuánto de ella se puede atribuir al Factor B (SS B). La suma residual de cuadrados se define así como la variabilidad en Y que no se puede atribuir a ninguno de nuestros dos factores. En otras palabras:
SS R = SS T − (SS A +SS B)
Por supuesto, hay una fórmula que puedes usar para calcular directamente el SS residual, pero creo que tiene más sentido conceptual pensarlo así. El objetivo de llamarlo residual es que es la variación sobrante, y la fórmula anterior lo deja claro. También debo señalar que, de acuerdo con la terminología utilizada en el capítulo de regresión, es común referirse a SS A +SS B como la varianza atribuible al “modelo ANOVA”, denotado SS M, por lo que a menudo decimos que la suma total de cuadrados es igual a la suma del modelo de cuadrados más la suma residual de cuadrados. Más adelante en este capítulo veremos que esto no es solo una similitud superficial: el ANOVA y la regresión son en realidad lo mismo bajo el capó.
En cualquier caso, probablemente valga la pena tomarse un momento para verificar que podamos calcular SS R usando esta fórmula, y verificar que sí obtenemos la misma respuesta que R produce en su tabla ANOVA. Los cálculos son bastante sencillos. Primero calculamos la suma total de cuadrados:
SS.tot <- sum( (clin.trial$mood.gain - grand.mean)^2 )
SS.tot
## [1] 4.845
y luego lo usamos para calcular la suma residual de cuadrados:
SS.res <- SS.tot - (SS.drug + SS.therapy)
SS.res
## [1] 0.9244444
Una vez más, obtenemos la misma respuesta.
¿Cuáles son nuestros grados de libertad?
Los grados de libertad se calculan de la misma manera que para el ANOVA unidireccional. Para cualquier factor dado, los grados de libertad son iguales al número de niveles menos 1 (es decir, R−1 para la variable de fila, Factor A y C−1 para la variable de columna, Factor B). Entonces, para el factor fármacos
obtenemos df=2, y para el factor de terapia
obtenemos df=1. Más adelante, cuando discutamos la interpretación del ANOVA como modelo de regresión (ver Sección 16.6) voy a dar una declaración más clara de cómo llegamos a este número, pero por el momento podemos usar la definición simple de grados de libertad, es decir, que los grados de libertad equivalen al número de cantidades que son observado, menos el número de restricciones. Entonces, para el factor drogas
, observamos 3 medias grupales separadas, pero éstas están restringidas por 1 gran media; y por lo tanto los grados de libertad son 2. Para los residuales, la lógica es similar, pero no exactamente la misma. El número total de observaciones en nuestro experimento es de 18. Las restricciones corresponden a la media 1 gran, los 2 grupos adicionales significan que el factor fármaco
introduce, y el 1 grupo adicional significa que el factor terapia
introduce, por lo que nuestros grados de libertad son 14. Como fórmula, esto es N−1− (R−1) − (C−1), lo que simplifica a N−R−C+1.
ANOVA factorial versus ANOVA de una vía
Ahora que hemos visto cómo funciona un ANOVA factorial, vale la pena tomarse un momento para compararlo con los resultados de los análisis unidireccionales, porque esto nos dará una muy buena idea de por qué es una buena idea ejecutar el ANOVA factorial. En el Capítulo 14 que ejecuté un ANOVA unidireccional que buscó ver si hay alguna diferencia entre fármacos, y un segundo ANOVA unidireccional para ver si había alguna diferencia entre las terapias. Como vimos en la Sección 16.1.1, las hipótesis nulas y alternativas probadas por los ANOVA unidireccionales son de hecho idénticas a las hipótesis probadas por el ANOVA factorial. Mirando aún más cuidadosamente las tablas de ANOVA, podemos ver que la suma de cuadrados asociados a los factores es idéntica en los dos análisis diferentes (3.45 para fármaco
y 0.92 para terapia
), al igual que los grados de libertad (2 para fármaco
, 1 para terapia
). ¡Pero no dan las mismas respuestas! Lo más notable es que cuando ejecutamos el ANOVA de una vía para terapia
en la Sección 14.11 no encontramos un efecto significativo (el valor de p fue de 0.21). Sin embargo, cuando observamos el efecto principal de la terapia
dentro del contexto del ANOVA bidireccional, sí obtenemos un efecto significativo (p=.019). Los dos análisis claramente no son los mismos.
¿Por qué sucede eso? La respuesta radica en entender cómo se calculan los residuos. Recordemos que toda la idea detrás de una prueba F es comparar la variabilidad que puede atribuirse a un factor particular con la variabilidad que no se puede explicar (los residuales). Si ejecuta un ANOVA unidireccional para terapia
y, por lo tanto, ignora el efecto del medicamento
, ¡el ANOVA terminará arrojando toda la variabilidad inducida por el fármaco en los residuos! Esto tiene el efecto de hacer que los datos parezcan más ruidosos de lo que realmente son, y el efecto de la terapia
que se encuentra correctamente como significativo en el ANOVA bidireccional ahora se vuelve no significativo. Si ignoramos algo que realmente importa (por ejemplo, la droga
) al tratar de evaluar la contribución de otra cosa (por ejemplo, la terapia
) entonces nuestro análisis se distorsionará. Por supuesto, está perfectamente bien ignorar variables que son genuinamente irrelevantes para el fenómeno de interés: si hubiéramos grabado el color de las paredes, y eso resultó ser no significativo en un ANOVA de tres vías (es decir, estado de ánimo.ganancia ~ fármaco + terapia + pared.color
), estaría perfectamente bien ignorarlo y solo reportar el ANOVA bidireccional más simple que no incluya este factor irrelevante. ¡Lo que no debes hacer es soltar variables que realmente marcan la diferencia!
¿Qué tipo de resultados capta este análisis?
El modelo ANOVA del que hemos estado hablando hasta ahora cubre una gama de diferentes patrones que podríamos observar en nuestros datos. Por ejemplo, en un diseño de ANOVA bidireccional, existen cuatro posibilidades: (a) solo importa el Factor A, (b) solo el Factor B importa, (c) tanto A como B importan, y (d) ni A ni B importan. Un ejemplo de cada una de estas cuatro posibilidades se traza en la Figura?? .