A continuación vamos a realizar inferencia en un conjunto de datos. Más concretamente, vamos a realizar análisis sobre la media y varianza de una y dos poblaciones.
Para comparar medias y varianzas se disponen de los comandos t.test y var.test, respectivamente:
t.test(var1, var2=NULL, alternative=“two.sided”, mu = 0, paired=FALSE, var.equal=FALSE, conf.level=0.95)
var.test(var1, var2, ratio=1, alternative = “two.sided”, conf.level=0.95)
Las opciones indicadas son las ofrecidas por defecto:
Si para t.test no se especifica el segundo conjunto de datos se trabajará con una única muestra.
El argumento alternative indica el tipo de contraste: bilateral si se especifica two.sided, unilateral con hipótesis alternativa del tipo mayor para greater y menor para less.
En el argumento mu se especifica el valor de la hipótesis nula en el caso de comparación de medias, mientras que para la comparación de varianzas se hace con el argumento ratio.
En el argumento paired se especifica si los datos están relacionados (pareados), True, o no, False.
En el argumento var.equal se especifica, en caso de dos muestras, si las varianzas son iguales (True) o no (False).
Con el argumento conf.level se indica el nivel de confianza a usar.
Al mismo tiempo, es importante manejar el concepto de p-valor. Se define el p-valor como el mínimo valor de significación al partir del cual se rechaza la hipótesi nula. Por tanto:
Si p-valor < 0.1, se rechaza la hipótesis nula al 10% de significación.
Si p-valor < 0.05, se rechaza la hipótesis nula al 5% de significación.
Si p-valor < 0.01, se rechaza la hipótesis nula al 1% de significación.
Como se verá a lo largo de este post, el p-valor resultará muy útil para interpretar los resultados obtenidos.
Conjunto de datos existentes en R
En un post anterior trataba la lectura de datos desde una hoja de cálculo o un archivo de texto plano. Otra posibilidad es usar los datos que vienen cargados en R. Para saber los datos que tenemos diponibles simplemente hay que ejecutar data() en la consola. O de forma más exhaustiva data(package = .packages(all.available = TRUE)).
Veamos algunos ejemplos que ilustren las órdenes anteriores usando los datos de la base cars (para más información help(cars)):
cars # quitar primera almohadilla para ver el conjunto de datos
Con la segunda orden puedo usar directamente el nombre de las variables: speed (velocidad en millas por hora) y dist (distancia de frenado en pies).
Inferencia para una población
Para contrastar si la velocidad media es igual, mayor o menor a 17 millas por hora ejecutaríamos:
t.test(speed, mu=17) # por defecto es a dos colas
One Sample t-test
data: speed
t = -2.1397, df = 49, p-value = 0.03739
alternative hypothesis: true mean is not equal to 17
95 percent confidence interval:
13.89727 16.90273
sample estimates:
mean of x
15.4
t.test(speed, mu=17, alternative="greater") # como alternativa la suposición que nos interesa
One Sample t-test
data: speed
t = -2.1397, df = 49, p-value = 0.9813
alternative hypothesis: true mean is greater than 17
95 percent confidence interval:
14.1463 Inf
sample estimates:
mean of x
15.4
t.test(speed, mu=17, alternative="less")
One Sample t-test
data: speed
t = -2.1397, df = 49, p-value = 0.01869
alternative hypothesis: true mean is less than 17
95 percent confidence interval:
-Inf 16.6537
sample estimates:
mean of x
15.4
Atendiendo al p-valor, rechazo que la media sea igual a 17 frente a la alternativa de que sea distinta (p-valor=0.03739<0.05), no rechazo que la media sea 17 frente a la alternativa de que sea mayor (p-valor=0.9813>0.05) y rechazo la hipótesis de que sea igual a 17 frente a la alternativa de que sea menor (p-valor=0.01869<0.05).
Adviértase que en todos los casos se nos proporciona el correspondiente intervalo de confianza, obteniéndose los mismos resultados, y el estimador puntual de la media (media muestral del conjunto de datos).
La decisión se podría haber tomado también a partir de la región de rechazo tradicional. Simplemente habría que comparar el valor experimental con el teórico. Así, por ejemplo, en el caso del contraste a dos colas, el valor teórico es:
qt(0.975, 49)
[1] 2.009575
Como el valor experimental (en valor absoluta), 2.1397, es mayor que el teórico, 2.009575, se rechaza la hipótesis nula de que la media es igual a 17. Como no podría ser de otra forma, ya que ya sea mediante esta opción o la anterior, siempre tenemos que llegar a la misma conclusión.
En los otros dos casos, los valores teóricos son, respectivamente:
qt(0.95, 49)
[1] 1.676551
qt(0.05, 49)
[1] -1.676551
Como -2.1397<1.676551, no se rechaza la hipótesis nula en el primer caso; mientras que como -2.1397<-1.676551, se rechaza la hipótesis nula en el segundo.
Inferencia para dos poblaciones
Empecemos viendo si speed y dist tienen la misma varianza. Observando las cuasivarianzas muestrales se debería verificar que la primera es menor que la segunda, por tal motivo nos interesamos por esta hipótesis:
var(speed)
[1] 27.95918
var(dist)
[1] 664.0608
var.test(speed, dist)
F test to compare two variances
data: speed and dist
F = 0.042103, num df = 49, denom df = 49, p-value < 2.2e-16
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
0.02389265 0.07419404
sample estimates:
ratio of variances
0.04210335
var.test(speed, dist, alternative="l") # con la inicial es suficiente
F test to compare two variances
data: speed and dist
F = 0.042103, num df = 49, denom df = 49, p-value < 2.2e-16
alternative hypothesis: true ratio of variances is less than 1
95 percent confidence interval:
0.00000000 0.06767227
sample estimates:
ratio of variances
0.04210335
Atendiendo a los resultados, en primer lugar se rechaza la hipótesis nula de que sean iguales frente a la alternativa de que sean distintas (p-valor=2.2·10^(-16)<0.05) y en segundo lugar también se rechaza frente a la alternativa de que sea menor (p-valor=2.2·10^(-16)<0.05). Adviértase que con el primer contraste sólo se podría afirmar (al nivel de significación usado) que las varianzas son distintas, sin embargo, con el segundo se podría decir que la primera varianza es menor que la segunda.
Al igual que antes también se obtienen los correspondientes intervalos de confianza, a partir de los cuales se obtendrían las mismas conclusiones (adviértase que el segundo es unilateral).
También se podrían haber usado la comparación de los valores experimentales con los teóricos para tomar una decisión en el contraste. En el caso del contraste a dos colas para el cociente de varianzas, los valores teóricos son:
qf(0.025,49,49)
[1] 0.5674762
qf(0.975,49,49)
[1] 1.762189
Entonces, como 0.042103 es menor que 0.5674762, entonces se rechaza la hipótesis nula de igualdad de varianza.
Para estudiar si speed y dist tienen la misma media habría que especificar que las varianzas son distintas:
t.test(speed, dist, var.equal=F)
Welch Two Sample t-test
data: speed and dist
t = -7.4134, df = 53.119, p-value = 9.628e-10
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-35.04152 -20.11848
sample estimates:
mean of x mean of y
15.40 42.98
Observando las medias muestrales (se obtienen como salida del comando anterior), como era esperable, rechazo la hipótesis de igualdad (p-valor=0.000000009628). Además, atendiendo al intervalo de confianza obtenido (siempre es negativo: (-35.04152, -20.11848)) se podría decir que la media de speed es menor que la de dist.
Usando la comparación del valor experimental con el teórico:
qt(0.975,53.119)
[1] 2.005641
Se rechaza la hipótesis nula de igualdad de medias ya que 7.4134 es mayor que 2.00487.
Esta conclusión se obtiene también especificando el siguiente contraste:
t.test(speed, dist, var.equal=F, alternative="l")
Welch Two Sample t-test
data: speed and dist
t = -7.4134, df = 53.119, p-value = 4.814e-10
alternative hypothesis: true difference in means is less than 0
95 percent confidence interval:
-Inf -21.35209
sample estimates:
mean of x mean of y
15.40 42.98
Como el p-valor=0.0000000004814 es menor que 0.05 se rechaza la hipótesis nula en favor de la laternativa de menor.
Inferencia para dos poblaciones (bis)
Consideremos que se dispone información para coches europeos y americanos, la cual es codificada como sigue:
europeo =c(0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 1, 1, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 0) # europeo igual a 1 si el coche es europeo y 0 en caso contrario
Usando la variable anterior vamos dividir en dos la muestra correspondiente a speed. Sus principales características se pueden calcular como:
factores =factor(europeo) # codifico 0 y 1 como factorestapply(speed, factores, mean) # help(tapply)
0 1
15.05882 16.12500
tapply(speed, factores, var) # en lugar de mean y var se puede usar cualquier función
0 1
32.23886 19.58333
Para ver si speed tiene la misma media en los coches europeos y americanos debemos ejecutar:
var.test(speed~europeo)
F test to compare two variances
data: speed by europeo
F = 1.6462, num df = 33, denom df = 15, p-value = 0.3049
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
0.6277125 3.7222238
sample estimates:
ratio of variances
1.64624
Two Sample t-test
data: speed by europeo
t = -0.66126, df = 48, p-value = 0.5116
alternative hypothesis: true difference in means between group 0 and group 1 is not equal to 0
95 percent confidence interval:
-4.307999 2.175646
sample estimates:
mean in group 0 mean in group 1
15.05882 16.12500
En primer lugar se tiene que las varianzas pueden considerarse iguales ya que no se rechaza la hipótesis nula (p-valor=0.3049>0.05). Considerando que el intervalo calculado, (0.6277125, 3.7222238), se llegaría a la misma conclusión ya que contiene al 1, por lo que el cociente puede considerarse igual a 1 y, por tanto, las varianzas serían iguales.
Atendiendo al resultado anterior, tendríamos en cuenta el último contraste (varianzas iguales). De nuevo no se rechaza la hipótesis de que las medias sean iguales (p-valor=0.5116>0.05), es decir, la velocidad media es la misma en los coches europeos y americanos al 5% de significación. Considerando el intervalo de confianza asociado, (-4.307999, 2.175646), se llega a la misma conclusión ya que el intervalo contiene al 0.
De esta forma se puede contrastar si dos medias son iguales dependiendo de cualquier factor con dos alternativas (poseer o no una característica de interés).
Generalización a más de dos poblaciones
Si el factor tiene más de dos alternativas, entonces se ha de realizar una ANOVA de un factor:
terna =rbinom(50,2,0.5) # genero aleatoriamente un factor 0, 1, 2table(terna)
$`0`
Min. 1st Qu. Median Mean 3rd Qu. Max.
4.00 13.00 16.00 16.12 19.00 24.00
$`1`
Min. 1st Qu. Median Mean 3rd Qu. Max.
4.00 12.25 15.50 15.45 19.75 25.00
$`2`
Min. 1st Qu. Median Mean 3rd Qu. Max.
7.00 10.00 12.00 14.18 18.00 24.00
tapply(speed, factores.terna, var)
0 1 2
23.98529 28.64069 35.96364
anova.speed=aov(speed~factores.terna) # hago el análisis de la varianzarequire(car) # necesario para el test de Levene de igualdad de varianzas
Loading required package: car
Warning: package 'car' was built under R version 4.2.2
Loading required package: carData
Warning: package 'carData' was built under R version 4.2.2
leveneTest(anova.speed, center=mean) # homogeneidad de varianzas
Levene's Test for Homogeneity of Variance (center = mean)
Df F value Pr(>F)
group 2 0.514 0.6014
47
summary(anova.speed) # resumen del ANOVA
Df Sum Sq Mean Sq F value Pr(>F)
factores.terna 2 25.1 12.57 0.439 0.647
Residuals 47 1344.9 28.61
En este caso no rechazo que las varianzas sean iguales (p-valor=0.4853>0.05) y que las tres medias son iguales (p.valor=0.884>0.05). En caso de rechazarla, sería interesante ver donde falla la igualdad. Con la siguiente orden se hacen todos los posibles pares:
Se tiene que una de las hipótesis que se han de cumplir para realizar inferencia es que las poblaciones se distribuyen como una normal. Por lo que es necesario comprobar si se verifica o no este aspecto.
¿Se distribuye la variable speed según una normal? Veamos distintas opciones:
# en todos los contrastes, H0: hay normalidadks.test(speed, pnorm, mean(speed), sd(speed)) # en este caso se especifica la media y varianza que ha de tener la normal
Warning in ks.test.default(speed, pnorm, mean(speed), sd(speed)): ties should
not be present for the Kolmogorov-Smirnov test
Asymptotic one-sample Kolmogorov-Smirnov test
data: speed
D = 0.068539, p-value = 0.9729
alternative hypothesis: two-sided
shapiro.test(speed)
Shapiro-Wilk normality test
data: speed
W = 0.97765, p-value = 0.4576
# install.packages("nortest") # ejecutarlo sólo una vezlibrary(nortest) # cargamos paquetelillie.test(speed)
Lilliefors (Kolmogorov-Smirnov) normality test
data: speed
D = 0.068539, p-value = 0.8068
Si el p-valor asociado a cada uno de los contrastes es mayor que 0.05 no se rechaza la hipótesis nula de normalidad. Es decir, no se rechaza que el conjunto de datos almacenado en speed procede de una normal. En este caso, todos los p-valores son mayores que 0.05 (0.9729, 0.4576 y 0.8068, respectivamente), luego no se rechaza la hipótesi nula de normalidad.
Aleatoriedad
Otra hipótesis a comprobar es que los datos son aleatorios. Mediante el siguiente comando se ejecuta el test de rachas:
# install.packages("tseries")library(tseries)
Warning: package 'tseries' was built under R version 4.2.2
Registered S3 method overwritten by 'quantmod':
method from
as.zoo.data.frame zoo
runs.test(as.factor(speed>median(speed))) # H0: hay aleatoriedad
Runs Test
data: as.factor(speed > median(speed))
Standard Normal = -6.8583, p-value = 6.966e-12
alternative hypothesis: two.sided
Puesto que el p-valor es inferior a 0.05, se rechaza la hipótesis nula de aleatoriedad. Es decir, hay evidencias para considerar los datos no aleatorios.