Tengo un conjunto de datos de muestra con 31 valores. Ejecuté una prueba t de dos colas usando R para probar si la media verdadera es igual a 10:
t.test(x=data, mu=10, conf.level=0.95)
Resultado:
t = 11.244, df = 30, p-value = 2.786e-12 alternative hypothesis: true mean is not equal to 10 95 percent confidence interval: 19.18980 23.26907 sample estimates: mean of x 21.22944
Ahora estoy tratando de hacer lo mismo manualmente:
t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = dt(t.value, df=length(lengths-1))
El valor t calculado usando este El método es el mismo que el resultado de la función t-test R. El valor p, sin embargo, resulta ser 3.025803e-12.
¿Alguna idea de lo que estoy haciendo mal?
¡Gracias!
EDIT
Aquí está el código R completo , incluido mi conjunto de datos:
# Raw dataset -- 32 observations data = c(21.75, 18.0875, 18.75, 23.5, 14.125, 16.75, 11.125, 11.125, 14.875, 15.5, 20.875, 17.125, 19.075, 25.125, 27.75, 29.825, 17.825, 28.375, 22.625, 28.75, 27, 12.825, 26, 32.825, 25.375, 24.825, 25.825, 15.625, 26.825, 24.625, 26.625, 19.625) # Student t-Test t.test(x=data, mu=10, conf.level=0.95) # Manually calculate p-value t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = dt(t.value, df=length(data) - 1)
Respuesta
Use pt
y convertirlo en dos colas.
> 2*pt(11.244, 30, lower=FALSE) [1] 2.785806e-12
Comentarios
- I creo que falta un detalle: ¿Cuándo usar lower = F? – Consulte mi pregunta aquí: stats.stackexchange.com/questions/133091/…
- El valor debe ser positivo, por lo que si usa esto con una variable que podría ser negativa, incluya
abs
. - Para una prueba de dos colas, ' busca la probabilidad de que el valor sea menor que -11.244 o mayor que +11.244. lower = F le dice a R que calcule la probabilidad de que el valor sea mayor que el primer parámetro. De lo contrario, le da la probabilidad de que el valor sea menor que el primer parámetro. Como tal, también puede hacer 2 * pt (-11.244, 30). Personalmente, suelo hacer 2 * pt (-abs (q), df = n-1) ya que R por defecto es lower = T.
Responder
Publiqué esto como comentario, pero cuando quise agregar un poco más en la edición, se volvió demasiado largo, así que lo moví aquí.
Editar : su estadística de prueba y df son correctos. La otra respuesta señala el problema con el cálculo del área de cola en la llamada a pt()
, y la duplicación para dos colas, que resuelve su diferencia. Sin embargo, dejaré mi discusión / comentario anterior porque hace puntos relevantes de manera más general sobre los valores p en colas extremas:
Es posible que no esté haciendo nada malo y aún así obtenga una diferencia, pero si publica un ejemplo reproducible, podría ser posible investigar más a fondo si tiene algún error (digamos en el df).
Estas cosas se calculan a partir de aproximaciones que pueden no ser particularmente precisas en la cola muy extrema .
Si las dos cosas no usan aproximaciones idénticas, es posible que no estén de acuerdo de cerca, pero esa falta de acuerdo no debería importar (para que el área de cola exacta sea un número significativo, las suposiciones requeridas tienen que mantener grados asombrosos de precisión). ¿Realmente tiene la normalidad exacta, la independencia exacta, la varianza exactamente constante?
No debe esperar necesariamente una gran precisión donde los números no significan nada de todos modos. ¿Hasta qué punto importa si el valor p aproximado calculado es $ 2 \ times 10 ^ {- 12} $ o $ 3 \ times 10 ^ {- 12} $? Ninguno de los números mide el valor p real de su situación real. Incluso si uno de los números representara el valor p real de su situación real, una vez que esté por debajo de aproximadamente $ 0.0001 $, ¿por qué le importaría cuál es ese valor en realidad?
Respuesta
La mejor forma de calcularlo manualmente es:
t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = 2*pt(-abs(t.value), df=length(data)-1)
Necesita el abs () porque, de lo contrario, corre el riesgo de obtener valores p mayores que $ 1 $ (cuando la media de los datos es mayor que la media dada).
Respuesta
Me gusta mucho la respuesta que proporcionó @Aaron, junto con los abs
comentarios. Encuentro que una confirmación útil es ejecutar
pt(1.96, 1000000, lower.tail = F) * 2
que produce 0.04999607
.
Aquí, estamos usando la propiedad conocida de que el 95% del área bajo la distribución normal ocurre en ~ 1.96 desviaciones estándar, por lo que la salida de ~ 0.05 da nuestro valor p. Usé 1000000 ya que cuando N es enorme, la distribución t es casi la misma que la distribución normal. Ejecutar esto me dio comodidad en la solución de @Aaron.