Nell’ambito dell’inferenza statistica, il test di ipotesi rappresenta uno strumento fondamentale per prendere decisioni informate sui parametri di una popolazione sulla base di dati campionari. Tuttavia, la sola esecuzione di un test di ipotesi, con la sua focalizzazione sul controllo dell’errore di primo tipo (α), non ci fornisce un quadro completo dell’efficacia del nostro procedimento. È qui che entra in gioco il concetto di potenza di un test statistico.
La potenza, definita come la probabilità di rifiutare correttamente un’ipotesi nulla falsa, è cruciale per valutare la sensibilità del nostro test nel rilevare un effetto reale, qualora esso esista. Un test con bassa potenza ha una probabilità elevata di commettere un errore di secondo tipo (β), ovvero di non riuscire a rifiutare un’ipotesi nulla falsa, conducendo a conclusioni errate sulla popolazione in studio.
Consideriamo specificamente il test Z per la media di una popolazione, uno strumento ampiamente utilizzato quando la deviazione standard della popolazione è nota. La procedura tradizionale per calcolare la potenza di questo test coinvolge la determinazione delle regioni di rifiuto definite dal livello di significatività (α) sotto l’ipotesi nulla (H₀), e successivamente il calcolo della probabilità che la statistica test cada in queste regioni sotto una specifica ipotesi alternativa (H₁). Questo processo, sebbene concettualmente chiaro, può risultare laborioso e richiedere l’accesso a tavole statistiche o software specializzato, specialmente quando si esplorano diversi scenari o si desidera visualizzare l’impatto di vari parametri (come la dimensione del campione o la dimensione dell’effetto) sulla potenza.
Proprio per superare queste limitazioni e rendere l’analisi della potenza più accessibile e intuitiva, abbiamo sviluppato un semplice ma efficace script in Python. Questo strumento automatizza i calcoli necessari per determinare la potenza di un test Z, sia per test a una coda che a due code. L’obiettivo è fornire ai ricercatori, agli studenti e ai professionisti uno strumento pratico per:
- Comprendere l’impatto dei parametri: Esplorare come variazioni nella dimensione del campione, nella dimensione dell’effetto (differenza tra la media sotto H₀ e la vera media), nella deviazione standard e nel livello di significatività influenzano la potenza del test.
- Pianificare studi più efficaci: Determinare la dimensione del campione necessaria per raggiungere un livello di potenza desiderato, evitando così studi sotto-dimensionati che potrebbero non riuscire a rilevare effetti importanti.
- Visualizzare la potenza: Ottenere una rappresentazione grafica della distribuzione sotto l’ipotesi nulla e sotto l’ipotesi alternativa, evidenziando le regioni di rifiuto e l’area corrispondente alla potenza del test.
Nel corso di questo articolo, guideremo il lettore attraverso l’utilizzo di questo script, illustrando passo per passo come inserire i parametri necessari e interpretare i risultati ottenuti. Attraverso esempi pratici, mostreremo come questo strumento possa facilitare la comprensione e l’applicazione del concetto di potenza nell’ambito del test di ipotesi Z, contribuendo a una pratica statistica più rigorosa e consapevole.
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import norm
def calcola_potenza(mu0, mu1, sigma, n, alpha, tipo_test):
se = sigma / np.sqrt(n)
if tipo_test == 'una coda':
z_alpha = norm.ppf(1 - alpha)
xbar_crit = mu0 + z_alpha * se
z_beta = (xbar_crit - mu1) / se
potenza = 1 - norm.cdf(z_beta)
elif tipo_test == 'due code':
z_alpha = norm.ppf(1 - alpha / 2)
xbar_crit_sup = mu0 + z_alpha * se
xbar_crit_inf = mu0 - z_alpha * se
z_beta_sup = (xbar_crit_sup - mu1) / se
z_beta_inf = (xbar_crit_inf - mu1) / se
potenza = norm.cdf(z_beta_inf) + (1 - norm.cdf(z_beta_sup))
else:
raise ValueError("Tipo di test non valido. Scegli 'una coda' o 'due code'.")
return {
'xbar_crit': xbar_crit if tipo_test == 'una coda' else (xbar_crit_inf, xbar_crit_sup),
'potenza': potenza
}
def plot_potenza(mu0, mu1, sigma, n, alpha, tipo_test, xbar_crit):
se = sigma / np.sqrt(n)
x = np.linspace(mu0 - 4*se, mu1 + 4*se, 1000)
y0 = norm.pdf(x, mu0, se)
y1 = norm.pdf(x, mu1, se)
plt.figure(figsize=(10, 6))
plt.plot(x, y0, label="Distribuzione sotto H₀", color='blue')
plt.plot(x, y1, label="Distribuzione sotto H₁", color='green')
if tipo_test == 'una coda':
x_alpha = x[x >= xbar_crit]
plt.fill_between(x_alpha, norm.pdf(x_alpha, mu0, se), color='red', alpha=0.3, label="Errore di I tipo (α)")
x_power = x[x >= xbar_crit]
plt.fill_between(x_power, norm.pdf(x_power, mu1, se), color='orange', alpha=0.4, label="Potenza (1 - β)")
plt.axvline(xbar_crit, color='black', linestyle='--', label="Soglia critica")
else:
xbar_crit_inf, xbar_crit_sup = xbar_crit
x_alpha_inf = x[x <= xbar_crit_inf] x_alpha_sup = x[x >= xbar_crit_sup]
plt.fill_between(x_alpha_inf, norm.pdf(x_alpha_inf, mu0, se), color='red', alpha=0.3)
plt.fill_between(x_alpha_sup, norm.pdf(x_alpha_sup, mu0, se), color='red', alpha=0.3, label="Errore di I tipo (α)")
x_power_inf = x[x <= xbar_crit_inf] x_power_sup = x[x >= xbar_crit_sup]
plt.fill_between(x_power_inf, norm.pdf(x_power_inf, mu1, se), color='orange', alpha=0.4)
plt.fill_between(x_power_sup, norm.pdf(x_power_sup, mu1, se), color='orange', alpha=0.4, label="Potenza (1 - β)")
plt.axvline(xbar_crit_inf, color='black', linestyle='--', label="Soglie critiche")
plt.axvline(xbar_crit_sup, color='black', linestyle='--')
plt.title("Grafico della Potenza del Test")
plt.xlabel("Valori della media campionaria X̄")
plt.ylabel("Densità")
plt.grid(True, linestyle='--', alpha=0.5)
plt.legend()
plt.show()
def main():
print("\n📊 Calcolo della Potenza di un Test Z\n")
mu0 = float(input("Inserisci μ₀ (media sotto H₀): ") or 50)
mu1 = float(input("Inserisci μ₁ (vera media): ") or 52)
sigma = float(input("Inserisci σ (deviazione standard): ") or 5)
n = int(input("Inserisci n (dimensione del campione): ") or 30)
alpha = float(input("Inserisci α (livello di significatività): ") or 0.05)
tipo_test = input("Tipo di test ('una coda' o 'due code'): ").strip().lower()
risultati = calcola_potenza(mu0, mu1, sigma, n, alpha, tipo_test)
xbar_crit = risultati['xbar_crit']
potenza = risultati['potenza']
print("\n📌 Risultati:")
if tipo_test == 'una coda':
print(f"Soglia critica X̄_crit = {xbar_crit:.2f}")
else:
print(f"Soglie critiche X̄_crit_inf = {xbar_crit[0]:.2f}, X̄_crit_sup = {xbar_crit[1]:.2f}")
print(f"Potenza del test = {potenza:.4f}")
plot_potenza(mu0, mu1, sigma, n, alpha, tipo_test, xbar_crit)
if __name__ == "__main__":
main()
📊 Calcolo della Potenza di un Test Z
Inserisci μ₀ (media sotto H₀): 3
Inserisci μ₁ (vera media): 2.5
Inserisci σ (deviazione standard): 2
Inserisci n (dimensione del campione): 120
Inserisci α (livello di significatività): 0.05
Tipo di test ('una coda' o 'due code'): due code
📌 Risultati:
Soglie critiche X̄_crit_inf = 2.64, X̄_crit_sup = 3.36
Potenza del test = 0.7819
Spiegazione del Grafico della Potenza del Test
Il grafico illustra visivamente la potenza del test Z .
Ecco una spiegazione dei suoi componenti:
Titolo: “Grafico della Potenza del Test”
Indica lo scopo del grafico, che è quello di visualizzare la potenza del test statistico.
Asse orizzontale (Valori della media campionaria X̄)
Questo asse rappresenta i possibili valori che la media campionaria potrebbe assumere.
Asse verticale (Densità)
Questo asse rappresenta la densità di probabilità per ciascuna distribuzione.
Curva blu (“Distribuzione sotto H₀”)
Questa curva mostra la distribuzione della media campionaria se l’ipotesi nulla (H₀) fosse vera (in questo caso, se la media della popolazione fosse effettivamente 3). È una distribuzione normale con media μ₀ = 3 e deviazione standard pari all’errore standard (σ / √n).
Curva verde (“Distribuzione sotto H₁”)
Questa curva mostra la distribuzione della media campionaria se l’ipotesi alternativa (H₁) fosse vera (in questo caso, se la vera media della popolazione fosse 2.5). È una distribuzione normale con media μ₁ = 2.5 e la stessa deviazione standard (errore standard).
Aree rosse (“Errore di I tipo (α)”)
Queste aree, situate nelle code della distribuzione sotto H₀, rappresentano la regione di rifiuto definita dal livello di significatività α (0.05). Per un test a due code, questa regione è divisa in due parti uguali, una in ciascuna coda. Se la media campionaria osservata cade in queste aree, rifiuteremmo l’ipotesi nulla anche se fosse vera. Le linee verticali tratteggiate nere delimitano queste aree.
Aree arancioni (“Potenza (1 – β)”)
Queste aree rappresentano la potenza del test. Sono le porzioni della distribuzione sotto H₁ che cadono all’interno della regione di rifiuto definita sotto H₀. In altre parole, è la probabilità che la media campionaria osservata porti al rifiuto dell’ipotesi nulla quando l’ipotesi alternativa è vera (cioè, quando la vera media è 2.5). L’area totale in arancione corrisponde al valore della potenza che hai calcolato (0.7819).
Linee verticali tratteggiate nere (“Soglie critiche”)
Queste linee indicano i valori critici della media campionaria (circa 2.64 e 3.36) che separano la regione di accettazione dall’area di rifiuto sotto l’ipotesi nulla. Se la media campionaria cade a sinistra della soglia inferiore o a destra della soglia superiore, rifiuteremo l’ipotesi nulla.
In Sintesi
Il grafico mostra:
- La sovrapposizione tra le due distribuzioni (sotto H₀ e H₁).
- La regione di rifiuto (in rosso) definita in base al livello di significatività α.
- La potenza del test (in arancione), che è la probabilità che la media campionaria, proveniente da una popolazione con media 2.5, cada nella regione di rifiuto definita per l’ipotesi nulla (media 3).
Questo grafico ti aiuta a visualizzare quanto è probabile che il tuo test rilevi correttamente una vera differenza nella media della popolazione (in questo caso, una media di 2.5 invece di 3). Una maggiore area arancione indica una maggiore potenza del test.
Questo codice è utile nella fase di pianificazione di uno studio o di un’analisi statistica in cui si intende eseguire un test di ipotesi sulla media di una popolazione. Permette di calcolare la potenza del test, ovvero la probabilità di rifiutare correttamente un’ipotesi nulla falsa. Questo è fondamentale per capire se un test ha una buona probabilità di rilevare un effetto reale, se presente.
Esempio di utilizzo
Scenario: Resistenza alla trazione dei bulloni
Un’azienda produce bulloni e afferma che la loro resistenza media alla trazione è di 80 MPa (MegaPascal). Un ingegnere di controllo qualità sospetta che la resistenza media sia in realtà inferiore. Decide di eseguire un test di ipotesi con un livello di significatività ([math]\alpha[/math]) del 5%. Prende un campione di 35 bulloni e la deviazione standard della resistenza alla trazione è nota essere di 5 MPa (basata su dati storici). L’ingegnere è particolarmente interessato a sapere la potenza del suo test nel caso in cui la vera resistenza media alla trazione fosse di 78 MPa.
Risoluzione manuale
Ipotesi:
- Ipotesi nulla ([math]H_0[/math]): [math]\mu = 80[/math] MPa (la resistenza media è 80 MPa)
- Ipotesi alternativa ([math]H_1[/math]): [math]\mu < 80[/math] MPa (la resistenza media è inferiore a 80 MPa) – Questo è un test a una coda (sinistra).
Livello di significatività:
[math]\alpha = 0.05[/math]
Statistica test:
Utilizzeremo un test Z perché la deviazione standard della popolazione ([math]\sigma = 5[/math] MPa) è nota. La statistica test è:
[math]Z = \frac{\bar{X} – \mu_0}{\sigma / \sqrt{n}}[/math]
Valore critico:
Per un test a una coda sinistra con [math]\alpha = 0.05[/math], il valore critico di Z è il quantile inferiore 0.05 della distribuzione normale standard, che è circa -1.645.
Regione di rifiuto:
Rifiuteremo [math]H_0[/math] se la statistica test Z è inferiore a -1.645.
Calcolo della potenza:
Vogliamo calcolare la probabilità di rifiutare [math]H_0[/math] quando la vera media ([math]\mu_1[/math]) è 78 MPa. Per fare ciò, dobbiamo trovare la probabilità che la media campionaria ([math]\bar{X}[/math]) sia inferiore al valore critico della media campionaria.
Il valore critico della media campionaria ([math]\bar{x}_{crit}[/math]) sotto [math]H_0[/math] è:
[math]\bar{x}_{crit} = \mu_0 + z_\alpha \cdot \frac{\sigma}{\sqrt{n}} = 80 + (-1.645) \cdot \frac{5}{\sqrt{35}} \approx 80 – 1.645 \cdot \frac{5}{5.916} \approx 80 – 1.645 \cdot 0.845 \approx 80 – 1.39 \approx 78.61[/math] MPa
Quindi, rifiuteremo [math]H_0[/math] se la media campionaria è inferiore a 78.61 MPa.
Ora, calcoliamo la probabilità di ottenere una media campionaria inferiore a 78.61 MPa quando la vera media è 78 MPa. Standardizziamo questo valore usando la media e la deviazione standard sotto [math]H_1[/math]:
[math]Z_\beta = \frac{\bar{x}_{crit} – \mu_1}{\sigma / \sqrt{n}} = \frac{78.61 – 78}{5 / \sqrt{35}} \approx \frac{0.61}{0.845} \approx 0.722[/math]
La potenza è la probabilità che la media campionaria sotto [math]H_1[/math] cada nella regione di rifiuto definita da [math]H_0[/math]. La regione di rifiuto è [math]\bar{X} < 78.61[/math]. Vogliamo [math]P(\bar{X} < 78.61 \mid \mu = 78)[/math].
Il valore Z corrispondente a [math]\bar{x}_{crit} = 78.61[/math] sotto la distribuzione con [math]\mu_1 = 78[/math] è:
[math]Z = \frac{78.61 – 78}{5 / \sqrt{35}} \approx 0.722[/math]
La potenza è [math]P(Z < 0.722)[/math]. Utilizzando la funzione di ripartizione cumulativa della normale standard (norm.cdf(0.722)):
[math]\text{Potenza} \approx 0.7649[/math]
Risoluzione con lo script:
Eseguiamo lo script con i seguenti input:
📊 Calcolo della Potenza di un Test Z
Inserisci μ₀ (media sotto H₀): 80
Inserisci μ₁ (vera media): 78
Inserisci σ (deviazione standard): 5
Inserisci n (dimensione del campione): 35
Inserisci α (livello di significatività): 0.05
Tipo di test ('una coda' o 'due code'): una coda
Output:
📌 Risultati:
Soglia critica X̄_crit = 78.61
Potenza del test = 0.7649Concetti Fondamentali del Test di Ipotesi sulla Media (Test Z)
1. Ipotesi e tipi di errore:
- Ipotesi Nulla (H₀): Un’affermazione sullo stato della popolazione che si cerca di confutare. Nell’esempio, [math]H_0: \mu = 80[/math] MPa.
- Ipotesi Alternativa (H₁): L’affermazione che si accetta se l’ipotesi nulla viene rifiutata. Nell’esempio, [math]H_1: \mu < 80[/math] MPa.
- Errore di Primo Tipo ([math]\alpha[/math]): La probabilità di rifiutare l’ipotesi nulla quando è vera (falso positivo). Solitamente fissato a 0.05 (5%).
- Errore di Secondo Tipo ([math]\beta[/math]): La probabilità di non rifiutare l’ipotesi nulla quando è falsa (falso negativo).
- Potenza ([math]1 – \beta[/math]): La probabilità di rifiutare correttamente l’ipotesi nulla quando è falsa. È la sensibilità del test nel rilevare un effetto reale.
2. Test Z per la Media (Deviazione Standard della Popolazione Nota):
Statistica Test:
Quando la deviazione standard della popolazione ([math]\sigma[/math]) è nota, si utilizza la statistica test Z:
[math]Z = \frac{\bar{X} – \mu_0}{\sigma / \sqrt{n}}[/math]
Dove:
- [math]\bar{X}[/math] è la media campionaria.
- [math]\mu_0[/math] è la media sotto l’ipotesi nulla.
- [math]\sigma[/math] è la deviazione standard della popolazione.
- [math]n[/math] è la dimensione del campione.
- [math]\sigma / \sqrt{n}[/math] è l’errore standard della media (SE).
Regione di Rifiuto:
La regione di rifiuto è l’insieme dei valori della statistica test che portano al rifiuto di [math]H_0[/math]. Questa regione è determinata dal livello di significatività ([math]\alpha[/math]) e dal tipo di test (una coda o due code).
- Test a una coda (sinistra): Rifiutiamo [math]H_0[/math] se [math]Z < -z_\alpha[/math].
- Test a una coda (destra): Rifiutiamo [math]H_0[/math] se [math]Z > z_\alpha[/math].
- Test a due code: Rifiutiamo [math]H_0[/math] se [math]|Z| > z_{\alpha/2}[/math].
Dove [math]z_\alpha[/math] è il valore critico della distribuzione normale standard corrispondente a un livello di significatività [math]\alpha[/math].
Spiegazione:
-
Quando è utile questo codice? Questo codice è utile quando si vuole valutare la probabilità che un test di ipotesi rilevi un effetto di una certa dimensione. Nell’esempio, l’ingegnere voleva sapere quanto fosse probabile rilevare una diminuzione della resistenza media dei bulloni a 78 MPa, dato un campione di 35 bulloni e un livello di significatività del 5%. Una potenza di circa 0.76 (76%) indica che il test ha una probabilità ragionevole di rilevare una differenza di 2 MPa nella resistenza media.
-
Esempio completo: L’esempio dei bulloni illustra come questo codice possa essere utilizzato in un contesto di controllo qualità. Prima di eseguire il test, l’ingegnere può utilizzare lo script per determinare se la dimensione del campione è adeguata per avere una buona probabilità di rilevare una differenza significativa nella resistenza media. Se la potenza fosse troppo bassa, potrebbe decidere di aumentare la dimensione del campione.
-
Confronto tra la soluzione manuale e lo script: La soluzione manuale, sebbene concettualmente corretta, richiede l’uso di tabelle statistiche o un calcolatore statistico per trovare i valori dei quantili della distribuzione normale. Lo script automatizza questo processo, fornendo rapidamente i risultati della potenza del test e la soglia critica. Questo rende il processo molto più efficiente e meno soggetto a errori di calcolo manuale.
In sintesi, questo codice è uno strumento prezioso per la power analysis, che è una parte essenziale della progettazione di esperimenti e studi statistici. Aiuta a prendere decisioni informate sulla dimensione del campione, il livello di significatività e la probabilità di ottenere risultati significativi.
