Eu tenho um conjunto de dados de amostra com 31 valores. Fiz um teste t bicaudal usando R para testar se a média verdadeira é 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 Agora estou tentando fazer a mesma coisa manualmente:
t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = dt(t.value, df=length(lengths-1)) O valor t calculado usando isso método é o mesmo que a saída da função R do teste t. O valor p, no entanto, acaba sendo 3,025803e-12.
Alguma ideia do que estou fazendo de errado?
Obrigado!
EDITAR
Aqui está o código R completo , incluindo meu conjunto de dados:
# 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) Resposta
Use pt e torná-lo bicaudal.
> 2*pt(11.244, 30, lower=FALSE) [1] 2.785806e-12 Comentários
- I acho que está faltando um detalhe: Quando usar menor = F? – Veja minha pergunta aqui: stats.stackexchange.com/questions/133091/…
- O valor precisa ser positivo, então, se usar isso com uma variável que pode ser negativa, envolva
abs. - Para um teste de duas caudas, você ' está procurando a probabilidade de que o valor seja menor que -11,244 ou maior que +11,244. inferior = F diz a R para calcular a probabilidade do valor ser maior do que o primeiro parâmetro. Caso contrário, fornece a probabilidade de o valor ser menor que o primeiro parâmetro. Como tal, você também pode fazer 2 * pt (-11,244, 30). Pessoalmente, eu geralmente faço 2 * pt (-abs (q), df = n-1) já que o padrão de R é inferior = T.
Resposta
Eu postei isso como um comentário, mas quando eu quis adicionar um pouco mais na edição, ficou muito longo, então mudei para cá.
Editar : sua estatística de teste e df estão corretos. A outra resposta indica o problema com o cálculo da área da cauda na chamada para pt(), e a duplicação para duas caudas, que resolve sua diferença. No entanto, deixarei minha discussão / comentário anterior porque faz pontos relevantes de forma mais geral sobre os valores-p em caudas extremas:
É possível que você não esteja fazendo nada de errado e ainda obtenha uma diferença, mas se postar um exemplo reproduzível, pode ser possível investigar mais a fundo se você tem algum erro (digamos no df).
Essas coisas são calculadas a partir de aproximações que podem não ser particularmente precisas na cauda extrema .
Se as duas coisas não usam aproximações idênticas, elas podem não concordar intimamente, mas essa falta de acordo não deve importar (para a área da cauda exata tão longe ser um número significativo, as suposições necessárias seriam tem que manter graus de precisão surpreendentes). Você realmente tem normalidade exata, independência exata, variância exatamente constante?
Você não deve necessariamente esperar grande precisão onde os números não significarão nada de qualquer maneira. Até que ponto importa se o valor p aproximado calculado é $ 2 \ vezes 10 ^ {- 12} $ ou $ 3 \ vezes 10 ^ {- 12} $? Nenhum dos números está medindo o valor p real de sua verdadeira situação. Mesmo que um dos números representasse o valor p real da sua verdadeira situação, uma vez que está abaixo de cerca de $ 0,0001 $, por que você se importaria com o que esse valor realmente era?
Resposta
A melhor maneira de calcular manualmente é:
t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = 2*pt(-abs(t.value), df=length(data)-1) Você precisa do abs () porque, de outra forma, você corre o risco de obter valores-p maiores que $ 1 $ (quando a média dos dados é maior do que a média fornecida)!
Resposta
Eu realmente gosto da resposta fornecida por @Aaron, junto com os abs comentários. Acho que uma confirmação útil é executar
pt(1.96, 1000000, lower.tail = F) * 2
que resulta em 0.04999607.
Aqui, estamos usando a propriedade bem conhecida de que 95% da área sob a distribuição normal ocorre em desvios-padrão de ~ 1,96, portanto, a saída de ~ 0,05 dá nosso valor p. Usei 1000000 já que quando N é enorme, a distribuição t é quase a mesma que a distribuição normal. Executá-la me deu conforto na solução de @Aaron.