Labolatorium 2 - test t-Studenta

Karolina Marek

2024-08-06

Test t-Studenta - wariant 1

Wariant 1 będzie dotyczył testu t-Studenta dla testowania jednej próby. Jest on używany w sytuacji gdy chcemy porównać średnią próby .

Na przykład:

  • Testowanie skuteczności nowego leku: Załóżmy, że średni czas reakcji na pewne leczenie jest znany (np. 30 minut). Chcemy sprawdzić, czy nowy lek skraca ten czas. Bierzemy próbkę pacjentów, podajemy im nowy lek i mierzymy czas reakcji, a następnie testujemy hipotezę, czy średnia z próby jest mniejsza niż 30 minut.

  • Porównanie średniej w ocenie produktów: Firma przeprowadzając na przykład ankietę zna średnią ocenę swojego starego produktu (np. 4 na 5). Następnie, chce sprawdzić, czy nowy produkt ma wyższą ocenę. W tym celu mogą przeprowadzić test dla średniej oceny nowego produktu.

Mała liczebność próby Test t-Studenta dla jednej próby jest szczególnie przydatny, gdy mamy małą liczbę obserwacji (zwykle \(n<30\)).

Nieznana wariancja Test t-Studenta jest stosowany, gdy wariancja populacji (\(\sigma^2\)) nie jest znana, co jest bardzo częstą sytuacją w badaniach empirycznych. W takich przypadkach używamy estymatora wariancji z próby (\(S^2\)) zamiast rzeczywistej wariancji populacji.

Zadanie 1

Wygenerujmy \(n=20\) obserwacji z rozkładu \(N(d,5)\), gdzie \(d=1\). Przetestujemy \(H_0: d\leq 0\) przeciwko \(H_1: d>0\). Wykonaj test t-studenta dwoma sposobami: pisząc cały test samodzielnie i używając gotowej funkcji.

Kroki:

  1. Określić jaki mamy test.

  2. Wycentrowanie obserwacji przy założeniu hipotezy zerowej.

Centrujemy \(Y=(Y_1, \dots, Y_n)= (X_1-d, \dots, X_n-d)\).

Teraz \(\sum_{i=1}^nY_i\sim N(0,n\sigma^2)\). Zatem \(\frac{1}{sqrt{n}}\sum_{i=1}^nY\sim N(0,\sigma^2)\).

  1. Obliczyć statystykę testową i określić jej rozkład.

Statystyka testowa \[U(Y)=\frac{\frac{1}{\sqrt{n}}\sum_{i=1}^nY_i}{\sqrt{\frac{1}{n-1}\sum_{i=1}^n(Y_i-\bar{Y})^2}}\sim t_{n-1}\]

  1. Wyznaczenie obszaru krytycznego.

Na wykładzie zostało pokazane, że poniższy test jest jednostajnie najmocniejszy wśród testów nieobciążonych.

\[\phi(X)=1,\quad U(Y)> c_{\alpha},\] gdzie \(c_{\alpha}=F^{-1}_{t_{n-1}}(1-\alpha)\).

  1. P-wartość

Skoro \(U(Y)=F^{-1}(1-p(Y))\), to p-wartość wyraża się poniższym wzorem \[p(Y)=1-F_{t_{n-1}}(U(Y)).\]

set.seed(5)
n<-20
m<-0
X<-rnorm(n,1,5)
Y=X-m
test_statistic<-sqrt(1/n)*sum(Y)/sqrt(sum((Y-mean(Y))^2)/(n-1)) 
p<-1-pt(test_statistic,n-1) 
test_statistic 
## [1] -0.3876304
p
## [1] 0.6487011

Duża p-wartość świadczy o tym, że nie mamy podstaw do odrzucenia hipotezy zerowej

W języku R jest również dostępna gotowa funkcja t.test().

t.test(X,mu=0,alternative="greater") 
## 
##  One Sample t-test
## 
## data:  X
## t = -0.38763, df = 19, p-value = 0.6487
## alternative hypothesis: true mean is greater than 0
## 95 percent confidence interval:
##  -2.19954      Inf
## sample estimates:
##  mean of x 
## -0.4027889

Zadanie 1 cd. Sporządźmy teraz wykres mocy w zależności od \(d\in[0,10]\) oraz \(s\in [1, 20]\).

W tym celu utworzymy siatkę \((d,s)\), gdzie \(d \in [0,10]\) oraz \(s\in[1,20]\). Następnie obliczymy moc empiryczną dla każdego punktu siatki i utworzymy wykres ciepła.

d <- seq(0,10, by=0.5)
s <- seq(1,20)
Moc<- matrix(0,nrow=20,ncol=21)   
mc = 1000
alpha=0.05
for (j in 1:length(d)){
  for (k in 1:length(s)) {
    M <- 0 
    for (i in 1:mc){
      X <- rnorm(n,j,k) 
      test_statistic <- sqrt(1/n) * (sum(X - 0)) / sqrt(sum((X - mean(X))^2) / (n-1)) 
      P=1-pt(test_statistic,n-1)    
      if(P< alpha){
        M <- 1+M
        }    
      }
    Moc[k,j]<-M/mc
    }
  }
heatmap(Moc,Rowv=NA,Colv=NA,scale="none", xlab="d", ylab="s")

Zauważmy, że im większa wariancja tym wolniej rośnie moc.

Test t-studenta - wariant 2

Wariant ten jest używany w kontekście testowania hipotezy dotyczącej średniej populacji \(m\), szczególnie gdy chcemy przetestować, czy średnia jest równa pewnej wartości \(m_0\). Czyli \(H_0: m=m_0\) przy alternatywie \(H_1: m\neq m_0\). Tutaj również przyjmujemy założenie, że dane pochodzą z rozkładu normalnego o nieznanych parametrach.

Zadanie 2

\(30\) samochodów użyto do zbadania dwóch rodzajów paliwa pod kątem spalania. Samochody zatankowały \(10\) litrów paliwa typu \(1\) i jeździły po torze do momentu, aż paliwo się skończyło. Zanotowano liczę przejechanych kilometrów. Później te same samochody zatankowały paliwo typu \(2\) i powtórzono eksperyment. Na podstawie zebranych wyników przetestuj czy te dwa typy paliwa różnią się między sobą. Dane znajdują się w pliku Baza 2.1.

Uwaga: Aby porównać średnie dla dwóch różnych typów paliwa musimy założyć, że mają one jednorodna wariancję.

Niech \(X=Y_i-Z_i\), gdzie \(Y_i\) będą dotyczyć danych dotyczących spalania paliwa typu \(1\), a \(Z_i\) paliwa typu \(2\).

Wtedy testujemy: \(H_0: m=0\), przeciwko \(H_1: m\neq 0\), gdzie \(m\) jest średnią dla próby \(X\).

Użyjemy tej samej statystyki testowej co w wariancie 1: \[T(X)=\frac{\sqrt{\frac{1}{n}}\sum_{j=1}^n(X_j-m_0)}{\sqrt{\frac{1}{n-1}\sum_{j=1}^n\left(X_j-\bar{X}\right)^2}}.\] W tym przypadku mamy test obustronny postaci: \[\phi(X)=1\quad, |T(X)|>c_{\alpha},\] gdzie \(c_{\alpha}=F^{-1}_{t_{n-1}}(1-\frac{1}{2})\).

Test ten jest jednostajnie najmocniejszy wśród testów nieobciążonych.

library(readxl)
spalanie <- read_excel("Baza 2.1.xlsx")
var(spalanie$"Typ 1")
## [1] 286.5844
var(spalanie$"Typ 2")
## [1] 387.0433
X<-spalanie$"Typ 1" - spalanie$"Typ 2"
n<-length(X)
test_statistic<-sqrt(1/n)*(sum(X-0))/sqrt(sum((X-mean(X))^2)/(n-1))
P<-2*(1-pt(abs(test_statistic),n-1)) 
mean(X) 
## [1] -2.80619
test_statistic
## [1] -0.6203171
P
## [1] 0.5398936

Zatem nie mamy podstaw do odrzucenia. Okazuje się, że 2 metry różnicy w praktyce nie jest statystycznie istotne.

Możemy również użyć gotowej funkcji t.test().

t.test(X,mu= 0, alternative = "two.side")
## 
##  One Sample t-test
## 
## data:  X
## t = -0.62032, df = 29, p-value = 0.5399
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -12.05840   6.44602
## sample estimates:
## mean of x 
##  -2.80619

Zadanie

Wygenerujmy \(n=20\) obserwacji z rozkładu \(N(d,5)\), gdzie \(d=1\). Oblicz empiryczną moc dla testu \(H_0: d\leq 0\) przeciwko \(H_1: d>0\). Wygeneruj histogram p-wartości.
Wygenerujmy \(n=20\) obserwacji z rozkładu \(N(d,5)\), gdzie \(d\in[-2,2]\). Narysuj krzywą mocy dla \(H_0: d= 0\) przeciwko \(H_1: d\neq0\).

Test t-Studenta - wariant 3

W tym wariancie rozważamy problem testowania hipotezy dotyczącej różnicy średnich dwóch niezależnych próbek, pochodzących z rozkładów normalnych i posiadających jednorodną wariancję.

\[X_1, \dots, X_m\sim N(m_x, \sigma^2)\]

\[Y_1, \dots, Y_n\sim N(m_y, \sigma^2)\] Testujemy: \[H_0: m_y\leq m_x\qquad H_1: m_y>m_x.\]

Ponieważ próby są różnej liczności, to nie możemy postąpić tak jak w poprzednim zadaniu. Dlatego tutaj będziemy porównywać różnicę średnich. Na początku wyznaczmy rozkład \(\bar{X}-\bar{Y}\).

\[E[\bar{X}-\bar{Y}]=E[\bar{X}]-E[\bar{Y}]=m_x-m_y\] \[Var[\bar{X}-\bar{Y}]=Var[\bar{X}]+Var[\bar{Y}]=\frac{\sigma^2}{n}+\frac{\sigma^2}{m}=\frac{(m+n)\sigma^2}{mn}\] Zatem, \[(\bar{X}-\bar{Y})\sqrt{\frac{mn}{(m+n)}}\sim N(0,\sigma^2 )\]

Statystyka testowa jest postaci \[U(X,Y)=\frac{\sqrt{\frac{mn}{(m+n)}}\left(\bar{Y}-\bar{X}\right)}{\sqrt{\frac{1}{m+n-2}\left(\sum_{i=1}^m(X_i-\bar{X})^2+\sum_{j=1}^n(Y_j-\bar{Y})^2\right)}}\sim t_{m+n-2}.\] Wobec tego, test postaci \[\phi(X,Y)=1, \quad U(X,Y)>c_{\alpha},\] gdzie \(c_{\alpha}=F^{-1}_{t_{m+n-2}}(1-\alpha)\) jest testem jednostajnie najmocniejszym wśród testów nieobciążonych.

P- wartość: \[p(X)=1-F_{t_{n+m-2}}\left(U(X,Y)\right).\]

Zadanie 3

Aby zbadać skuteczność pewnego leku na cukrzycę przeprowadzono badanie na \(40\) zdrowych i \(34\) chorych na cukrzycę osobach. Osoby te przez tydzień jadły to samo. Osoby chore przyjmowały lek. Po tygodniu zbadano wszystkim poziom cukru we krwi na czczo. Wyniki zebrano w pliku Baza 2.2. Przetestować hipotezę o skuteczności leku przeciwko alternatywie o nieskuteczności leku.

Baza<- read_excel("Baza 2.2.xlsx")
summary(Baza) # w zbiorze danych znajdują się wartości NA
##      Zdrowi           Chorzy      
##  Min.   : 78.00   Min.   : 71.00  
##  1st Qu.: 83.75   1st Qu.: 95.25  
##  Median : 94.50   Median :102.00  
##  Mean   : 94.22   Mean   :101.44  
##  3rd Qu.:100.50   3rd Qu.:109.75  
##  Max.   :120.00   Max.   :120.00  
##                   NA's   :6
X<-na.exclude(Baza$"Zdrowi")
Y<-na.exclude(Baza$"Chorzy") 
nx<-length(X) 
ny<-length(Y)
test_statistic<-sqrt((nx*ny)/(nx+ny))*((mean(Y)-mean(X))/sqrt((sum((X-mean(X))^2)+sum((Y-mean(Y))^2))/(nx+ny-2))) 
p<-1-pt(test_statistic,nx+ny-2) 
mean(Y)-mean(X) 
## [1] 7.216176
p 
## [1] 0.002969341

Odrzucamy hipotezę zerową na rzecz alternatywy. Pomimo zastosowania leku, poziom cukry był istotnie różny.

Możemy również użyć gotowej funkcji t.test().

t.test(Y,X, var.equal = TRUE,alternative ="greater")
## 
##  Two Sample t-test
## 
## data:  Y and X
## t = 2.8353, df = 72, p-value = 0.002969
## alternative hypothesis: true difference in means is greater than 0
## 95 percent confidence interval:
##  2.975291      Inf
## sample estimates:
## mean of x mean of y 
##  101.4412   94.2250

Test t-studenta - wariant 4

Dla dwóch prób o rozkładzie normalnym i jednorodnej wariancji (tak, jak w wariancie 3), testujemy: \(H_0: m_y=m_x\), przeciwko \(m_x\neq m_y\).

Używamy tej samej statystyki co dla wariantu 3, tj \[U(X,Y)=\frac{\sqrt{\frac{mn}{(m+n)}}\left(\bar{Y}-\bar{X}\right)}{\sqrt{\frac{1}{m+n-2}\left(\sum_{i=1}^m(X_i-\bar{X})^2+\sum_{j=1}^n(Y_j-\bar{Y})^2\right)}}\sim t_{m+n-2}.\]

Test \(\phi(X,Y)=1\) dla \(|U(X,Y)|>c_{\alpha}\), gdzie \(c_{\alpha}=F^{-1}_{t_{m+n_2}}(1-\alpha/2)\) jest testem jednostajnie najmocniejszym wśród testów nieobciążonych.

P-wartość: \(p(X,Y)=2\left(1-F_{t_{m+n_2}}\left(U(X,Y)\right)\right).\)

Zadanie 4

Przetestować hipotezę o braku różnicy między grupami zdrowych i chorych przeciwko alternatywie, że jest istotna różnica dla zbioru danych Baza 2.2.xlsx.

Baza<- read_excel("Baza 2.2.xlsx")
X<-na.exclude(Baza$"Zdrowi")
Y<-na.exclude(Baza$"Chorzy") 
nx<-length(X) 
ny<-length(Y) 
test_statistic<-sqrt((nx*ny)/(nx+ny))*((mean(Y)-mean(X))/sqrt((sum((X-mean(X))^2)+sum((Y-mean(Y))^2))/(nx+ny-2)))
p<-2*(1-pt(abs(test_statistic),nx+ny-2))
abs(mean(Y)-mean(X)) 
## [1] 7.216176
p
## [1] 0.005938681

Odrzucamy hipotezę zerową na rzecz alternatywy.

Używając gotowej t.test().

t.test(Y,X,var.equal=TRUE,alternative="two.sided")
## 
##  Two Sample t-test
## 
## data:  Y and X
## t = 2.8353, df = 72, p-value = 0.005939
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##   2.142611 12.289742
## sample estimates:
## mean of x mean of y 
##  101.4412   94.2250

Różnica pomiędzy wariantem 3 a wariantem 4

Zauważmy, że w zadaniu 3 sprawdzamy, czy lek nie powoduje wzrostu poziomu cukru, co jest testem jednostronnym. Natomiast, zadanie 4 koncentruje się na sprawdzeniu, czy istnieje jakakolwiek różnica w poziomie cukru między zdrowymi a chorymi, co jest testem dwustronnym.

Warto również pamiętać, o odporności testu t Studenta na naruszenie założenia normalności.

Zadanie

Aby zbadać skuteczność pewnego leku na cukrzycę przeprowadzono badanie na \(40\) zdrowych i \(34\) chorych na cukrzycę osobach. Osoby te przez tydzień jadły to samo. Osoby chore przyjmowały lek. Po tygodniu zbadano wszystkim poziom cukru we krwi na czczo. Wyniki zebrano w pliku Baza 2.2. Przetestować hipotezę o skuteczności leku przeciwko alternatywie o nieskuteczności leku.

Definiujemy zmienne \(X\) i \(Y\)na odwrót i używamy tej samej statystyki testowej. Wtedy, \(H_0: m_x\leq m_y\) przeciw \(H_1: m_x>m_y\). Zbuduj test lewostronny. Wyznacz p wartość i empiryczną moc testu.

Baza<- read_excel("Baza 2.2.xlsx")
summary(Baza) # w zbiorze danych znajdują się wartości NA
##      Zdrowi           Chorzy      
##  Min.   : 78.00   Min.   : 71.00  
##  1st Qu.: 83.75   1st Qu.: 95.25  
##  Median : 94.50   Median :102.00  
##  Mean   : 94.22   Mean   :101.44  
##  3rd Qu.:100.50   3rd Qu.:109.75  
##  Max.   :120.00   Max.   :120.00  
##                   NA's   :6
X<-na.exclude(Baza$"Chorzy")
Y<-na.exclude(Baza$"Zdrowi") 
nx<-length(X) 
ny<-length(Y)
test_statistic<-sqrt((nx*ny)/(nx+ny))*((mean(Y)-mean(X))/sqrt((sum((X-mean(X))^2)+sum((Y-mean(Y))^2))/(nx+ny-2)))