Ustal, czy proces rozproszenia ciężkiego ogona znacznie się poprawił

Obserwuję czasy przetwarzania procesu przed i po zmianie, aby dowiedzieć się, czy proces poprawił się przez zmianę. Proces poprawił się, jeśli czas przetwarzania został skrócony. Rozkład czasu przetwarzania jest gruby, dlatego porównanie na podstawie średniej nie jest rozsądne. Zamiast tego chciałbym wiedzieć, czy prawdopodobieństwo zaobserwowania krótszego czasu przetwarzania po zmianie znacznie przekracza 50%.

Niech będzie zmienną losową dla czasu przetwarzania po zmianie, a poprzednią. Jeśli jest znacznie powyżej , to powiedziałbym, że proces się poprawił.Y P ( X < Y ) $X$$Y$$P(X < Y)$$0.5$

Teraz mam obserwacje z i obserwacje z . Obserwowana prawdopodobieństwa jest .x i X m y j$n$$x_i$$X$$m$$y_j$$Y$P = $P(X < Y)$$\hat p = \frac{1}{n m} \sum_i \sum_j 1_{x_i < y_j}$

Co mogę powiedzieć o dana uwagi i ?x i y $P(X < Y)$$x_i$$y_j$

Twoje oszacowanie jest równe statystyce Mann-Whitney podzielonej przez (dzięki, Glen!), A zatem jest równoważne statystyce sumy rang Wilcoxona (znanej również jako statystyka Wilcoxona-Manna-Whitneya) : , gdzie jest wielkością próby (zakładając brak powiązań). Możesz zatem użyć tabel / oprogramowania testu Wilcoxona i przekształcić je z powrotem w aby uzyskać przedział ufności lub wartość .$\hat{p}$$U$$mn$$W$$W = U + {n(n+1)\over{2}}$$n$$y$$U$$p$

Niech będzie wielkością próbki , = . Następnie asymptotycznie$m$$x$$N$$m+n$

$W^* = \frac{W-\frac{m(N+1)}{2}}{\sqrt{\frac{mn(N+1)}{12}}} \sim \text{N}(0,1)$

Źródło: Hollander i Wolfe , nieparametryczne metody statystyczne, w przybliżeniu str. 117, ale zapewne zaprowadzi Cię większość książek o statystykach nieparametrycznych.

@jbowman zapewnia (ładne) standardowe rozwiązanie problemu oszacowania który jest znany jako model wytrzymałości na stres .$\theta=P(X<Y)$

Inną nieparametryczną alternatywę zaproponowano w Baklizi i Eidous (2006) dla przypadku, w którym i są niezależne. Jest to opisane poniżej.$X$$Y$

Z definicji mamy to

gdzie jest CDF i jest gęstością . Następnie, przy użyciu próbek i można uzyskać jądra estymatory o i i konsekwentnie i estymator X f Y Y X Y F X f Y$F_X$$X$$f_Y$$Y$$X$$Y$$F_X$$f_Y$$\theta$

Jest to zaimplementowane w następującym kodzie R przy użyciu jądra Gaussa.

# Optimal bandwidth
h = function(x){
n = length(x)
return((4*sqrt(var(x))^5/(3*n))^(1/5))
}

# Kernel estimators of the density and the distribution
kg = function(x,data){
hb = h(data)
k = r = length(x)
for(i in 1:k) r[i] = mean(dnorm((x[i]-data)/hb))/hb
return(r )
} 

KG = function(x,data){
hb = h(data)
k = r = length(x)
for(i in 1:k) r[i] = mean(pnorm((x[i]-data)/hb))
return(r )
} 

# Baklizi and Eidous (2006) estimator
nonpest = function(dat1B,dat2B){
return( as.numeric(integrate(function(x) KG(x,dat1B)*kg(x,dat2B),-Inf,Inf)$value))  
}

# Example when X and Y are Cauchy
datx = rcauchy(100,0,1)
daty =  rcauchy(100,0,1)

nonpest(datx,daty)

Aby uzyskać przedział ufności dla , możesz pobrać próbkę ładowania tego estymatora w następujący sposób.$\theta$

# bootstrap
B=1000
p = rep(0,B)

for(j in 1:B){
dat1 =  sample(datx,length(datx),replace=T)
dat2 =  sample(daty,length(daty),replace=T)
p[j] = nonpest(dat1,dat2)
}

# histogram of the bootstrap sample
hist(p)

# A confidence interval (quantile type)
c(quantile(p,0.025),quantile(p,0.975))

Można również rozważyć inne rodzaje interwałów ładowania.

$X_i-Y_i$$P(X_i-Y_i<0) = p$$I\{X_i-Y_i<0\}$$i=1,2,..,n$$X$$X_i < Y_i$$n$ $p=P(X_i-Y_i<0)$$X/n$