Rangi
Rozważmy nieparametryczne problemy testowania, które będą niezmiennicze względem grupy ciągłych i ściele rosnących (niekoniecznie liniowych) przekształceń.
Przypomnijmy z wykładu.
Niech \((X_1, \dots, X_n)\) będzie próbą prostą o ciągłej dystrybuancie \(F\) oraz niech \[G=\{g: \mathbb{R}^n\rightarrow \mathbb{R}^n: g(\mathbf{X})=\left(h(x_1), \dots, h(x_n)\right), \text{ h - ciągła i ściśle rosnąca funckja}\}\] będzie grupą przekształceń.
Rangą i-tej obserwacji wektora \((x_1, \dots, x_n)\) nazywamy liczbę współrzędnych tego wektora mniejszych lub równych \(x_i\), tzn: \[r_i(x_1, \dots, x_n)=\sum_{j=1}^n\mathbb{I}_{(-\infty, x_i]}(x_j).\]
Przykład: \(Rang(9,6,1,4,5)=\)
Odwzorowanie Rang jest maksymalnym niezmiennikiem względem grupy przekształceń G (pokazane na wykładzie).
Zadanie 1
Wygenerować \(20\) obserwacji \(X_1, \dots, X_{20}\) z rozkładu \(N(0,1)\). Narysować punkty \((X_i, R_i)\).
set.seed(321)
n <- 20
X <- rnorm(n)
R <- numeric(n)
for (i in 1:n) {
R[i] <- sum(X <= X[i])
}
plot(X,R)
Zadanie 2
Wygenerować \(20\) obserwacji \(X_1, \dots, X_{20}\) z rozkładu Poissona \(P(5)\). Narysować punkty \((X_i, R_i)\).
## [1] 9 9 3 3 4 4 5 4 5 7 5 4 7 2 5 3 6 4 4 6
Zadanie 2 - cd
Narysować wykres jeszcze raz uwzględniając poprawki na powtarzające się rangi.
Dodajemy do rang pewną losowość aby uniknąć powtórzeń.
set.seed(321)
n <- 20
X <- rpois(n,5)
R <- numeric(n)
for (i in 1:n) {
R[i] <- sum(X <= X[i])
}
U <- runif(n, 0, 1)
C <- R + U
RC <- numeric(n)
for (i in 1:n) {
RC[i] <- sum(C <= C[i])
}
plot(X,RC)
points(X,R,pch=3,col="red")
Możemy również użyć gotowej funkcji rank()
.
## [1] 20 20 4 4 10 10 14 10 14 18 14 10 18 1 14 4 16 10 10 16
## [1] 20 19 2 4 9 8 11 10 12 18 13 6 17 1 14 3 15 5 7 16
Test Manna- Whitneya
Przypomnimy, test t-Studenta (wariant 3 i 4) mieliśmy: \(X_1, \dots, X_m\sim N(m_x, \sigma^2)\) oraz \((Y_1, \dots, Y_n)\sim(N(m_y, \sigma^2))\) i testowaliśmy:
(wariant 3) \(H_0: m_y\leq m_x\), przeciwko \(H_1: m_y> m_x\),
(wariant 4) \(H_0: m_y = m_x\), przeciwko \(H_1: m_y\neq m_x\).
Co gdy mamy sytuację, że nie możemy założyć normalności rozkładów, ani jednorodności wariancji?
Zatem mamy:
\((X_1, \dots, X_m)\), \((Y_1, \dots, Y_n)\) niezależne próby proste o dystrybuantach odpowiednio \(F_x\), \(F_Y\) takich, że \[\exists_{\Delta: \mathbb{R}\rightarrow \mathbb{R}}F_Y(t)=F_X\left(t\right)-\Delta(t),\] gdzie \(\Delta\) jest stałego znaku oraz \(F_X\) jest ciągła.
Zauważmy, że nieznanymi parametrami są \(\theta=(\Delta,F_X)\) - zbiór takich funkcji nie należy do skończenie wymiarowej przestrzeni Euklidesowej - nieparametryczny problem testowania.
Przypomnijmy, że odwzorowanie \[Rang(x_1, \dots, x_m, y_1, \dots, y_n)=(r_1, \dots, r_{m+n})\] jest maksymalnym niezmiennikiem względem grupy przekształceń \[G=\left\{g: \mathbb{R}^{m+n}\rightarrow \mathbb{R}^{m+n}: g\left(\mathbf{x,y}\right)=\left(h(x_1), \dots, h(x_m), h(y_1), \dots, h(y_n)\right), \text{ h- ciągła i ściśle rosnąca}\right\}.\] Zatem, test będzie oparty na wektorze rang \((R_1, \dots, R_{m+n})\) odpowiadającym łącznej próbie \((X_1, \dots, X_m, Y_1, \dots, Y_n)\).
Idea testu Manna-Whitney opiera się na porównaniu sum rang dla \(X\) i \(Y\):
\(\bar{R}_X=\frac{1}{m}\sum_{i=1}^m R_i\) - rangi odpowiadające \(X\),
\(\bar{R}_Y=\frac{1}{n}\sum_{i=m+1}^{m+n} R_i\) - rangi odpowiadające \(Y\).
Zauważmy, że jeżeli:
\(\Delta\leq 0\), wtedy \(\mathbf{E}(\bar{R}_Y)\leq \mathbf{E}(\bar{R}_X)\),
\(\Delta> 0\), wtedy \(\mathbf{E}(\bar{R}_Y) >\mathbf{E}(\bar{R}_X)\).
Wobec tego,
(dla wariantu 3) tj. \(H_0: \Delta\leq 0\), przeciwko \(H_1: \Delta>0\), będziemy odrzucać hipotezę zerową, gdy \(\bar{R}_Y-\bar{R}_X\) będzie duże,
(dla wariantu 4) tj. \(H_0: \Delta = 0\), przeciwko \(H_1: \Delta\neq0\), będziemy odrzucać hipotezę zerową, gdy \(\left\|\bar{R}_Y-\bar{R}_X\right\|\) będzie duże,
Wykonując standaryzację, otrzymujemy statystykę testową Manna-Whitney: \[W=\frac{\bar{R}_Y-\bar{R}_X}{\sqrt{V_{\Delta=0}\left(\bar{R}_Y-\bar{R}_X\right)}}=\frac{\bar{R}_Y-\frac{(m+n+1)}{2}}{\sqrt{\frac{m(m+n+1)}{12m}}}.\] Gdy \(F_X=F_Y\) zachodzi \(W\rightarrow N(0,1)\) dla \(n\rightarrow \infty\) oraz \(\min\{m,n\}\rightarrow \infty\).
Ostatecznie, dla $H_0: $, przeciwko alternatywie \(\Delta>0\) mamy \(\phi(W)=1\) dla \(W>c_{\alpha}\), gdzie \(c_{\alpha}=F^{-1}_{N(0,1)}(1-\alpha)\) i jest testem asymptotycznie zgodnym.
Asymptotyczna p-wartość: \(p(W)=1-F_{N(0,1)}(W)\).
Natomiast dla $H_0: = 0 $, przeciwko alternatywie \(\Delta \neq 0\) mamy \(\phi(W)=1\) dla \(|W|>c_{\alpha}\), gdzie \(c_{\alpha}=F^{-1}_{N(0,1)}(1-\frac{\alpha}{2})\) i jest testem asymptotycznie zgodnym.
Asymptotyczna p-wartość: \(p(W)=2\left(1-F_{N(0,1)}(W)\right)\).
Zadanie 3
Pacjenci u których zdiagnozowano pewną chorobę są leczeni jedną z dwóch terapii. Zmierzono czas życia od momentu rozpoczęcia terapii dla \(58\) pacjentów leczonych terapią \(1\) oraz dla \(52\) pacjentów leczonych terapią \(2\). Czy terapie się różnią, a jeśli tak, to która jest lepsza? (Dane są w pliku Baza 6.1)
library(readxl)
Baza<- read_excel("Baza 6.1.xlsx")
Z<-Baza$"Czas"
G<- Baza$"Leczenie"
X <- Z[G == 1]
Y <- Z[G == 2]
m <- length(X)
n <- length(Y)
# Zwrócić uwagę na przesunięcie dystrybuant
par(mfrow=c(2,2))
hist(X,prob=TRUE)
hist(Y,prob=TRUE)
# Wykresy dystrybuant empirycznych
plot(ecdf(X))
lines(ecdf(Y),col="blue")
par(mfrow=c(1,1))
R<-rank(Z)
# Wybór rang dla grupy Y
RY <- R[G == 2]
W<-(mean(RY)-(m+n+1)/2)/sqrt(m*(m+n+1)/(12*n))
p<-2*(1-pnorm(abs(W),0,1))
p
## [1] 0.02400052
## [1] 0.9879997
Zatem, terapie różnią się miedzy sobą rozkładem. Możemy również stwierdzić, że $ $, wobec tego terapia druga jest lepsza (krótsza) od terapii pierwszej.
Możemy również użyć gotowej funkcji wilcox.test
.
##
## Wilcoxon rank sum test with continuity correction
##
## data: X and Y
## W = 1885, p-value = 0.02419
## alternative hypothesis: true location shift is not equal to 0
##
## Wilcoxon rank sum test with continuity correction
##
## data: X and Y
## W = 1885, p-value = 0.9881
## alternative hypothesis: true location shift is less than 0
Obliczanie p-wartości testu Manna-Whitneya metodą MC
Przypomnijmy, że gdy \(\Delta=0\), to wektor rang ma rozkład jednostajny na zbiorze permutacji \((1, \dots, m+n)\). Zatem równoważnie możemy wygenerować zmienne \(\mathbf{U}=U_1, \dots, U_{m+n}\sim(U(0,1))\) i obliczyć \(Rang(U)\).
Kroki:
Generujemy \(mc\) (np. \(mc=1000\)) prób prostych \(U_1, \dots, U_{m+n}\sim U[0,1]\).
Dla każdej z prób z kroku \(1\) obliczamy rangi \(R_1, \dots, R_{m+n}\).
Dla każdego wektora rang z korku \(2\) obliczamy statystykę Manna-Whitney \(W_1, \dots, W_{mc}\).
Dla \(H_0: \Delta\leq 0\), przeciwko \(H_1: \Delta>0\) mamy \[\hat{p}=\frac{1}{mc}\sum_{i=1}^{mc}\mathbb{I}(W_i\geq W).\] Dla \(H_0: \Delta = 0\), przeciwko \(H_1: \Delta\neq 0\) mamy \[\hat{p}=\frac{1}{mc}\sum_{i=1}^{mc}\mathbb{I}\left(|W_i|\geq |W|\right).\] Statystyka \(W\) jest obliczona dla oryginalnych danych.
Zadanie 3 - ciąg dalszy
Obliczyć dokładną p wartość metodą MC.
set.seed(123)
mc <- 100000
WMC <- numeric(mc)
for (i in 1:mc){
UMC <- runif(m + n)
RMC <- rank(UMC)
RYMC <- RMC[(m + 1):(m + n)]
WMC[i] <- (mean(RYMC) - (m + n + 1) / 2) / sqrt(m * (m + n + 1) / (12 * n))
}
pMC <- mean(abs(WMC) >= abs(W))
pMC
## [1] 0.0229
## [1] 0.98855
Zadanie
Niech \((X_1, \dots, X_{20})\) - i.i.d. z rozkładu \(N(m_x=0, 1)\) oraz \((Y_1, \dots, Y_{30})\) - i.i.d z rozkładu \(N(m_y, 1)\).
Sporządź wykres mocy empirycznej dla \(m_y\in [-1,2]\) dla testu t-Studenta. (wariant 4: \(H_0: m_X=m_y\), przeciwko \(H_1: m_X\neq m_y\)).
Dla tych samych danych obliczyć p-wartości testu t-Studenta i testu Manna-Whitney. Narysuj wykres różnić p-wartości.
Dodać do wykresu z podpunktu a) moc empiryczną testu Manna-Whitney.