Curso-Taller MCAyGL

Análisis Estadístico bivariado

Facultad de Matemáticas, Universidad Autónoma de Guerrero

Librerías que se ocupan en esta sesión de trabajo:

library("agricolae")   # Funciones: kruskal, 
library("knitr")       # función kable
library("ggcorrplot")  # función  ggcorrplot
library("dplyr")       # Funciones: rename, mutate, summarise, grop_by, etc
library("vcd")         # Función assocstats
library("vcdExtra")    # Función CMHtest
library("tidyverse")    # %>%
library("nortest")       # pearson.test

1 Introducción

El análisis bivariado es el estudio de dos variables, \(X\) y \(Y\). Dicho estudio se realiza a partir de \(x_{1},...,x_{n}\) una muestra de la variable aleatoria \(X\) y \(y_{1},...,y_{m}\) una muestra de la variable aleatoria \(Y\).

En un análisis bivariado nos debemos:

• Conocer el objetivo de la investigación

• Identificar la escala de las variables

• Identificar la herramienta estadística a utilizar.

Tipos de variables:

• Variables categóricas (cualitativas):

• son las variables nominales (dicotómicas, politómicas) cuyas categorías no tienen orden;

• son las variables ordinales cuyas categorías presentan un orden.

• Variables cuantitativas: variables de escala métrica.

2 Análisis descriptivo bivariado.

2.1 Datos iris

dat_iris<-iris%>%
          rename(Long.Sepal=Sepal.Length, Ancho.Sepal=Sepal.Width, 
                 Long.Petalo=Petal.Length, Ancho.Petalo=Petal.Width,
                 Especies=Species)

->    Long.Sepal Ancho.Sepal Long.Petalo Ancho.Petalo Especies
->  1        5.1         3.5         1.4          0.2   setosa
->  2        4.9         3.0         1.4          0.2   setosa
->  3        4.7         3.2         1.3          0.2   setosa
->  4        4.6         3.1         1.5          0.2   setosa
->  5        5.0         3.6         1.4          0.2   setosa
->  6        5.4         3.9         1.7          0.4   setosa

res<-dat_iris%>%
           select(Especies, Long.Petalo)%>%
           group_by(Especies)%>%
           dplyr::summarise(Total=sum(Long.Petalo))%>%
           mutate(Porcentaje=(Total/sum(Total))*100)

kable(res,digits = 3, caption = "Resultados de la variable Long.Petalo")

Resultados de la variable Long.Petalo

Especies	Total	Porcentaje
setosa	73.1	12.968
versicolor	213.0	37.786
virginica	277.6	49.246

res_mean<-dat_iris%>%
                group_by(Especies)%>%      # Agrupas las datos por categoría
                summarise_all(list(mean))  # Calcula el promedio de las variables por categoría.

kable(res_mean, caption = "Promedio de las variables")

Promedio de las variables

Especies	Long.Sepal	Ancho.Sepal	Long.Petalo	Ancho.Petalo
setosa	5.006	3.428	1.462	0.246
versicolor	5.936	2.770	4.260	1.326
virginica	6.588	2.974	5.552	2.026

res_sd<-dat_iris%>%
                group_by(Especies)%>%
                summarise_all(list(sd))

kable(res_sd, caption = "Desviación estándar de las variables")

Desviación estándar de las variables

Especies	Long.Sepal	Ancho.Sepal	Long.Petalo	Ancho.Petalo
setosa	0.3524897	0.3790644	0.1736640	0.1053856
versicolor	0.5161711	0.3137983	0.4699110	0.1977527
virginica	0.6358796	0.3224966	0.5518947	0.2746501

ggplot(dat_iris, aes(Long.Petalo, fill=Especies)) +
  geom_histogram()+
  facet_wrap(~Especies)+
  theme_bw()+
  labs(x="Longitud del pétado",y="Valores")

ggplot(dat_iris, aes(y=Long.Petalo, fill=Especies)) +
      geom_boxplot()+
      facet_wrap(~Especies)+
      labs(x="Longitud del pétado",y="Valores")+
      theme_bw()+
      theme(legend.position ="none")

Gráfico de cajas y bigotes:

ggplot(dat_iris, aes(x=Especies, y=Long.Petalo, color=Especies))+
geom_boxplot()+                                      # Gráfico de cajas y bigotes
theme_classic()+                                     # Fondo del gráfico blanco
theme(text = element_text(size = 10),                # Tamaño de letra de todas las etiquetas
      axis.text.x = element_text(angle =60,hjust=1), # ángulo de las etiquetas del eje x
      legend.position = "none"                       # Se omite leyenda de color
          )

2.2 Datos corn

corn: Data from a completely randomized design where four different methods of growing corn resulted in various yields per acre on various plots of ground where the four methods were tried.

data(corn)   
str(corn)

->  'data.frame':   34 obs. of  3 variables:
->   $ method     : int  1 1 1 1 1 1 1 1 1 2 ...
->   $ observation: int  83 91 94 89 89 96 91 92 90 91 ...
->   $ rx         : num  11 23 28.5 17 17 31.5 23 26 19.5 23 ...

View(corn)

head(corn)

->    method observation   rx
->  1      1          83 11.0
->  2      1          91 23.0
->  3      1          94 28.5
->  4      1          89 17.0
->  5      1          89 17.0
->  6      1          96 31.5

res_method<-corn%>%
            group_by(method)%>%
            dplyr::summarise(Promedio=mean(rx), Total=sum(rx))%>%
            mutate(Porcentaje=(Total/sum(Total))*100)

kable(res_method,digits = 3, caption = "Resultados del cultivo de maíz.")

Resultados del cultivo de maíz.

method	Promedio	Total	Porcentaje
1	21.833	196.5	33.025
2	15.300	153.0	25.714
3	29.571	207.0	34.790
4	4.812	38.5	6.471

2.3 Datos cotton

data(cotton)
str(cotton)

->  'data.frame':   96 obs. of  5 variables:
->   $ site   : Factor w/ 2 levels "Lima","Pisco": 1 1 1 1 1 1 1 1 1 1 ...
->   $ block  : Factor w/ 6 levels "I","II","III",..: 1 2 3 4 5 6 1 2 3 4 ...
->   $ lineage: int  1 1 1 1 1 1 1 1 1 1 ...
->   $ epoca  : int  1 1 1 1 1 1 2 2 2 2 ...
->   $ yield  : int  12 13 11 13 11 15 15 20 17 18 ...

dat_cotton<-cotton
colnames(dat_cotton)<-c("Sitio","Bloque","Linaje","Epoca","Rendimiento") # Se renombran las columnas

View(dat_cotton)

head(dat_cotton)

->    Sitio Bloque Linaje Epoca Rendimiento
->  1  Lima      I      1     1          12
->  2  Lima     II      1     1          13
->  3  Lima    III      1     1          11
->  4  Lima     IV      1     1          13
->  5  Lima      V      1     1          11
->  6  Lima     VI      1     1          15

res_sitio<-dat_cotton%>%
           group_by(Sitio)%>%
           dplyr::summarise(Promedio=mean(Rendimiento), Total=sum(Rendimiento))%>%
           mutate(Porcentaje=(Total/sum(Total))*100)

kable(res_sitio, digits = 3, caption = "Rendimiento del algodón por sítio")

Rendimiento del algodón por sítio

Sitio	Promedio	Total	Porcentaje
Lima	13.188	633	42.144
Pisco	18.104	869	57.856

res_bloque<-dat_cotton%>%
           group_by(Bloque)%>%
           dplyr::summarise(Promedio=mean(Rendimiento), Total=sum(Rendimiento))%>%
           mutate(Porcentaje=(Total/sum(Total))*100)

kable(res_bloque,digits = 3, caption = "Rendimiento del algodón por bloque.")

Rendimiento del algodón por bloque.

Bloque	Promedio	Total	Porcentaje
I	13.000	208	13.848
II	16.562	265	17.643
III	14.812	237	15.779
IV	16.812	269	17.909
V	16.188	259	17.244
VI	16.500	264	17.577

res_linaje<-dat_cotton%>%
             group_by(Linaje)%>%
             dplyr::summarise(Promedio=mean(Rendimiento), Total=sum(Rendimiento))%>%
             mutate(Porcentaje=(Total/sum(Total))*100)

kable(res_linaje, digits = 3, caption = "Rendimiento del algodón por linaje.")

Rendimiento del algodón por linaje.

Linaje	Promedio	Total	Porcentaje
1	16.542	397	26.431
2	16.958	407	27.097
3	15.333	368	24.501
4	13.750	330	21.971

res_epoca<-dat_cotton%>%
              group_by(Epoca)%>%
              dplyr::summarise(Promedio=mean(Rendimiento), Total=sum(Rendimiento))%>%
              mutate(Porcentaje=(Total/sum(Total))*100)

kable(res_epoca, digits = 3, caption = "Rendimiento del algodón por época.")

Rendimiento del algodón por época.

Epoca	Promedio	Total	Porcentaje
1	13.375	642	42.743
2	17.917	860	57.257

ggplot(dat_cotton, aes(Rendimiento, fill=Sitio)) +
  geom_histogram()+
  facet_wrap(~Sitio)+
  theme_bw()+
  labs(x="Rendimento",y="Valores")

ggplot(dat_cotton, aes(Rendimiento, fill=Bloque)) +
  geom_histogram()+
  facet_wrap(~Bloque)+
  theme_bw()+
  labs(x="Rendimento",y="Valores")

ggplot(dat_cotton, aes(Rendimiento, fill=Linaje)) +
  geom_histogram()+
  facet_wrap(~Linaje)+
  theme_bw()+
  labs(x="Rendimento",y="Valores")

ggplot(dat_cotton, aes(Rendimiento, fill=Epoca)) +
  geom_histogram()+
  facet_wrap(~Epoca)+
  theme_bw()+
  labs(x="Rendimento",y="Valores")

3 Análisis estadístico bivariado

Tipos de errores que se cometen en una prueba de hipótesis:

Tabla	\(H_{0}\) es verdadera	\(H_{0}\) es falsa
Se acepta \(H_{0}\)	Decisión correcta	Error tipo II
Se rechaza \(H_{0}\)	Error tipo I	Decisión correcta

3.1 Comparación de medias en dos poblaciones

• Si el tamaño de muestra menor a \(50\) y la distribución de la variable es normal, entonces se usan las pruebas paramétricas en caso contrario se usa una prueba no paramétrica.

Supongamos que \(X_{1},...,X_{n}\) es una muestra aleatoria de una población con distribución normal \(N(\mu_{X},\sigma_{X}^2)\), y \(Y_{1},...,Y_{m}\) es una muestra aleatoria de una población con distribución normal \(N(\mu_{Y},\sigma_{Y}^2)\); los juegos de hipótesis están dados por:

1. \(H_{0}: \mu_{X}-\mu_{Y}=\delta_{0}\) vs \(H_{0}: \mu_{X}- \mu_{Y}\neq \delta_{0}\)

2. \(H_{0}: \mu_{X}-\mu_{Y} \leq \delta_{0}\) vs \(H_{0}: \mu_{X}-\mu_{Y} > \delta_{0}\)

3. \(H_{0}: \mu_{X}-\mu_{Y}\geq \delta_{0}\) vs \(H_{0}: \mu_{X}-\mu_{Y} < \delta_{0}\)

Se rechaza \(H_{0}\) si p-valor \(<\alpha\)

• Si \(\sigma_{X}=\sigma_{Y}=\sigma\), con \(\sigma\) conocida, el estadístico de prueba está dado por \[ Z=\frac{(\bar{X}-\bar{Y})-(\bar{\mu}_{X}-\bar{\mu}_{Y})}{\sigma\sqrt{1/n+1/m}}\sim N(0,1) \]
• Si Si \(\sigma_{X}=\sigma_{Y}=\sigma\), con \(\sigma\) desconocida , el estadístico de prueba está dado por \[ Z=\frac{(\bar{X}-\bar{Y})-(\bar{\mu}_{X}-\bar{\mu}_{Y})}{S_{p}\sqrt{(1/n+1/m)}}\sim t(n+m-2) \] donde \[ S_{p}^2=\frac{S_{1}^2(n-1)+S_{2}^2(m-1)}{n+m-1} \] En este caso la hipótesis de dos colas se rechaza cuando \[ -t_{\alpha/2,n+m-2}<t<t_{\alpha/2,n+m-2} \]

Ejemplo con la base de datos :

res_lima<-dat_cotton%>%
           select(Sitio,Rendimiento)%>%
           filter(Sitio=="Lima")  
  
res_pisco<-dat_cotton%>%
           select(Sitio,Rendimiento)%>%
           filter(Sitio=="Pisco") 

res<-t.test(res_lima$Rendimiento,res_pisco$Rendimiento)
res

->  
->      Welch Two Sample t-test
->  
->  data:  res_lima$Rendimiento and res_pisco$Rendimiento
->  t = -5.7252, df = 93.038, p-value = 1.256e-07
->  alternative hypothesis: true difference in means is not equal to 0
->  95 percent confidence interval:
->   -6.622016 -3.211317
->  sample estimates:
->  mean of x mean of y 
->   13.18750  18.10417

En este caso se rechaza \(H_{0}\), ya que 1.2557065^{-7} \(<\alpha\), lo cual quiere decir que las medias son diferentes; el intervalo de confianza confirma lo mismo -6.6220165, -3.2113168

Normalidad

\(H_{0}\) Los datos tienen distribución normal

\(H_{1}\) Los datos tienen distribución normal

Rechaza \(H_{0}\) si p-valor \(<\alpha\)

pearson.test(res_lima$Rendimiento) # Prueba para evaluar normalidad

->  
->      Pearson chi-square normality test
->  
->  data:  res_lima$Rendimiento
->  P = 7, p-value = 0.4289

3.1.1 Muestras relacionadas

Un tipo de análisis bivariado se presenta cuando las muestras están relacionadas (pareadas), las cuales se presentan cuando las unidades experimentales se miden antes y después del tratamiento.

Si se tienen \(n\) pares de unidades experimentales para comparar dos tratamientos, y se asignas los dos tratamientos aleatoriamente, uno a cada miembro del par, se obtienen dos muestras relacionadas. La información quedaría dada por \((X_{1},Y_{1}),...,(X_{n},Y_{n})\), donde \(X_{i}\) es la respuesta al primer tratamiento en el par \(i\)-ésimo, y \(Y_{i}\) es la respuesta al segundo tratamiento. Para el análisis estadístico se tiene \(D_{i}=X_{i}-Y_{i}\), \(i=1,...,n\). La independencia entre los pares garantiza que \(D_{1},...,D_{n}\) es una muestra aleatoria de una variable, con media \(\mu_{D}\) y varianza \(\sigma^2_D\). De acuerdo con esto, si \(\mu_{D}=0\), entonces las medias de los tratamiento son iguales , si \(\mu_{D}>0\) indica que la media del primer tratamiento es mayor que la del segundo, y si \(\mu_{D}<0\) entonces la media del segundo tratamiento es mayor. Se asume que la distribución de las muestras relacionadas en normal.

Los juegos de hipótesis están dadas por:

1. \(H_{0}: \mu_{D}=d_{0}\) vs \(H_{0}: \mu_{D}\neq d_{0}\)

2. \(H_{0}: \mu_{D}\leq d_{0}\) vs \(H_{0}: \mu_{D}> d_{0}\)

3. \(H_{0}: \mu_{D}\geq d_{0}\) vs \(H_{0}: \mu_{D}< d_{0}\)

Donde \(d_{0}\) es una constante elegida por el investigador, que generalmente es el cero, \(k=0\).

El estadístico de prueba está dado por

\[ t_{0}=\frac{\sqrt{n}(\bar{D}-k)}{S_{D}} \] donde \(\bar{D}\) y \(S_{D}\) son la media y la varianza de las \(D_{i}\); se asume que \(t_{0}\) tiene una distribución \(t\) de Student con \(n-1\) grados de libertad, cuando \(\mu_{D}=d_{0}\).

Se rechaza \(H_{0}\) si p-valor\(<\alpha\).

Ejemplo: niveles de andrógenos en un grupo se siervos (Walpole, 1992; pag. 328).

Ciervo = 1:15 
Antes = c(2.76,5.18,2.68,3.05,4.1,7.05,6.6,4.79,7.39,7.3,11.78,3.9,26,67.48,17.04)
Despues = c(7.02,3.10,5.44,3.99,5.21,10.26,13.91,18.53,7.91,4.85,11.1,3.74,94.03,94.03,41.70)

dat_pareados = data.frame(Ciervo,Antes, Despues)

attach(dat_pareados)

En ambos casos se cumple el supuesto de normalidad.

Los juegos de hipótesis:

Nivel de significancia \(\alpha=0.05\)

\(H_{0}: \mu_{D}=k\) vs \(H_{0}: \mu_{D}\neq k\)

Rechaza \(H_{0}\) si p-valor\(<\alpha\).

res<-t.test(Despues,Antes, mu = 0,alternative ="two.sided",conf.level = 0.95, paired = TRUE)
res

->  
->      Paired t-test
->  
->  data:  Despues and Antes
->  t = 2.0646, df = 14, p-value = 0.058
->  alternative hypothesis: true mean difference is not equal to 0
->  95 percent confidence interval:
->   -0.3823535 20.0783535
->  sample estimates:
->  mean difference 
->            9.848

Note que 0.0579983\(>\alpha\), por lo tanto no se rechaza \(H_{0}\); note que el intervalo de confianza incluye el cero \(res\)conf.int$, entonces \(\mu_{D}=0\).

Para probar independencia en datos pareados están las pruebas de McNemar y Q de Cochran.

3.2 Comparaciones de medianas

Comparaciones de medianas con la preba de Mann-Whitney:

Sean \(X_{1},....,X_{n}\) y \(Y_{1},....,Y_{n}\) dos muestras aleatorias independientes de poblaciones, cuyos modelos probabilísticos tienen como medianas \(M_{X}\) y \(M_{Y}\) respectivamente.

Juegos de hipótesis:

1. \(H_{0}: M_{X}=M_{Y}\) vs \(H_{1}: M_{X}\neq M_{Y}\)

2. \(H_{0}: M_{X}\leq M_{Y}\) vs \(H_{1}: M_{X} > M_{Y}\)

3. \(H_{0}: M_{X} \geq M_{Y}\) vs \(H_{1}: M_{X} < M_{Y}\)

Se rechaza \(H_{0}\) si p-valor\(<\alpha\)

LA wilcox.test() realiza un test de Mann–Whitney–Wilcoxon entre dos muestras cuando se indica que paired = False y además genera el intervalo de confianza para la diferencia de localización.

muestraX <- c(1.1, 3.4, 4.3, 2.1, 7.0 , 2.5 )
muestraY <- c(7.0, 8.0, 3.0, 5.0, 6.2 , 4.4 )

res<-wilcox.test(x = muestraX, y = muestraY, alternative = "two.sided", mu = 0,
            paired = FALSE, conf.int = 0.95)
res

->  
->      Wilcoxon rank sum test with continuity correction
->  
->  data:  muestraX and muestraY
->  W = 6.5, p-value = 0.07765
->  alternative hypothesis: true location shift is not equal to 0
->  95 percent confidence interval:
->   -4.9000274  0.4000139
->  sample estimates:
->  difference in location 
->               -2.390035

No se rechaza \(H_{0}\), ya que 0.0776483 \(>\alpha\), note que el intervalo tambíen con firma la decisión : -4.9000274, 0.4000139

Prueba de Kruskal-Wallis se utiliza para compara las medianas de \(t\) muestras aleatorias independientes.

3.3 Asociación de variables

En la asociación de variables se busca saber si dadas dos variables, una variable explica a la otra.

3.3.1 Variables cuantitativas

Dadas dos variables, \(X\) y \(Y\), el coeficiente de correlación de Pearson, está dado por

\[ \rho=\frac{Cov(X,Y)}{\sigma_{X}\sigma_{Y}} \] donde \(\sigma_{X}\) es la desviación estándar de \(X\), y \(\sigma_{Y}\) es la desviación estándar de \(Y\).

• \(\rho\) mide la relación lineal que existe entre \(X\) y \(Y\).

• El coeficiente de correlación está entre: \(-1 \leq ρ \leq 1\).

• El coeficiente de correlación no depende de la escala de \(X\) y \(Y\).

• El \(ρ\) es invariante bajo trasformaciones de escala y localidad.

• A mayor asociación lineal entre las v.a \(X\) y \(Y\), más cercano estará \(ρ\) de \(-1\) ó \(1\).

• Este coeficiente de correlación se utiliza para v.a cuantitativas, es decir, las v.a discretas y continuas.

• Si \(X\) y \(Y\) son v.a con distribución normal y \(ρ=0\), entonces \(X\) y \(Y\) son independientes.

Prueba de hipótesis:

Nivel de significancia \(\alpha=0.05\)

Juego de hipótesis:

\(H_{0}: r=0\) vs \(H_{1}: r\neq 0\)

Rechaza \(H_{0}\) si p-valor\(<\alpha\)

->    Long.Sepal Ancho.Sepal Long.Petalo Ancho.Petalo Especies
->  1        5.1         3.5         1.4          0.2   setosa
->  2        4.9         3.0         1.4          0.2   setosa
->  3        4.7         3.2         1.3          0.2   setosa
->  4        4.6         3.1         1.5          0.2   setosa
->  5        5.0         3.6         1.4          0.2   setosa
->  6        5.4         3.9         1.7          0.4   setosa

attach(dat_iris)
res<-cor.test(Long.Petalo,Long.Sepal, methods="pearson")
res

->  
->      Pearson's product-moment correlation
->  
->  data:  Long.Petalo and Long.Sepal
->  t = 21.646, df = 148, p-value < 2.2e-16
->  alternative hypothesis: true correlation is not equal to 0
->  95 percent confidence interval:
->   0.8270363 0.9055080
->  sample estimates:
->        cor 
->  0.8717538

Se rechaza \(H_{0}\), porque 1.0386674^{-47} \(<\alpha\)

correlaciones<-cor(dat_iris[,-5]) # Correlaciones

p.mat<-cor_pmat(dat_iris[,-5]) # p-valores


ggcorrplot(correlaciones,
           hc.order = TRUE, 
           type = "lower", 
           lab = T, 
           insig = "pch", 
           p.mat = p.mat,
           legend.title = "Correlación",
           outline.col = "white",
           sig.level = 0.05,
           lab_size = 4.5,
           title="Correlaciones: Municipio de Victoria")

Correlación no significativa:

res<-cor.test(Ancho.Sepal,Long.Sepal, methods="pearson")
res

->  
->      Pearson's product-moment correlation
->  
->  data:  Ancho.Sepal and Long.Sepal
->  t = -1.4403, df = 148, p-value = 0.1519
->  alternative hypothesis: true correlation is not equal to 0
->  95 percent confidence interval:
->   -0.27269325  0.04351158
->  sample estimates:
->         cor 
->  -0.1175698

No se rechaza \(H_{0}\), ya que 0.1518983\(>\alpha\)

Un intervalo de confianza para el coeficiente de correlación de Pearson, está dado por

3.3.2 Variables categóricas

Dadas dos variables categóricas \(X\) y \(Y\), se pueden tener la siguiente combinación de variables:

• \(Y\) Nominal y \(X\) Nominal

• \(Y\) Nominal y \(X\) Ordinal

• \(Y\) Ordinal y \(X\) Ordinal

En este caso se utilizan estadísticos de prueba no paramétricos para probar la existencia de asociación entre las dos variables categóricas:

Juego de hipótesis:

\(H_{0}\) las variables son independientes vs \(H_{1}\) las variables no son independientes

Se rechaza \(H_{0}\) si \(p-valor<\alpha\).

Estadísticos de prueba:

Variable \(X\)	Variable \(Y\)	Prueba
Nominal	Nominal	Chi-cuadrado, Prueba de Fisher
Ordinal	Nominal	Chi-cuadrado, Prueba de Fisher
Ordinal	Ordinal	si número de categorías en ambas ordinales es \(\leq 5\), utilizar el Coeficiente de Spearman.
ordinal	ordinal	si número de categorías de una o ambas ordinales es \(<5\), utilizar el coeficiente de Kendall tau.
Ordinal	Nominal	Chi-cuadrada, Prueba exacta de Fisher (\(n<100\))

La prueba de Fisher se usa cuando las variables categóricas tienen dos niveles, es decir, cuanto se tiene una tabla de \(2*2\).

La prueba de Chi-cuadrada es buena cuando el número de observaciones esperado en los niveles (categorías) es mayor o igual a \(5\), y el número de observaciones totales en las variables es mayor a \(30\).

Para medir la fuerza de asociación entre dos variables categóricas están los estadísticos de prueba phi, \(\phi\), \(0\leq \phi\leq 1\), Coeficiente V de Cramer, \(0\leq V \leq 1\), y Coeficiente de Contingencia C de Pearson, \(0\leq C \leq \sqrt{2}/2\).

Grado de asociación para tablas de \(2*2\):

• Baja \(0.25-0.50\)

• Moderada \(0.51-0.75\),

• Alta \(0.76-1.0\)

En tablas de \(p*q\) donde \(p>2\) y \(q>2\) es mejor utilizar el Coeficiente V de Cramer.

Ejemplo con Nominal-Nominal:

attach(dat_cotton)  # Se agrega la base de datos dat_cotton

tapply(Rendimiento,list(Sitio), sum)

->   Lima Pisco 
->    633   869

tapply(Rendimiento,list(Epoca), sum)

->    1   2 
->  642 860

tapply(Rendimiento,list(Sitio), mean)

->      Lima    Pisco 
->  13.18750 18.10417

tapply(Rendimiento,list(Epoca), mean)

->         1        2 
->  13.37500 17.91667

tabla <- xtabs(Rendimiento ~Sitio+Epoca, data =dat_cotton)

tabla

->         Epoca
->  Sitio     1   2
->    Lima  268 365
->    Pisco 374 495

fisher.test(tabla, alternative = "two.sided")

->  
->      Fisher's Exact Test for Count Data
->  
->  data:  tabla
->  p-value = 0.792
->  alternative hypothesis: true odds ratio is not equal to 1
->  95 percent confidence interval:
->   0.7856214 1.2017842
->  sample estimates:
->  odds ratio 
->   0.9718066

assocstats(tabla)

->                        X^2 df P(> X^2)
->  Likelihood Ratio 0.073322  1  0.78656
->  Pearson          0.073305  1  0.78658
->  
->  Phi-Coefficient   : 0.007 
->  Contingency Coeff.: 0.007 
->  Cramer's V        : 0.007

class(tabla)

->  [1] "xtabs" "table"

dat<-data.frame(tabla)
dat

->    Sitio Epoca Freq
->  1  Lima     1  268
->  2 Pisco     1  374
->  3  Lima     2  365
->  4 Pisco     2  495

View(dat)

ggplot(dat, aes(x=Sitio, y=Freq, fill=Epoca)) +
       geom_bar(stat="identity") +
       labs(y="Rendimiento total",title="Rendimiento de algodón en dos localidades de Perú")+
       theme()

Ejemplo con Nomial-Nomial:

data(yacon)
attach(yacon)

View(yacon)
levels(locality)

->  [1] "CAJ" "LIM" "OXA"

tabla <- xtabs(glucose~locality+dose, data=yacon)

tabla

->          dose
->  locality     F0   F150    F80
->       CAJ 241.01 333.94 304.23
->       LIM 647.98 492.21 573.08
->       OXA 141.81 128.58 146.96

tapply(glucose,list(locality), sum)

->      CAJ     LIM     OXA 
->   879.18 1713.27  417.35

tapply(glucose,list(dose), sum)

->       F0    F150     F80 
->  1030.80  954.73 1024.27

dat1<-as.data.frame(tabla)
dat1

->    locality dose   Freq
->  1      CAJ   F0 241.01
->  2      LIM   F0 647.98
->  3      OXA   F0 141.81
->  4      CAJ F150 333.94
->  5      LIM F150 492.21
->  6      OXA F150 128.58
->  7      CAJ  F80 304.23
->  8      LIM  F80 573.08
->  9      OXA  F80 146.96

View(dat1)


ggplot(dat1, aes(x=locality, y=Freq, fill=dose)) +
       geom_bar(stat="identity") +
       labs(x="Localidad")+
       theme(axis.text.x=element_text(angle=0,hjust=1,size=9))

chisq.test(tabla) # se rechaza H_0

->  
->      Pearson's Chi-squared test
->  
->  data:  tabla
->  X-squared = 34.627, df = 4, p-value = 5.54e-07

assocstats(tabla) 

->                      X^2 df   P(> X^2)
->  Likelihood Ratio 34.812  4 5.0780e-07
->  Pearson          34.627  4 5.5404e-07
->  
->  Phi-Coefficient   : NA 
->  Contingency Coeff.: 0.107 
->  Cramer's V        : 0.076

CMHtest(tabla)

->  Cochran-Mantel-Haenszel Statistics for locality by dose 
->  
->                   AltHypothesis   Chisq Df       Prob
->  cor        Nonzero correlation  4.1683  1 4.1187e-02
->  rmeans  Row mean scores differ 11.5222  2 3.1476e-03
->  cmeans  Col mean scores differ 17.1903  2 1.8500e-04
->  general    General association 34.6159  4 5.5706e-07

res<-chisq.test(tabla)

res$observed

->          dose
->  locality     F0   F150    F80
->       CAJ 241.01 333.94 304.23
->       LIM 647.98 492.21 573.08
->       OXA 141.81 128.58 146.96

res$expected

->          dose
->  locality       F0     F150      F80
->       CAJ 301.1026 278.8822 299.1952
->       LIM 586.7628 543.4614 583.0457
->       OXA 142.9345 132.3864 142.0291

res$residuals

->          dose
->  locality          F0        F150         F80
->       CAJ -3.46309210  3.29692307  0.29107528
->       LIM  2.52721433 -2.19847622 -0.41272252
->       OXA -0.09406033 -0.33082014  0.41375240