Objetivos

Continuar utilizando las funciones del conjunto de librerías tidyverse para la manipulación y análisis de datos, y la librería flextable para la presentación de resultados.
Utilizar BinomDiffCI y BinomRatioCI del paquete DescTools para calcular, respectivamente, intervalos de confianza para la diferencia y el cociente de proporciones utilizando muestras independientes. En el caso de muestras emparejadas, utilizar diffpropci.mp del paquete propCIs para obtener un intervalo de confianza para la diferencia de proporciones.
Utilizar t.test para calcular intervalos de confianza para la diferencia de medias de muestras independientes o emparejadas.
Utilizar var.test para calcular intervalos de confianza para el cociente de varianzas o desviaciones típicas de muestras independientes, y Var.test del paquete PairedData cuando la muestras son emparejadas.
Seguir utilizando las funciones do() (de tidyverse) y tidy (de la librería broom) para integrar el cálculo de intervalos de confianza en cadenas de comandos de tidyverse enlazados por tuberías.
Seguir practicando con la construcción de tablas resumen con estadísticos e intervalos de confianza para varias variables.

Preliminares: carga de las librerías que utilizaremos en esta práctica.

library(readxl)
library(flextable)
library(janitor)
library(DescTools)
library(broom)
library(PropCIs)
library(PairedData)
library(tidyverse)

En esta práctica utilizaremos de nuevo la función formatInterval() que creamos en la práctica anterior para dar formato a los intervalos de confianza. Como ya señalamos en dicha práctica si se va a utilizar esta función, es preciso incluirla en el documento Rmd antes de llamarla por primera vez. Esta función se ha modificado ligeramente para poder etiquetar con un nombre al estimador (por defecto el nombre es estimate).

# La función recibe como argumentos, el estimador del parámetro y los límites 
# inferior y superior del intervalo de confianza. Por defecto se redondean los
# valores a 4 dígitos y se considera que el nivel de confianza el 0.95. El valor
# del nivel de confianza se usa solamente para poner nombre a la variable que 
# designa el intervalo, como suele ser habitual en publicaciones científicas
# La función devuelve un data.frame con el estimador y el intervalo formateados
# con el número de dígitos indicado. La opción `label` permite especificar una
# etiqueta para el nombre del estimador
formatInterval <- function(est, inf, sup, digits=4, 
                           conf.level=0.95, label="Estimate"){
  CI <- sprintf("(%s, %s)",round(inf,digits),round(sup,digits))
  estimate <- round(est,digits)
  df <- data.frame(estimate,CI, stringsAsFactors = FALSE)
  names(df)=c(label,paste("CI",100*conf.level,"%",sep=""))
  return(df)
}

1. Muestras independientes

Utilizaremos los datos del archivo tortugas.csv.

tortugas <- read_csv2("tortugas.csv")

1.1. Intervalo de confianza para la diferencia de proporciones

Ejercicio 1

Construir e interpretar un intervalo de confianza al 95% para la diferencia entre las proporciones de nidos afectados por cangrejos en las playas de Ervatao y Ponta Cosme.

Para construir un intervalo de confianza para la diferencia de proporciones cuando la muestra es grande podemos utilizar la función prop.test() del paquete base de R, o la función BinomDiffCI() de la librería DescTools. Para aplicar cualquiera de estas funciones necesitamos saber cuántos nidos se revisaron en cada playa y en cuántos había cangrejos. Para ello, aprovechando que la variable cangrejos vale 1 cuando el nido tiene cangrejos y 0 en caso contrario, podemos usar la sintaxis:

tortugas %>%
  filter(!is.na(cangrejos)) %>% 
  group_by(playa) %>% 
  summarize(nE=sum(cangrejos),n=n())

## # A tibble: 4 x 3
##   playa             nE     n
##   <chr>          <dbl> <int>
## 1 Calheta           97   142
## 2 Ervatao          159   290
## 3 Ponta Cosme      258   356
## 4 Porto Ferreiro    81   112

Podíamos haber utilizado también tabyl():

tortugas %>% 
  filter(!is.na(cangrejos)) %>% 
  tabyl(playa,cangrejos) %>% 
  adorn_totals("col")

##           playa   0   1 Total
##         Calheta  45  97   142
##         Ervatao 131 159   290
##     Ponta Cosme  98 258   356
##  Porto Ferreiro  31  81   112

Así pues, en Ervatao se revisaron 290 nidos de los cuales tenían cangrejos 159; y en Ponta Cosme se revisaron 356, de los que 258 tenían cangrejos. Las proporciones estimadas de nidos con cangrejos en cada playa son entonces:

\[\hat{p}_{Ervatao}=\frac{159}{290}=0.5483\] \[\hat{p}_{Ponta Cosme}=\frac{258}{356}=0.7247\]

Es decir, en Ponta Cosme un 72.47% de los nidos tenían cangrejos frente a un 54.83% en Ervatao. La diferencia estimada es entonces de un 72.47-54.83 = 17.64%

Para obtener el intervalo de confianza para esta diferencia utilizando prop.test():

prop.test(x=c(159,258),n=c(290,356),conf.level=0.95)

## 
##  2-sample test for equality of proportions with continuity
##  correction
## 
## data:  c(159, 258) out of c(290, 356)
## X-squared = 20.979, df = 1, p-value = 4.644e-06
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.25328397 -0.09960251
## sample estimates:
##    prop 1    prop 2 
## 0.5482759 0.7247191

Para obtener el intervalo usando BinomDiffCI():

BinomDiffCI(159,290,258,356,conf.level=0.95)

##             est     lwr.ci     upr.ci
## [1,] -0.1764432 -0.2490401 -0.1019968

Como vemos, los dos intervalos son parecidos (se diferencian en milésimas), si bien el intervalo calculado por BinomDiffCI(), debido a Agresti-Caffo, es más preciso. Por tanto, podemos afirmar que, si bien en las muestras la diferencia de proporciones ha sido del 17.6%, con un 95% de confianza podemos asegurar que en las poblaciones la diferencia está entre un 10.2% y un 24.9%; o dicho de otra forma, en la playa de Ponta Cosme hay entre un 10.2% y un 24.9% más de nidos afectados por Cangrejos que en la playa de Ervatao.

Mejora de la presentación del resultado

Para mejorar la presentación, podemos proceder del mismo modo que ya hicimos en la práctica 8 para construir un intervalo de confianza para una proporción, con la única diferencia de que ahora usamos BinomDiffCI() y debemos agrupar por playa para obtener el número de cangrejos afectados en cada playa:

tortugas %>% 
  filter(!is.na(cangrejos)&(playa=="Ponta Cosme"|playa=="Ervatao")) %>%
  group_by(playa) %>% 
  summarise(nE=sum(cangrejos),n=n()) %>% 
  do(data.frame(BinomDiffCI(.$nE[1],.$n[1],.$nE[2],.$n[2],conf.level=0.95))) %>% 
  do(formatInterval(.$est,.$lwr.ci,.$upr.ci,label="difference")) %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

difference	CI95%
-0.1764	(-0.249, -0.102)

Ejercicio 2

Se desea explorar si la diferencia entre las proporciones de nidos con cangrejos en Ervatao y Ponta Cosme se ha mantenido estable a lo largo del tiempo. Construye un intervalo de confianza para la diferencia entre ambas proporciones cada uno de los años en que se realizaron campañas.

Podemos explorar cada año cuantos nidos tenían cangrejos y en qué proporción:

porAño <- tortugas %>% 
  filter(!is.na(cangrejos)&(playa=="Ponta Cosme"|playa=="Ervatao")) %>%
  group_by(Año,playa) %>% 
  summarise(nE=sum(cangrejos),n=n(), pE=nE/n) 
porAño %>% flextable() %>% autofit()

Año	playa	nE	n	pE
1999	Ervatao	25	60	0.4166667
1999	Ponta Cosme	34	47	0.7234043
2000	Ervatao	39	59	0.6610169
2000	Ponta Cosme	39	58	0.6724138
2001	Ervatao	27	41	0.6585366
2001	Ponta Cosme	47	59	0.7966102
2002	Ervatao	28	49	0.5714286
2002	Ponta Cosme	39	54	0.7222222
2003	Ervatao	17	37	0.4594595
2003	Ponta Cosme	34	51	0.6666667
2004	Ervatao	23	44	0.5227273
2004	Ponta Cosme	65	87	0.7471264

Utilizando ahora pivot_wider() podemos calcular la diferencia en proporciones entre playas cada año:

porAño %>% 
  select(Año,playa,pE) %>% 
  pivot_wider(names_from = playa,values_from = pE) %>% 
  mutate(dif=Ervatao-`Ponta Cosme`) %>% 
  flextable() %>% 
  colformat_num(2:4,digits=4) %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Año	Ervatao	Ponta Cosme	dif
1999	0.4167	0.7234	-0.3067
2000	0.6610	0.6724	-0.0114
2001	0.6585	0.7966	-0.1381
2002	0.5714	0.7222	-0.1508
2003	0.4595	0.6667	-0.2072
2004	0.5227	0.7471	-0.2244

Como vemos, todos los años la proporción de nidos con cangrejos ha sido mayor en Ponta Cosme que en Ervatao. Para calcular un intervalo de confianza para la diferencia entre proporciones cada año, podemos utilizar la misma sintaxis que empleamos en el ejercicio anterior simplemente añadiendo el año como variable de agrupación. En primer lugar agrupamos por año y playa para obtener el total de nidos y el total de nidos afectados por cangrejos:

tortugas %>% 
  filter(!is.na(cangrejos)&(playa=="Ponta Cosme"|playa=="Ervatao")) %>%
  group_by(Año,playa) %>% 
  summarise(nE=sum(cangrejos),n=n())

## # A tibble: 12 x 4
## # Groups:   Año [6]
##      Año playa          nE     n
##    <dbl> <chr>       <dbl> <int>
##  1  1999 Ervatao        25    60
##  2  1999 Ponta Cosme    34    47
##  3  2000 Ervatao        39    59
##  4  2000 Ponta Cosme    39    58
##  5  2001 Ervatao        27    41
##  6  2001 Ponta Cosme    47    59
##  7  2002 Ervatao        28    49
##  8  2002 Ponta Cosme    39    54
##  9  2003 Ervatao        17    37
## 10  2003 Ponta Cosme    34    51
## 11  2004 Ervatao        23    44
## 12  2004 Ponta Cosme    65    87

Ahora deberemos volver a agrupar por año para que construya los intervalos de confianza a partir de la tabla anterior:

tortugas %>% 
  filter(!is.na(cangrejos)&(playa=="Ponta Cosme"|playa=="Ervatao")) %>%
  group_by(Año,playa) %>% 
  summarise(nE=sum(cangrejos),n=n()) %>%
  group_by(Año) %>% 
  do(data.frame(BinomDiffCI(.$nE[1],.$n[1],.$nE[2],.$n[2],conf.level=0.95))) %>% 
  do(formatInterval(.$est,.$lwr.ci,.$upr.ci,label="Difference")) %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Año	Difference	CI95%
1999	-0.3067	(-0.4712, -0.1186)
2000	-0.0114	(-0.1796, 0.1577)
2001	-0.1381	(-0.3114, 0.0399)
2002	-0.1508	(-0.3259, 0.0346)
2003	-0.2072	(-0.4007, 0.003)
2004	-0.2244	(-0.3904, -0.0492)

1.2. Intervalo de confianza para el cociente de proporciones

Ejercicio 3

Construir e interpretar un intervalo de confianza al 95% para en cociente entre las proporciones de nidos afectados por cangrejos en las playas de Ervatao y Ponta Cosme.

Ya hemos visto que la proporción muestral de nidos afectados por cangrejos en Ervatao es del 58.83% y en Ponta Cosme es el 72.47%. Su cociente es por tanto:

\[\frac{0.5483}{0.7247}=0.7565\]

lo que significa que la proporción de nidos afectados por cangrejos en Ervatao es un 75.66% de la proporción de los afectados en Ponta Cosme. Para obtener un intervalo de confianza para este cociente usamos la función BinomRatioCI():

BinomRatioCI(159,290,258,356,conf.level=0.95)

##            est    lwr.ci    upr.ci
## [1,] 0.7565357 0.6692946 0.8551485

Podemos mejorar la presentación usando la misma sintaxis de antes, simplemente cambiando BinomDiffCI() por BinomRatioCI():

tortugas %>% 
  filter(!is.na(cangrejos)&(playa=="Ponta Cosme"|playa=="Ervatao")) %>%
  group_by(playa) %>% 
  summarise(nE=sum(cangrejos),n=n()) %>% 
  do(data.frame(BinomRatioCI(.$nE[1],.$n[1],.$nE[2],.$n[2],conf.level=0.95))) %>% 
  do(formatInterval(.$est,.$lwr.ci,.$upr.ci, label="Prop.Ratio")) %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Prop.Ratio	CI95%
0.7565	(0.6693, 0.8551)

1.3. Intervalo de confianza para la diferencia de medias

Ejercicio 4

Construir un intervalo de confianza al 99% para la diferencia entre las profundidades medias de los nidos excavados en la playa de Ervatao y los excavados en la playa de Ponta Cosme.

La profundidad media de los nidos en cada playa puede obtenerse mediante:

tortugas %>% 
  filter(playa=="Ervatao"|playa=="Ponta Cosme") %>% 
  group_by(playa) %>% 
  summarize(profMedia=mean(profNido))

## # A tibble: 2 x 2
##   playa       profMedia
##   <chr>           <dbl>
## 1 Ervatao          58.9
## 2 Ponta Cosme      45.1

Vemos así que en Ervatao la profundidad media es 58.9 cm, mayor que los 45.1 cm de media que se observan en Ponta Cosme. La diferencia entre ambas es, pues, de $58.9-45.1=13.8$ cm

Para construir un intervalo de confianza para la diferencia de medias utilizamos la función t.test. En la sintaxis básica de esta función hay que indicar el conjunto de datos sobre el que se realiza el test, que en este caso es el formado por los nidos de las dos playas indicadas:

t.test(profNido~playa,
       data=subset(tortugas,playa=="Ervatao"|playa=="Ponta Cosme"),
       conf.level = 0.99)

## 
##  Welch Two Sample t-test
## 
## data:  profNido by playa
## t = 32.353, df = 914.65, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 99 percent confidence interval:
##  12.67380 14.87141
## sample estimates:
##     mean in group Ervatao mean in group Ponta Cosme 
##                  58.88425                  45.11165

Esta salida nos muestra el valor medio en cada playa, así como el intervalo de confianza al 99%, que es $[12.67, 14.87]$. Este intervalo nos informa de que, aunque en la muestra disponible la diferencia observada entre las medias ha sido de 13.8 cm, lo más que podemos llegar a asegurar con una confianza del 99% es que la diferencia media entre las profundidades de los nidos en ambas playas está entre 12.67 y 14.87 cm.

Podemos integrar este resultado en una cadena de comandos de tidyverse, para mejorar automáticamente la presentación final del resultado. Observemos primero los nombres de las variables generadas tras calcular el t.test a través de la función tidy:

tortugas %>% 
  filter(playa=="Ponta Cosme"|playa=="Ervatao") %>%
  do(tidy(t.test(.$profNido~.$playa,conf.level=0.99)))

## # A tibble: 1 x 10
##   estimate estimate1 estimate2 statistic   p.value parameter conf.low
##      <dbl>     <dbl>     <dbl>     <dbl>     <dbl>     <dbl>    <dbl>
## 1     13.8      58.9      45.1      32.4 1.10e-153      915.     12.7
## # … with 3 more variables: conf.high <dbl>, method <chr>,
## #   alternative <chr>

Como podemos ver, el estimador de la diferencia se llama estimate y los límites del intervalo son conf.low y conf.high. Debemos introducir estos nombres al llamar a la función formatInterval():

tortugas %>% 
  filter(playa=="Ponta Cosme"|playa=="Ervatao") %>%
  do(tidy(t.test(.$profNido~.$playa,conf.level=0.99))) %>% 
  do(formatInterval(.$estimate,.$conf.low,.$conf.high,conf.level = 0.99, label="Mean Difference")) %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Mean Difference	CI99%
13.7726	(12.6738, 14.8714)

Ejercicio 5

Calcula las medias de las siguientes variables en las playas de Ervatao y Ponta Cosme y construye un intervalo de confianza al 95% para sus diferencias. Muestra todos los resultados en una única tabla.

Distancia hasta la orilla
Número de huevos
Peso de la tortuga nidificante

En principio vamos a resolver cada apartado por separado:

Distancia hasta la orilla

Como el problema nos pide la media en cada playa, podemos observar que t.test ya nos devuelve esos valores:

tortugas %>% 
  filter(playa=="Ervatao"|playa=="Ponta Cosme") %>%
  do(tidy(t.test(.$distancia~.$playa,conf.level=0.95)))

## # A tibble: 1 x 10
##   estimate estimate1 estimate2 statistic p.value parameter conf.low
##      <dbl>     <dbl>     <dbl>     <dbl>   <dbl>     <dbl>    <dbl>
## 1    -2.77      23.6      26.4     -5.80 8.95e-9      908.    -3.71
## # … with 3 more variables: conf.high <dbl>, method <chr>,
## #   alternative <chr>

La media en la primera playa se ha guardado como estimate1 y en la segunda playa como estimate2. Podemos incluir estos valores en la cadena de comandos para obtener el intervalo, y así mostramos no sólo la diferencia de medias y el intervalo de confianza, sino también las medias en cada playa. Incluso podemos añadir el nombre de la variable para mayor claridad:

mediaDist <- tortugas %>% 
  filter(playa=="Ponta Cosme"|playa=="Ervatao") %>%
  do(tidy(t.test(.$distancia~.$playa,conf.level=0.95)))%>% 
  do(tibble(`Variable`= "Distancia media a la orilla",
            `Ervatao`=.$estimate1,
            `Ponta Cosme `=.$estimate2,
            formatInterval(.$estimate,.$conf.low,.$conf.high,
                           conf.level = 0.95, label="Mean Difference")))

mediaDist %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Variable	Ervatao	Ponta Cosme	Mean Difference	CI95%
Distancia media a la orilla	23.59093	26.36104	-2.7701	(-3.7068, -1.8334)

Número de huevos

mediaHuevos <- tortugas %>% 
  filter(playa=="Ponta Cosme"|playa=="Ervatao") %>%
  do(tidy(t.test(.$Huevos~.$playa,conf.level=0.95)))%>% 
  do(tibble(`Variable`= "Número medio de Huevos",
            `Ervatao`=.$estimate1,
            `Ponta Cosme `=.$estimate2,
            formatInterval(.$estimate,.$conf.low,.$conf.high,
                           conf.level = 0.95, label="Mean Difference")))
mediaHuevos %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Variable	Ervatao	Ponta Cosme	Mean Difference	CI95%
Número medio de Huevos	83.95943	82.80321	1.1562	(-0.4882, 2.8007)

Peso de la tortuga nidificante

mediaPeso <- tortugas %>% 
  filter(playa=="Ponta Cosme"|playa=="Ervatao") %>%
  do(tidy(t.test(.$peso~.$playa,conf.level=0.95)))%>% 
  do(tibble(`Variable`= "Peso medio",
            `Ervatao`=.$estimate1,
            `Ponta Cosme `=.$estimate2,
            formatInterval(.$estimate,.$conf.low,.$conf.high,
                           conf.level = 0.95, label="Mean Difference"))) 
mediaPeso %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Variable	Ervatao	Ponta Cosme	Mean Difference	CI95%
Peso medio	60.03866	60.86546	-0.8268	(-1.5107, -0.1429)

Podemos usar bind-rows() para juntar los tres resultados anteriores en una única tabla:

tablaMedias <- bind_rows(mediaDist,mediaHuevos,mediaPeso)
tablaMedias %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Variable	Ervatao	Ponta Cosme	Mean Difference	CI95%
Distancia media a la orilla	23.59093	26.36104	-2.7701	(-3.7068, -1.8334)
Número medio de Huevos	83.95943	82.80321	1.1562	(-0.4882, 2.8007)
Peso medio	60.03866	60.86546	-0.8268	(-1.5107, -0.1429)

Con una sintaxis un poco más elaborada, a la vez que más simple, podemos generar la tabla con las tres variables con una única cadena de comandos:

tortugas %>% 
  filter(playa=="Ponta Cosme"|playa=="Ervatao") %>%
  select(playa,distancia,Huevos,peso) %>% 
  pivot_longer(2:4,names_to = "Variable", values_to = "valor") %>% 
  group_by(Variable) %>%
  do(tidy(t.test(.$valor~.$playa,conf.level=0.95)))%>% 
  do(tibble(`Ervatao`=.$estimate1,
            `Ponta Cosme `=.$estimate2,
            formatInterval(.$estimate,.$conf.low,.$conf.high,
                           conf.level = 0.95, label="Mean Difference"))) %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Variable	Ervatao	Ponta Cosme	Mean Difference	CI95%
distancia	23.59093	26.36104	-2.7701	(-3.7068, -1.8334)
Huevos	83.95943	82.80321	1.1562	(-0.4882, 2.8007)
peso	60.03866	60.86546	-0.8268	(-1.5107, -0.1429)

Ejercicio 6

Explora si la diferencia entre las profundidades medias de los nidos excavados en la playa de Ervatao y los excavados en la playa de Ponta Cosme ha variado durante los años en que se han realizado las campañas. Calcula un intervalo de confianza para esta diferencia cada año.

Igual que hicimos en el ejercicio 2 con las proporciones, podemos en primer lugar calcular la profundidad media por playa y año, y su diferencia:

tortugas %>% 
  filter(playa=="Ervatao"|playa=="Ponta Cosme") %>% 
  group_by(Año,playa) %>% 
  summarize(profMedia=mean(profNido)) %>% 
  pivot_wider(names_from = playa, values_from = profMedia) %>% 
  mutate(diferencia=Ervatao-`Ponta Cosme`) %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  colformat_num(2:4,digits=2) %>% 
  autofit()

Año	Ervatao	Ponta Cosme	diferencia
1999	58.37	45.01	13.36
2000	58.90	45.39	13.51
2001	58.58	45.23	13.35
2002	58.75	45.41	13.34
2003	59.09	44.35	14.73
2004	59.66	45.16	14.50

Como vemos, la profundidad media de los nidos se mantiene bastante constante a lo largo de las distintas campañas, siendo siempre mayor en Ervatao que en Ponta Cosme. Para calcular un intervalo de confianza para cada año usamos la misma sintaxis del ejercicio 4, pero añadiendo una agrupación por año:

tortugas %>% 
  filter(playa=="Ponta Cosme"|playa=="Ervatao") %>%
  group_by(Año) %>% 
  do(tidy(t.test(.$profNido~.$playa,conf.level=0.99))) %>% 
  do(formatInterval(.$estimate,.$conf.low,.$conf.high,conf.level = 0.99)) %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit() %>% 
  set_caption("Diferencia media en la profundidad de los nidos de las playas de Ervatao y Ponta Cosme")

Diferencia media en la profundidad de los nidos de las playas de Ervatao y Ponta Cosme
Año	Estimate	CI99%
1999	13.3625	(10.7777, 15.9473)
2000	13.5107	(11.0018, 16.0196)
2001	13.3502	(10.12, 16.5803)
2002	13.3420	(10.4296, 16.2545)
2003	14.7312	(12.2192, 17.2431)
2004	14.4991	(11.854, 17.1442)

1.4. Intervalo de confianza para el cociente de varianzas o desviaciones típicas

Ejercicio 7

Calcula la desviación típica de la distancia desde el nido hasta la orilla en las playas de Ervatao y Ponta Cosme. Calcula el cociente entre ambas y construye un intervalo de confianza al 95% para dicho cociente. ¿Muestra ese intervalo alguna evidencia de que dicha variabilidad difiere entre estas dos playas?

Para calcular la desviación típica en cada playa procedemos de modo similar a como hemos hecho en el ejercicio 4 con las medias:

tortugas %>% 
  filter(playa=="Ervatao"|playa=="Ponta Cosme") %>% 
  group_by(playa) %>% 
  summarize(sdDistancia=sd(distancia))

## # A tibble: 2 x 2
##   playa       sdDistancia
##   <chr>             <dbl>
## 1 Ervatao            6.92
## 2 Ponta Cosme        7.52

La variabilidad en Ponta Cosme es ligeramente mayor que el Ervatao. El cociente entre ambas desviaciones típicas es:

\[\frac{6.92}{7.52}=0.92\]

Podemos obtener un intervalo de confianza para dicho cociente mediante la función var.test():

vt <- var.test(distancia~playa,
         data=subset(tortugas,playa=="Ervatao"|playa=="Ponta Cosme"),
         conf.level=0.95)
vt

## 
##  F test to compare two variances
## 
## data:  distancia by playa
## F = 0.84776, num df = 418, denom df = 497, p-value = 0.07996
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
##  0.7056736 1.0199654
## sample estimates:
## ratio of variances 
##          0.8477589

En realidad el intervalo que hemos obtenido es para las varianzas; para convertirlo en un intervalo para las desviaciones típicas solamente hay que hacer la raiz cuadrada; téngase en cuenta que el resultado de var.test() es una lista, y el elemento que contiene al intervalo de confianza es conf.int:

sqrt(vt$conf.int)

## [1] 0.8400438 1.0099334
## attr(,"conf.level")
## [1] 0.95

Así pues, con un 95% de confianza podemos decir que el cociente entre las desviaciones típicas poblacionales de las distancias de los nidos a la orilla en las dos playas es un valor comprendido entre 0.84 y 1.0099; por tanto, a este nivel de confianza no es descartable la posibilidad de que dicho cociente sea 1, lo que significa que no es descartable la posibilidad de que las desviaciones típicas poblacionales sean iguales en ambas playas, y que la diferencia que se ha observado en sus valores muestrales sea debida exclusivamente al error muestral (la variabilidad natural que se observa al tomar muestras aleatorias).

Mediante el uso de do() y tidy() podemos obtener directamente el intervalo formateado y listo para presentación. En este caso, como el intervalo es para las desviaciones típicas, debemos tener la precaución de calcular las raíces cuadradas de los términos:

tortugas %>% 
  filter(playa=="Ervatao"|playa=="Ponta Cosme") %>%
  do(tidy(var.test(.$distancia~.$playa,conf.level=0.95))) %>% 
  do(formatInterval(sqrt(.$estimate),sqrt(.$conf.low),sqrt(.$conf.high))) %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Estimate	CI95%
0.9207	(0.84, 1.0099)

2. Muestras emparejadas

Para los ejercicios de esta sección utilizaremos el archivo luzNidos.csv que puede descargarse del campus virtual o directamente desde internet mediante:

luzNidos <- read_csv2("http://estadistica-dma.ulpgc.es/estadFCM/datos/luzNidos.csv")

En este archivo se recogen datos de 62 tortugas con las que se realizó una experiencia controlada para determinar si la profundidad del nido y/o el tamaño de la puesta están relacionados con la iluminación de la playa. Para ello, se identificaron todas las tortugas que anidaron más de una vez en las playas de Ervatao y Ponta Cosme en el año 2010; de éstas, se seleccionaron 62 que anidaron dos veces en zonas de similares características (tipo de arena, humedad, temperatura), pero una vez con la playa iluminada por la luna llena y otra vez con la playa oscura por no haber luna en el cielo.

El archivo contiene las siguientes variables:

id: número de identificación de la tortuga, de 1 a 62
playa: la playa donde se hizo el nido
luna: variable con dos valores, “sí” y “no”, que indica si el nido se hizo con luna llena o no.
profNido: profundidad del nido
Huevos: tamaño de la puesta (número de huevos)

2.1. Intervalo de confianza para la diferencia de proporciones

Ejercicio 8

Construir e interpretar un intervalo de confianza al 95% para la diferencia entre la proporción de tortugas que han anidado en Ervatao con luna llena y la proporción de las que han anidado en la misma playa cuando no hay luna.

Podemos calcular dichas proporciones fácilmente mediante:

luzNidos %>% tabyl(luna,playa) %>%
  adorn_totals("col") %>% 
  adorn_percentages() %>% 
  adorn_ns() %>% 
  flextable() %>% 
  autofit()

luna	Ervatao	Ponta Cosme	Total
no	0.2258065 (14)	0.7741935 (48)	1 (62)
si	0.3387097 (21)	0.6612903 (41)	1 (62)

Es decir, con luna llena anidaron en Ervatao 21 de las 62 tortugas (un 33.87%); sin luna anidaron 14 de las 64 (un 22.58%). Por tanto la diferencia es de un 33.87-22.58=11.29%. Para construir un intervalo de confianza para esta diferencia, debemos tener en cuenta que son las mismas 62 tortugas las que se han observado dos veces; por tanto estas muestras están emparejadas y para construir un intervalo de confianza para la diferencia debemos construir una tabla que nos indique cuántas de estas tortugas volvieron a la misma playa y cuántas anidaron en playas distintas. Para ello creamos un nuevo data.frame con 62 filas (una por cada tortuga) y dos columnas, una con la playa de anidamiento con luna llena y otra con la playa de anidamiento sin luna:

luzNidos %>% 
  select(id,luna,playa) %>% 
  pivot_wider(values_from = playa, names_from = luna)

## # A tibble: 62 x 3
##       id si          no         
##    <dbl> <chr>       <chr>      
##  1     1 Ponta Cosme Ponta Cosme
##  2     2 Ponta Cosme Ervatao    
##  3     3 Ponta Cosme Ervatao    
##  4     4 Ponta Cosme Ponta Cosme
##  5     5 Ervatao     Ervatao    
##  6     6 Ervatao     Ervatao    
##  7     7 Ponta Cosme Ponta Cosme
##  8     8 Ervatao     Ponta Cosme
##  9     9 Ponta Cosme Ponta Cosme
## 10    10 Ponta Cosme Ponta Cosme
## # … with 52 more rows

Ahora podemos usar tabyl para cruzar estas variables:

luzNidos %>% 
  select(id,luna,playa) %>% 
  pivot_wider(values_from = playa, names_from = luna) %>% 
  tabyl(si,no)

##           si Ervatao Ponta Cosme
##      Ervatao       4          17
##  Ponta Cosme      10          31

Así, de las 21 tortugas que anidaron en Ervatao con luna llena, solo 4 anidaron también en Ervatao sin luna; las restantes 17 anidaron en Ponta Cosme; asimismo de las 41 que anidaron con luna llena en Ponta Cosme, 31 volvieron a anidar en esta playa sin luna y las otras 10 se fueron a Ervatao.

El intervalo de confianza para la diferencia entre as proporciones de las que anidaron en Ervatao con y sin luna llena se calcula mediante la función diffpropci.mp del paquete PropCIs, donde debemos colocar solamente los números de tortugas que cambiaron de playa (17 de Ervatao a Ponta Cosme y 10 de Ponta Cosme a Ervatao), y el total de tortugas (62):

diffpropci.mp(10,17,62,conf.level=0.95)$conf.int

## [1] -0.05044342  0.26919342
## attr(,"conf.level")
## [1] 0.95

2.2. Intervalo de confianza para la diferencia de medias

Ejercicio 9

Construir e interpretar un intervalo de confianza al 95% para la diferencia entre la profundidad media de los nidos excavados con luna llena y los excavados cuando no hay luna.

La profundidad media de los nidos según haya o no luna llena puede obtenerse mediante:

luzNidos %>% 
  group_by(luna) %>% 
  summarize(mean(profNido))

## # A tibble: 2 x 2
##   luna  `mean(profNido)`
##   <chr>            <dbl>
## 1 no                56.3
## 2 si                51.3

Así pues, en esta muestra de 62 tortugas, cuando hay luna llena la profundidad media de los nidos es de 51.3 cm; cuando no hay luna dicha profundidad es de 56.3 cm, es decir, hay una diferencia media de 5 cm según haya o no luna. Para obtener un intervalo de confianza utilizamos t.test() con la opción paired=TRUE (NOTA: el uso de esta sintaxis requiere que los datos de profundidad del nido estén correctamente ordenados en el mismo orden en el grupo con luna y en el grupo sin luna, ya que en caso contrario el emparejamiento no se haría bien pues la función t.test en ningún caso comprueba el número de identificación para hacer las parejas):

t.test(profNido~luna,data=luzNidos,paired=TRUE,conf.level=0.95)

## 
##  Paired t-test
## 
## data:  profNido by luna
## t = 10.64, df = 61, p-value = 1.578e-15
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  4.072079 5.956953
## sample estimates:
## mean of the differences 
##                5.014516

Por tanto, con un 95% de confianza, la diferencia media de profundidad entre los nidos excavados con y sin luna llena se sitúa entre 4.07 y 5.96 cm. Esto significa que podemos asegurar con un 95% de confianza que en la oscuridad las tortugas hacen nidos en promedio más profundos que cuando hay luz.

Si no hubiésemos utilizado la opción paired=TRUE el intervalo habría sido diferente (e incorrecto):

t.test(profNido~luna,data=luzNidos,conf.level=0.95)

## 
##  Welch Two Sample t-test
## 
## data:  profNido by luna
## t = 3.7062, df = 121.71, p-value = 0.0003179
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  2.336061 7.692971
## sample estimates:
## mean in group no mean in group si 
##         56.30323         51.28871

Sintaxis alternativa para el intervalo para la diferencia de medias emparejadas.

En muchas ocasiones los datos de las muestras emparejadas están dispuestos en columnas. Por ejemplo, los datos de profundidad de los nidos según que haya luna o no podrían estar en un archivo csv como el de la imagen:

Leemos el archivo:

profNidos <- read_csv2("archivoNidos.csv")

Para comparar las medias de las profundidades de los nidos con y sin luna simplemente hay que especificar los nombres de las dos variables en el comando t.test(); en este caso el nombre de la variable debe especificarse junto con el data.frame que la contiene, separado con el símbolo $:

t.test(profNidos$conLuna, profNidos$sinLuna, conf.level=0.95)

## 
##  Welch Two Sample t-test
## 
## data:  profNidos$conLuna and profNidos$sinLuna
## t = -3.7062, df = 121.71, p-value = 0.0003179
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -7.692971 -2.336061
## sample estimates:
## mean of x mean of y 
##  51.28871  56.30323

Como vemos, hemos obtenido exactamente el mismo resultado de antes. Este procedimiento tiene además la garantía de que, si los datos se han introducido correctamente en el archivo, los valores están correctamente emparejados.

Comprobación del emparejamiento

Si recibimos los datos en la primera de las formas, esto es, una variable que contiene la identificación de la tortuga, una segunda variable con el grupo al que pertenece (con luna o sin luna) y una tercera variable con la profundidad de los nidos, pero no estamos seguros de que estén ordenados correctamente para llevar a cabo el t.test, podemos utilizar la siguiente sintaxis para ordenar los datos en columnas debidamente emparejadas y guardarlos en el data.fame nidos:

nidos <- luzNidos %>% select(id,luna,profNido) %>% 
  pivot_wider(names_from = luna,values_from = profNido)
nidos

## # A tibble: 62 x 3
##       id    si    no
##    <dbl> <dbl> <dbl>
##  1     1  50    53.1
##  2     2  54.3  58.8
##  3     3  49.1  52.7
##  4     4  47.9  60.8
##  5     5  65.5  66.2
##  6     6  49.5  54.4
##  7     7  54.8  59.5
##  8     8  53.4  57.3
##  9     9  46    49.2
## 10    10  56.7  57.9
## # … with 52 more rows

Obsérvese que de esta forma cada identificador de tortuga se acompaña de la profundidad del nido que se realizó con luna y el que se realizó sin luna. Ahora, para hacer el t.test podemos utilizar la sintaxis siguiente:

t.test(nidos$si,nidos$no,conf.level=0.95)

## 
##  Welch Two Sample t-test
## 
## data:  nidos$si and nidos$no
## t = -3.7062, df = 121.71, p-value = 0.0003179
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -7.692971 -2.336061
## sample estimates:
## mean of x mean of y 
##  51.28871  56.30323

También podemos combinar los dos chunks anteriores y dar formato al intervalo:

luzNidos %>% select(id,luna,profNido) %>% 
  pivot_wider(names_from = luna,values_from = profNido) %>% 
  do(tidy(t.test(.$si,.$no,conf.level=0.95))) %>% 
  do(formatInterval(.$estimate,.$conf.low,.$conf.high,label="Mean difference")) %>% 
  flextable() %>% 
  fontsize(size=14, part="all") %>% 
  autofit()

Mean difference	CI95%
-5.0145	(-7.693, -2.3361)

2.3. Intervalo de confianza para el cociente de varianzas o desviaciones típicas.

Ejercicio 10

Construir e interpretar un intervalo de confianza al 95% para el cociente entre la desviación típica del tamaño de la puesta cuando hay luna llena y cuando no. A partir de dicho intervalo, ¿podría decirse que la variabilidad en el tamaño de la puesta es diferente segíun que haya o no luna llena?

La desviación típica del tamaño de la puesta según haya o no luna llena puede obtenerse mediante:

luzNidos %>% 
  group_by(luna) %>% 
  summarize(sd(Huevos))

## # A tibble: 2 x 2
##   luna  `sd(Huevos)`
##   <chr>        <dbl>
## 1 no            7.68
## 2 si            8.58

Así pues, en esta muestra de 62 tortugas, cuando hay luna llena la variabilidad del tamaño de la puesta (sd=8.58) es mayor que cuando no hay luna (sd=7.68). El cociente entre ambas desviaciones típicas es 8.58/7.68 = 1.117. Como se trata de muestras emparejadas, para obtener un intervalo de confianza para este cociente utilizamos la función Var.test() del paquete PairedData, indicando la opción paired=TRUE. Esta función requiere que se le pasen los valores emparejados en dos variables; por ello, primero usamos pivot_wider() para generar un nuevo data.frame con esa característica:

tamPuesta <- luzNidos %>% select(id,luna,Huevos) %>% 
  pivot_wider(names_from = luna,values_from = Huevos)
tamPuesta

## # A tibble: 62 x 3
##       id    si    no
##    <dbl> <dbl> <dbl>
##  1     1    75    70
##  2     2    75    72
##  3     3    93    90
##  4     4    85    87
##  5     5    88    86
##  6     6    90    71
##  7     7    89    77
##  8     8    87   101
##  9     9    75    87
## 10    10    78    78
## # … with 52 more rows

con lo que hemos generado una variable llamada “si” y otra “no”; son estas dos variables las que comparamos mediante Var.test() (esta función no admite el argumento data y por tanto hay que especificar las variables con el nombre del dataframe al que pertenecen, separadas con el símbolo $):

Var.test(tamPuesta$si,tamPuesta$no,paired=TRUE,conf.level=0.95)

## 
##  Paired Pitman-Morgan test
## 
## data:  x and y
## t = 0.86757, df = 60, p-value = 0.3891
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
##  0.7508828 2.0796888
## sample estimates:
## variance of x variance of y 
##      73.62374      58.91592

Por tanto, con un 95% de confianza, el cociente entre las varianzas del tamaño de la puesta con y sin luna llena está entre 0.75 y 2.08.

Si queremos el intervalo para el cociente de desviaciones típicas, calculamos la raíz cuadrada del intervalo anterior:

sqrt(Var.test(tamPuesta$si,tamPuesta$no,
              paired=TRUE,conf.level=0.95)$conf.int)

## [1] 0.8665349 1.4421126
## attr(,"conf.level")
## [1] 0.95

Como el valor 1 cae dentro de este intervalo, no es descartable que las desviaciones típicas (o las varianzas) poblacionales sean iguales.

Práctica 9. Intervalos de Confianza para dos muestras con R

Intervalos para medias, varianzas y proporciones en muestras independientes y emparejadas