Ho un set di dati di esempio con 31 valori. Ho eseguito un test t a due code utilizzando R per verificare se la media vera è uguale a 10:
t.test(x=data, mu=10, conf.level=0.95) Risultato:
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 Ora sto cercando di fare la stessa cosa manualmente:
t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = dt(t.value, df=length(lengths-1)) Il valore t calcolato utilizzando questo è lo stesso delloutput della funzione t-test R. Il valore p, tuttavia, risulta essere 3.025803e-12.
Qualche idea su cosa sto sbagliando?
Grazie!
EDIT
Ecco il codice R completo , compreso il mio set di dati:
# 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) Answer
Usa pt e renderlo a due code.
> 2*pt(11.244, 30, lower=FALSE) [1] 2.785806e-12 Commenti
- I penso che manchi un dettaglio: Quando usare lower = F? – Consulta la mia domanda qui: stats.stackexchange.com/questions/133091/…
- Il valore deve essere positivo, quindi se lo utilizzi con una variabile che potrebbe essere negativa, includi
abs. - Per un test a due code, devi ' sta cercando la probabilità che il valore sia inferiore a -11,244 o superiore a +11,244. inferiore = F dice a R di calcolare la probabilità che il valore sia maggiore del primo parametro. Altrimenti, ti dà la probabilità che il valore sia inferiore al primo parametro. In quanto tale, potresti anche fare 2 * pt (-11.244, 30). Personalmente, di solito eseguo 2 * pt (-abs (q), df = n-1) poiché R per impostazione predefinita è inferiore = T.
Risposta
Lho pubblicato come commento ma quando volevo aggiungerne un po di più in modifica, è diventato troppo lungo, quindi lho spostato qui.
Modifica : la statistica del test e df sono corretti. Laltra risposta rileva il problema con il calcolo dellarea della coda nella chiamata a pt(), e il raddoppio per due code, che risolve la tua differenza. Tuttavia lascio la mia discussione / commento precedente perché rende i punti rilevanti più in generale sui valori p nelle code estreme:
È possibile che tu non stia facendo nulla di sbagliato e che continui a ottenere una differenza, ma se pubblichi un esempio riproducibile potrebbe essere possibile indagare ulteriormente se hai qualche errore (ad esempio nel df).
Queste cose sono calcolate da approssimazioni che potrebbero non essere particolarmente accurate nella coda estrema .
Se le due cose non utilizzano approssimazioni identiche, potrebbero non essere strettamente concordi, ma tale mancanza di accordo non dovrebbe importare (affinché larea di coda esatta fino a quel punto sia un numero significativo, i presupposti richiesti sarebbero devono mantenere livelli di accuratezza sbalorditivi). Hai davvero la normalità esatta, lindipendenza esatta, la varianza esattamente costante?
Non dovresti necessariamente aspettarti una grande precisione là dove i numeri non significano nulla comunque. In che misura è importante se il valore p approssimativo calcolato è $ 2 \ times 10 ^ {- 12} $ o $ 3 \ times 10 ^ {- 12} $? Nessuno dei due numeri misura il valore p effettivo della tua situazione reale. Anche se uno dei numeri rappresentasse il valore p reale della tua situazione reale, una volta che è inferiore a $ 0,0001 $, perché ti interesserebbe quale fosse effettivamente quel valore?
Risposta
Il modo migliore per calcolarla manualmente è:
t.value = (mean(data) - 10) / (sd(data) / sqrt(length(data))) p.value = 2*pt(-abs(t.value), df=length(data)-1) Ti serve il funzione abs () perché altrimenti corri il rischio di ottenere valori p maggiori di $ 1 $ (quando la media dei dati è maggiore della media data)!
Risposta
Mi piace molto la risposta fornita da @Aaron, insieme ai abs commenti. Trovo che una comoda conferma sia lesecuzione
pt(1.96, 1000000, lower.tail = F) * 2
che restituisce 0.04999607.
Qui stiamo usando la nota proprietà che il 95% dellarea sotto la distribuzione normale si verifica a ~ 1,96 deviazioni standard, quindi loutput di ~ 0,05 fornisce il nostro valore p. Ho usato 1000000 poiché quando N è enorme, la distribuzione t è quasi la stessa della distribuzione normale. Lesecuzione di questo mi ha dato conforto nella soluzione di @Aaron.