K Sample Location Problem

Author

Subhendu Ghosh, Swarnadeep Datta

K SAMPLE LOCATION PROBLEM

Introduction

Suppose that \(F_1(x)\), \(F_2(x)\),…,\(F_k(x)\) be the distribution functions corresponding to k(>2) populations where, in general, \[ F_i(x) = F(x - \theta_i)\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,i=1,2,...,k\] where \(F_i(x)\) is an unknown continuous distribution and the differences between the k populations are through the values of the location parameters \(\theta_1\), \(\theta_2\),…, \(\theta_k\). Here our concern is to check whether all the distributions have same location or not, i.e, to test if, \[\theta_1= \theta_2=...= \theta_k\]

Real Life Examples

Suppose we are are conducting a survey to compare the average income of people of different districts in a particular state. Here the number of districts in that particular state is our ‘k’ and the average income of people in a district is our parameter(location parameter) \(\theta\) of interest.
We are to check whether the districts have similar average income ,i.e, \[\theta_1= \theta_2=...= \theta_k\]

Parametric Version

If we had the information that the parent distributions of the random variables \(X_1\),\(X_2\),…,\(X_k\) are normal, i.e, \[X_i\,\,\,\overset{ind}{\sim}\,\,\,\,N(\theta_i,\sigma^2)\,\,\,\, \forall\,\,\,i\,=\,1(1)k\] then we can do the ANOVA test to meet the problem in question.

What if the distributions are not Normal?

When we collect the data from practical sources we don’t know the underlying distribution.
In that case we first try to check whether the parent distribution is normal or not, but if we don’t reach to any convincing conclusion to accept the fact the parent distribution is normal we use to treat the given problem as a non-parametric testing problem, and use the Kruskal-Wallis test to check the location shift of the distributions in question,…

KRUSKAL-WALLIS TEST

Set up

Suppose we are given that, \[ X_1,_1,\, X_1,_2,...,X_1,n_1\,\,\,\,\,\,\,\,iid\,\,\, F_1(x) \] \[X_2,_1,\, X_2,_2,...,X_2,n_2\,\,\,\,\,\,\,\,iid\,\,\, F_2(x)\] \[...............................................\] \[X_k,_1,\, X_k,_2,...,X_k,n_k\,\,\,\,\,\,\,\,iid\,\,\, F_k(x)\] where, \(F_i(x)\,\,\, are\,\,\, assumed\,\,\,to\,\,\,be \,\,\,continuous\,\,\,\, \forall\,\,\,i\,\,=\,\,1(1)k\)

Testing Problem

Consider, \[ F_i(x) = F(x - \theta_i)\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,i=1,2,...,k\] Now we are to test whether the distributions have same location or not. so our appropriate null hypothesis is \[ H_0 :\,\,\,\,\,\theta_1= \theta_2=...= \theta_k=\theta\] against the alternative \[ H_1 :\,\,\,\,\theta_i\neq\,\theta_j\,\,\,\,\,\,\,for\,\,some\,\,\,i\,\neq\,j\] Note that if the locations are not equal then we reject the null hypothesis in favour of the alternate

Test Statistic

At first combine the the samples from k different populations into a combined sample and rank the observations in the according to the combined sample
Let us define \(R_{ij}\) be the rank of \(X_{ij}\) in the combined arrangements of N observations
Define \(\bar{R}_i\,\, :=\,\,\,\dfrac{1}{n_i}\sum_{j=1}^{ n_i}R_{ij},\,\,\,\,\forall\,\,i\,=\,1(1)n_i\) ,and

\(\bar{R} = \dfrac{\sum_{i=1}^{k}\sum_{j=1}^{ n_i}R_{ij}}{N} = \dfrac{N+1}{2}\)

\[\mathcal{H} = \dfrac{12}{N(N+1)}\sum_{r=1}^{n}n_i\left(\bar{R}_i - \frac{N+1}{2}\right)^2 \]

Alternative Form

The Test Statistic can alternatively be written as:-
\[\mathcal{H}\,\,\,\,=\,\,\,\dfrac{12}{N(N+1)}\sum_{r=1}^{n}n_i\bar{R}_i^2\,\,\,\, - \,\,\,\,3(N+1)\]

Distribution Free under \(H_0\)

Note that under Null hypothesis all the distributions have same location, thus the rank vector follows uniform distribution over the set of all permutations of 1 to N.
Hence the Test Statistic which depends only upon the ranks is distribution free. We can calculate the Expectation and variance of the Test statistic under null hypothesis
We shall verify this claim by checking whether data generated from different distributions have same null distribution of the test Statistic or not,..

\(\mathcal{N}(0,10)\;\;\;n_1=25,n_2=20,n_3=22,n_4=30\)

H<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rnorm(n1,0,10)
  
  n2<-20
  s2<-rnorm(n2,0,10)
  
  n3<-22
  s3<-rnorm(n3,0,10)
  
  n4<-30
  s4<-rnorm(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}

Hist<-hist(H,breaks = 50, col = "orchid2",main = "",probability = T)
curve(dchisq(x,df=3),0,20,add=TRUE)
box()

\(\mathcal{C}(0,10)\;\;\;n_1=25,n_2=20,n_3=22,n_4=30\)

H<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rcauchy(n1,0,10)
  
  n2<-20
  s2<-rcauchy(n2,0,10)
  
  n3<-22
  s3<-rcauchy(n3,0,10)
  
  n4<-30
  s4<-rcauchy(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}

Hist<-hist(H,breaks = 50, col = "orchid2",main = "",probability = T)
curve(dchisq(x,df=3),0,20,add=TRUE)
box()

\(\mathcal{logistic}(0,10)\;\;\;n_1=25,n_2=20,n_3=22,n_4=30\)

H<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rlogis(n1,0,10)
  
  n2<-20
  s2<-rlogis(n2,0,10)
  
  n3<-22
  s3<-rlogis(n3,0,10)
  
  n4<-30
  s4<-rlogis(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}

Hist<-hist(H,breaks = 50, col = "orchid2",main = "",probability = T)
curve(dchisq(x,df=3),0,20,add=TRUE)
box()

\(\mathcal{t}(3),\;\mathcal{N}(0,10),\;\mathcal{C}(0,10),\;\mathcal{logistic}(0,10)\)

\[n_1=25,n_2=20,n_3=22,n_4=30\]

H<-NULL
for(k in 1:10000)
{
 n1<-25
  s1<-rt(n1,3)
  
  n2<-20
  s2<-rcauchy(n2,0,10)
  
  n3<-22
  s3<-rnorm(n3,0,10)
  
  n4<-30
  s4<-rlogis(n4,0,10)
   
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}

Hist<-hist(H,breaks = 50, col = "orchid2",main = "",probability = T)
curve(dchisq(x,df=3),0,20,add=TRUE)
box()

\(\mathcal{exp}(1)\;\;\;n_1=25,n_2=20,n_3=22,n_4=30\)

H<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rexp(n1,1)
  
  n2<-20
  s2<-rexp(n2,1)
  
  n3<-22
  s3<-rexp(n3,1)
  
  n4<-30
  s4<-rexp(n4,1)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}

Hist<-hist(H,breaks = 50, col = "orchid2",main = "",probability = T)
curve(dchisq(x,df=3),0,20,add=TRUE)
box()

\(\mathcal{exp}(1)\;\;\;n_1=2,n_2=3,n_3=2,n_4=3\)

H<-NULL
for(k in 1:10000)
{
  n1<-2
  s1<-rexp(n1,1)
  
  n2<-3
  s2<-rexp(n2,1)
  
  n3<-2
  s3<-rexp(n3,1)
  
  n4<-3
  s4<-rexp(n4,1)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}

Hist<-hist(H,breaks = 50, col = "orchid2",main = "",probability = T)
box()

\(\mathcal{N}(1)\;\;\;n_1=2,n_2=3,n_3=2,n_4=3\)

H<-NULL
for(k in 1:10000)
{
  n1<-2
  s1<-rnorm(n1,0,1)
  
  n2<-3
  s2<-rnorm(n2,0,1)
  
  n3<-2
  s3<-rnorm(n3,0,1)
  
  n4<-3
  s4<-rnorm(n4,0,1)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}

Hist<-hist(H,breaks = 50, col = "orchid2",main = "",probability = T)
box()

Formal Tests

We shall do a pearsonian chi square test to check whether the obtained distribution are same or not \[\mathcal{N}(0,10)\,\,\,\ vs\,\,\,\,\mathcal{C}(0,10)\;\;n_1=25,n_2=20,n_3=22,n_4=30\]

H<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rnorm(n1,0,10)
  
  n2<-20
  s2<-rnorm(n2,0,10)
  
  n3<-22
  s3<-rnorm(n3,0,10)
  
  n4<-30
  s4<-rnorm(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}
c1=0
c2=0
c3=0
c4=0
c5=0
c6=0
c7=0

for(k in 1:10000){
if(H[k]<=2) {
  c1=c1+1
}
else if (H[k]>2 && H[k]<=4){
  c2=c2+1
}
else if (H[k]>4 && H[k]<=6){
  c3=c3+1
}
  else if (H[k]>6 && H[k]<=8){
    c4=c4+1
  }
  else if (H[k]>8 && H[k]<=10){
    c5=c5+1
  }
  else if (H[k]>10 && H[k]<=12){
    c6=c6+1
  }
  else if (H[k]>12){
    c7=c7+1
  }
  
}

H<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rcauchy(n1,0,10)
  
  n2<-20
  s2<-rcauchy(n2,0,10)
  
  n3<-22
  s3<-rcauchy(n3,0,10)
  
  n4<-30
  s4<-rcauchy(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}
d1=0
d2=0
d3=0
d4=0
d5=0
d6=0
d7=0

for(k in 1:10000){
if(H[k]<=2) {
  d1=d1+1
}
else if (H[k]>2 && H[k]<=4){
  d2=d2+1
}
else if (H[k]>4 && H[k]<=6){
  d3=d3+1
}
  else if (H[k]>6 && H[k]<=8){
    d4=d4+1
  }
  else if (H[k]>8 && H[k]<=10){
    d5=d5+1
  }
  else if (H[k]>10 && H[k]<=12){
    d6=d6+1
  }
  else if (H[k]>12){
    d7=d7+1
  }
  
}
chi<-(c1-d1)*(c1-d1)/d1+(c2-d2)*(c2-d2)/d2+(c3-d3)*(c3-d3)/d3+(c4-d4)*(c4-d4)/d4+(c5-d5)*(c5-d5)/d5+(c6-d6)*(c6-d6)/d6+(c7-d7)*(c7-d7)/d7
chi
[1] 16.35223
qchisq(0.01,6,lower.tail = FALSE)
[1] 16.81189

Formal Tests

We shall do a pearsonian chi square test to check whether the obtained distribution are same or not \[\mathcal{N}(0,10)\,\,\,\ vs\,\,\,\,\mathcal{C}(0,10)\;\;n_1=2,n_2=3,n_3=2,n_4=3\]

H<-NULL
for(k in 1:10000)
{
  n1<-2
  s1<-rnorm(n1,0,10)
  
  n2<-3
  s2<-rnorm(n2,0,10)
  
  n3<-2
  s3<-rnorm(n3,0,10)
  
  n4<-3
  s4<-rnorm(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}
c1=0
c2=0
c3=0
c4=0


for(k in 1:10000){
if(H[k]<=2) {
  c1=c1+1
}
else if (H[k]>2 && H[k]<=4){
  c2=c2+1
}
else if (H[k]>4 && H[k]<=6){
  c3=c3+1
}
  else if (H[k]>6){
    c4=c4+1
  }
 
  
}

H<-NULL
for(k in 1:10000)
{
  n1<-2
  s1<-rcauchy(n1,0,10)
  
  n2<-3
  s2<-rcauchy(n2,0,10)
  
  n3<-2
  s3<-rcauchy(n3,0,10)
  
  n4<-3
  s4<-rcauchy(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}
d1=0
d2=0
d3=0
d4=0


for(k in 1:10000){
if(H[k]<=2) {
  d1=d1+1
}
else if (H[k]>2 && H[k]<=4){
  d2=d2+1
}
else if (H[k]>4 && H[k]<=6){
  d3=d3+1
}
  else if (H[k]>6 ){
    d4=d4+1
  }
  
  
}
chi<-(c1-d1)*(c1-d1)/d1+(c2-d2)*(c2-d2)/d2+(c3-d3)*(c3-d3)/d3+(c4-d4)*(c4-d4)/d4
chi
[1] 6.458716
qchisq(0.01,3,lower.tail = FALSE)
[1] 11.34487

Small Sample Distribution

Here we are interested to find the small sample distribution of the Kruskal-Wallis statistic. Note that when the sample sizes \(n_i\) are small the kruskal-Wallis statistic doesn’t follow any theoritical distribution. Thus in such cases we have to draw small sample tables

\(n_1=1,n_2=2,n_3=2,n_4=2\)

library(combinat)

Attaching package: 'combinat'
The following object is masked from 'package:utils':

    combn
require(kableExtra)
Loading required package: kableExtra
Warning: package 'kableExtra' was built under R version 4.3.3
n1=1

n2=2

n3=2

n4=2

N=n1+n2+n3+n4

h = lapply(permn(1:N),FUN = function(r,n1=1,n2=2,n3=2,n4=2){
  N=n1+n2+n3+n4 

    r1<-sum(r[1:n1])/n1

      r2<-mean(r[(n1+1):(n1+n2)])

    r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  
    r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
    ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)})

h = unlist(h)

out = round(table(h)/length(h),3)

values = names(out)

probability = as.vector(out)

cdf = cumsum(probability)

send = cbind(values, probability, cdf)

kbl(send,format = 'html',caption = 'small sample Table') %>% kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
small sample Table
values probability cdf
0 0.01 0.01
0.214285714285712 0.019 0.029
0.321428571428569 0.019 0.048
0.535714285714285 0.019 0.067
0.642857142857142 0.019 0.086
0.75 0.019 0.105
0.857142857142854 0.01 0.115
1.07142857142857 0.019 0.134
1.17857142857143 0.019 0.153
1.28571428571428 0.019 0.172
1.39285714285714 0.019 0.191
1.5 0.019 0.21
1.60714285714285 0.038 0.248
1.92857142857143 0.048 0.296
2.14285714285714 0.019 0.315
2.25 0.038 0.353
2.35714285714285 0.019 0.372
2.46428571428571 0.057 0.429
2.78571428571428 0.038 0.467
2.89285714285714 0.057 0.524
3.10714285714285 0.038 0.562
3.21428571428571 0.019 0.581
3.42857142857143 0.01 0.591
3.64285714285714 0.019 0.61
3.75 0.057 0.667
3.96428571428571 0.019 0.686
4.07142857142857 0.038 0.724
4.17857142857143 0.019 0.743
4.5 0.057 0.8
4.71428571428571 0.019 0.819
4.82142857142857 0.057 0.876
5.03571428571428 0.057 0.933
5.35714285714285 0.029 0.962
5.67857142857143 0.038 1

\(n_1=2,n_2=3,n_3=2,n_4=3\)

library(combinat)
require(kableExtra)

n1=2

n2=3

n3=2

n4=3

N=n1+n2+n3+n4

h = lapply(permn(1:N),FUN = function(r,n1=1,n2=2,n3=2,n4=2){
  N=n1+n2+n3+n4 

    r1<-sum(r[1:n1])/n1

      r2<-mean(r[(n1+1):(n1+n2)])

    r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  
    r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
    ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)})

h = unlist(h)

out = table(h)/length(h)

values = names(out)

probability = as.vector(out)

cdf = cumsum(probability)

send = cbind(values, probability, cdf)

kbl(send,format = 'html',caption = 'small sample Table') %>% kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
small sample Table
values probability cdf
0 7.93650793650794e-05 7.93650793650794e-05
0.214285714285712 0.000158730158730159 0.000238095238095238
0.321428571428569 0.000158730158730159 0.000396825396825397
0.535714285714285 0.000158730158730159 0.000555555555555556
0.642857142857142 0.000158730158730159 0.000714285714285714
0.75 0.000158730158730159 0.000873015873015873
0.857142857142854 7.93650793650794e-05 0.000952380952380952
1.07142857142857 0.000158730158730159 0.00111111111111111
1.17857142857143 0.000158730158730159 0.00126984126984127
1.28571428571428 0.000158730158730159 0.00142857142857143
1.39285714285714 0.000158730158730159 0.00158730158730159
1.5 0.000158730158730159 0.00174603174603175
1.60714285714285 0.000317460317460317 0.00206349206349206
1.82142857142857 7.93650793650794e-05 0.00214285714285714
1.92857142857143 0.000396825396825397 0.00253968253968254
2.03571428571428 7.93650793650794e-05 0.00261904761904762
2.14285714285714 0.000317460317460317 0.00293650793650794
2.25 0.000396825396825397 0.00333333333333333
2.35714285714285 0.000158730158730159 0.00349206349206349
2.46428571428571 0.000476190476190476 0.00396825396825397
2.57142857142857 0.000317460317460317 0.00428571428571429
2.67857142857143 7.93650793650794e-05 0.00436507936507937
2.78571428571428 0.000317460317460317 0.00468253968253968
2.89285714285714 0.000714285714285714 0.0053968253968254
3.10714285714285 0.000555555555555556 0.00595238095238095
3.21428571428571 0.000238095238095238 0.00619047619047619
3.42857142857143 0.000396825396825397 0.00658730158730159
3.53571428571428 7.93650793650794e-05 0.00666666666666667
3.64285714285714 0.000238095238095238 0.0069047619047619
3.75 0.000793650793650794 0.0076984126984127
3.85714285714285 0.000317460317460317 0.00801587301587302
3.96428571428571 0.000396825396825397 0.00841269841269841
4.07142857142857 0.000555555555555556 0.00896825396825397
4.17857142857143 0.000476190476190476 0.00944444444444444
4.28571428571428 0.000158730158730159 0.0096031746031746
4.39285714285714 0.000238095238095238 0.00984126984126984
4.5 0.000793650793650794 0.0106349206349206
4.60714285714285 0.000317460317460317 0.010952380952381
4.71428571428571 0.000396825396825397 0.0113492063492063
4.82142857142857 0.000952380952380952 0.0123015873015873
4.92857142857143 0.000396825396825397 0.0126984126984127
5.03571428571428 0.000555555555555556 0.0132539682539683
5.14285714285714 0.000396825396825397 0.0136507936507937
5.35714285714285 0.000634920634920635 0.0142857142857143
5.46428571428571 0.000873015873015873 0.0151587301587302
5.57142857142857 0.000238095238095238 0.0153968253968254
5.67857142857143 0.00111111111111111 0.0165079365079365
5.78571428571428 0.000317460317460317 0.0168253968253968
5.89285714285714 0.000317460317460317 0.0171428571428571
6 0.000476190476190476 0.0176190476190476
6.10714285714285 0.000714285714285714 0.0183333333333333
6.21428571428571 0.000396825396825397 0.0187301587301587
6.32142857142857 0.000952380952380952 0.0196825396825397
6.42857142857143 0.000555555555555556 0.0202380952380952
6.53571428571428 0.000476190476190476 0.0207142857142857
6.64285714285714 0.000238095238095238 0.020952380952381
6.75 0.000634920634920635 0.0215873015873016
6.85714285714285 0.000714285714285714 0.0223015873015873
6.96428571428571 0.000952380952380952 0.0232539682539683
7.07142857142857 0.000634920634920635 0.0238888888888889
7.17857142857143 0.000714285714285714 0.0246031746031746
7.28571428571428 0.000634920634920635 0.0252380952380952
7.39285714285714 0.00142857142857143 0.0266666666666667
7.5 0.000317460317460317 0.026984126984127
7.60714285714285 0.000476190476190476 0.0274603174603175
7.71428571428571 0.00111111111111111 0.0285714285714286
7.82142857142857 0.000238095238095238 0.0288095238095238
7.92857142857143 0.00111111111111111 0.0299206349206349
8.03571428571428 0.00142857142857143 0.0313492063492063
8.14285714285714 0.000396825396825397 0.0317460317460317
8.25 0.00142857142857143 0.0331746031746032
8.35714285714285 0.000873015873015873 0.034047619047619
8.46428571428572 0.000714285714285714 0.0347619047619048
8.57142857142857 0.00126984126984127 0.036031746031746
8.67857142857142 0.000317460317460317 0.0363492063492063
8.78571428571428 0.000714285714285714 0.0370634920634921
8.89285714285714 0.00103174603174603 0.0380952380952381
9 0.00126984126984127 0.0393650793650794
9.10714285714285 0.000714285714285714 0.0400793650793651
9.21428571428572 0.000952380952380952 0.041031746031746
9.32142857142857 0.000634920634920635 0.0416666666666667
9.42857142857142 0.000555555555555556 0.0422222222222222
9.53571428571428 0.000714285714285714 0.0429365079365079
9.64285714285714 0.00126984126984127 0.0442063492063492
9.75 0.00142857142857143 0.0456349206349206
9.85714285714285 0.00111111111111111 0.0467460317460317
9.96428571428572 0.00253968253968254 0.0492857142857143
10.0714285714286 0.000873015873015873 0.0501587301587302
10.1785714285714 0.000793650793650794 0.050952380952381
10.2857142857143 0.00126984126984127 0.0522222222222222
10.3928571428571 0.000714285714285714 0.0529365079365079
10.5 0.00174603174603175 0.0546825396825397
10.6071428571429 0.00158730158730159 0.0562698412698413
10.7142857142857 0.000793650793650794 0.0570634920634921
10.8214285714286 0.00166666666666667 0.0587301587301587
10.9285714285714 0.00142857142857143 0.0601587301587302
11.0357142857143 0.00103174603174603 0.0611904761904762
11.1428571428571 0.00142857142857143 0.0626190476190476
11.25 0.00134920634920635 0.063968253968254
11.3571428571429 0.000634920634920635 0.0646031746031746
11.4642857142857 0.0023015873015873 0.0669047619047619
11.5714285714286 0.00150793650793651 0.0684126984126984
11.6785714285714 0.00206349206349206 0.0704761904761905
11.7857142857143 0.00142857142857143 0.0719047619047619
11.8928571428571 0.00134920634920635 0.0732539682539683
12 0.000873015873015873 0.0741269841269841
12.1071428571429 0.00103174603174603 0.0751587301587302
12.2142857142857 0.00150793650793651 0.0766666666666667
12.3214285714286 0.00166666666666667 0.0783333333333333
12.4285714285714 0.00142857142857143 0.0797619047619048
12.5357142857143 0.0030952380952381 0.0828571428571429
12.6428571428571 0.000714285714285714 0.0835714285714286
12.75 0.00150793650793651 0.0850793650793651
12.8571428571429 0.00238095238095238 0.0874603174603175
12.9642857142857 0.000476190476190476 0.0879365079365079
13.0714285714286 0.00253968253968254 0.0904761904761905
13.1785714285714 0.00293650793650794 0.0934126984126984
13.2857142857143 0.000793650793650794 0.0942063492063492
13.3928571428571 0.00293650793650794 0.0971428571428571
13.5 0.00182539682539683 0.098968253968254
13.6071428571429 0.000873015873015873 0.0998412698412698
13.7142857142857 0.00246031746031746 0.102301587301587
13.8214285714286 0.00166666666666667 0.103968253968254
13.9285714285714 0.00158730158730159 0.105555555555556
14.0357142857143 0.00285714285714286 0.108412698412698
14.1428571428571 0.00253968253968254 0.110952380952381
14.25 0.0023015873015873 0.113253968253968
14.3571428571429 0.00261904761904762 0.115873015873016
14.4642857142857 0.000952380952380952 0.116825396825397
14.5714285714286 0.00158730158730159 0.118412698412698
14.6785714285714 0.0019047619047619 0.12031746031746
14.7857142857143 0.00301587301587302 0.123333333333333
14.8928571428571 0.00158730158730159 0.124920634920635
15 0.00285714285714286 0.127777777777778
15.1071428571429 0.0030952380952381 0.130873015873016
15.2142857142857 0.000952380952380952 0.131825396825397
15.3214285714286 0.0019047619047619 0.133730158730159
15.4285714285714 0.00246031746031746 0.136190476190476
15.5357142857143 0.00103174603174603 0.137222222222222
15.6428571428571 0.00166666666666667 0.138888888888889
15.75 0.00412698412698413 0.143015873015873
15.8571428571429 0.00198412698412698 0.145
15.9642857142857 0.00468253968253968 0.14968253968254
16.0714285714286 0.00277777777777778 0.152460317460317
16.1785714285714 0.00111111111111111 0.153571428571429
16.2857142857143 0.00222222222222222 0.155793650793651
16.3928571428571 0.0023015873015873 0.158095238095238
16.5 0.00285714285714286 0.160952380952381
16.6071428571429 0.00468253968253968 0.165634920634921
16.7142857142857 0.00126984126984127 0.166904761904762
16.8214285714286 0.00285714285714286 0.169761904761905
16.9285714285714 0.00174603174603175 0.171507936507936
17.0357142857143 0.00317460317460317 0.17468253968254
17.1428571428571 0.00174603174603175 0.176428571428571
17.25 0.00253968253968254 0.178968253968254
17.3571428571429 0.0026984126984127 0.181666666666667
17.4642857142857 0.00222222222222222 0.183888888888889
17.5714285714286 0.00246031746031746 0.186349206349206
17.6785714285714 0.00444444444444444 0.190793650793651
17.7857142857143 0.00150793650793651 0.192301587301587
17.8928571428571 0.00214285714285714 0.194444444444444
18 0.00420634920634921 0.198650793650794
18.1071428571429 0.000952380952380952 0.199603174603175
18.2142857142857 0.00380952380952381 0.203412698412698
18.3214285714286 0.00396825396825397 0.207380952380952
18.4285714285714 0.0030952380952381 0.21047619047619
18.5357142857143 0.00444444444444444 0.214920634920635
18.6428571428571 0.00341269841269841 0.218333333333333
18.75 0.00126984126984127 0.219603174603175
18.8571428571429 0.0030952380952381 0.222698412698413
18.9642857142857 0.00182539682539683 0.22452380952381
19.0714285714286 0.0019047619047619 0.226428571428571
19.1785714285714 0.00365079365079365 0.230079365079365
19.2857142857143 0.00380952380952381 0.233888888888889
19.3928571428571 0.0023015873015873 0.236190476190476
19.5 0.0030952380952381 0.239285714285714
19.6071428571429 0.00206349206349206 0.241349206349206
19.7142857142857 0.00222222222222222 0.243571428571429
19.8214285714286 0.00206349206349206 0.245634920634921
19.9285714285714 0.00285714285714286 0.248492063492063
20.0357142857143 0.00301587301587302 0.251507936507937
20.1428571428571 0.00301587301587302 0.25452380952381
20.25 0.00571428571428571 0.260238095238095
20.3571428571429 0.00214285714285714 0.262380952380952
20.4642857142857 0.00206349206349206 0.264444444444444
20.5714285714286 0.00317460317460317 0.267619047619048
20.6785714285714 0.00174603174603175 0.269365079365079
20.7857142857143 0.0030952380952381 0.272460317460317
20.8928571428571 0.00634920634920635 0.278809523809524
21 0.00174603174603175 0.280555555555556
21.1071428571429 0.00611111111111111 0.286666666666667
21.2142857142857 0.00341269841269841 0.290079365079365
21.3214285714286 0.00253968253968254 0.292619047619048
21.4285714285714 0.00349206349206349 0.296111111111111
21.5357142857143 0.00198412698412698 0.298095238095238
21.6428571428571 0.00253968253968254 0.300634920634921
21.75 0.00547619047619048 0.306111111111111
21.8571428571429 0.0026984126984127 0.308809523809524
21.9642857142857 0.00325396825396825 0.312063492063492
22.0714285714286 0.00238095238095238 0.314444444444444
22.1785714285714 0.00317460317460317 0.317619047619048
22.2857142857143 0.00222222222222222 0.31984126984127
22.3928571428571 0.00246031746031746 0.322301587301587
22.5 0.00404761904761905 0.326349206349206
22.6071428571429 0.00317460317460317 0.32952380952381
22.7142857142857 0.00452380952380952 0.334047619047619
22.8214285714286 0.00507936507936508 0.339126984126984
22.9285714285714 0.00222222222222222 0.341349206349206
23.0357142857143 0.0019047619047619 0.343253968253968
23.1428571428571 0.00380952380952381 0.347063492063492
23.25 0.00198412698412698 0.349047619047619
23.3571428571429 0.0046031746031746 0.353650793650794
23.4642857142857 0.00515873015873016 0.358809523809524
23.5714285714286 0.00222222222222222 0.361031746031746
23.6785714285714 0.005 0.366031746031746
23.7857142857143 0.00365079365079365 0.36968253968254
23.8928571428571 0.00158730158730159 0.371269841269841
24 0.00492063492063492 0.376190476190476
24.1071428571429 0.00349206349206349 0.37968253968254
24.2142857142857 0.0019047619047619 0.381587301587302
24.3214285714286 0.00507936507936508 0.386666666666667
24.4285714285714 0.00452380952380952 0.391190476190476
24.5357142857143 0.00333333333333333 0.39452380952381
24.6428571428571 0.00404761904761905 0.398571428571429
24.75 0.00206349206349206 0.400634920634921
24.8571428571429 0.0023015873015873 0.402936507936508
24.9642857142857 0.00301587301587302 0.405952380952381
25.0714285714286 0.00436507936507937 0.41031746031746
25.1785714285714 0.00396825396825397 0.414285714285714
25.2857142857143 0.00285714285714286 0.417142857142857
25.3928571428571 0.00571428571428571 0.422857142857143
25.5 0.00198412698412698 0.42484126984127
25.6071428571429 0.00301587301587302 0.427857142857143
25.7142857142857 0.00452380952380952 0.432380952380952
25.8214285714286 0.00119047619047619 0.433571428571429
25.9285714285714 0.00492063492063492 0.438492063492063
26.0357142857143 0.00571428571428571 0.444206349206349
26.1428571428571 0.0030952380952381 0.447301587301587
26.25 0.0053968253968254 0.452698412698413
26.3571428571429 0.00341269841269841 0.456111111111111
26.4642857142857 0.00182539682539683 0.457936507936508
26.5714285714286 0.00373015873015873 0.461666666666667
26.6785714285714 0.00222222222222222 0.463888888888889
26.7857142857143 0.00349206349206349 0.467380952380952
26.8928571428571 0.00444444444444444 0.471825396825397
27 0.00420634920634921 0.476031746031746
27.1071428571429 0.00484126984126984 0.480873015873016
27.2142857142857 0.00357142857142857 0.484444444444444
27.3214285714286 0.0026984126984127 0.487142857142857
27.4285714285714 0.00182539682539683 0.488968253968254
27.5357142857143 0.00277777777777778 0.491746031746032
27.6428571428571 0.00373015873015873 0.49547619047619
27.75 0.0046031746031746 0.500079365079365
27.8571428571429 0.00341269841269841 0.503492063492064
27.9642857142857 0.00595238095238095 0.509444444444444
28.0714285714286 0.00246031746031746 0.511904761904762
28.1785714285714 0.00325396825396825 0.51515873015873
28.2857142857143 0.00571428571428571 0.520873015873016
28.3928571428571 0.00198412698412698 0.522857142857143
28.5 0.00388888888888889 0.526746031746032
28.6071428571429 0.00579365079365079 0.532539682539683
28.7142857142857 0.00222222222222222 0.534761904761905
28.8214285714286 0.00603174603174603 0.540793650793651
28.9285714285714 0.00253968253968254 0.543333333333333
29.0357142857143 0.00134920634920635 0.54468253968254
29.1428571428571 0.00404761904761905 0.548730158730159
29.25 0.00277777777777778 0.551507936507937
29.3571428571429 0.00158730158730159 0.553095238095238
29.4642857142857 0.00626984126984127 0.559365079365079
29.5714285714286 0.00349206349206349 0.562857142857143
29.6785714285714 0.0030952380952381 0.565952380952381
29.7857142857143 0.00476190476190476 0.570714285714286
29.8928571428571 0.0023015873015873 0.573015873015873
30 0.00325396825396825 0.576269841269841
30.1071428571429 0.00253968253968254 0.578809523809524
30.2142857142857 0.00547619047619048 0.584285714285714
30.3214285714286 0.00349206349206349 0.587777777777778
30.4285714285714 0.00404761904761905 0.591825396825397
30.5357142857143 0.00603174603174603 0.597857142857143
30.6428571428571 0.00198412698412698 0.59984126984127
30.75 0.00293650793650794 0.602777777777778
30.8571428571429 0.00396825396825397 0.606746031746032
30.9642857142857 0.00150793650793651 0.608253968253968
31.0714285714286 0.00420634920634921 0.612460317460317
31.1785714285714 0.00531746031746032 0.617777777777778
31.2857142857143 0.00246031746031746 0.620238095238095
31.3928571428571 0.00484126984126984 0.625079365079365
31.5 0.005 0.630079365079365
31.6071428571429 0.00166666666666667 0.631746031746032
31.7142857142857 0.0026984126984127 0.634444444444444
31.8214285714286 0.0026984126984127 0.637142857142857
31.9285714285714 0.00357142857142857 0.640714285714286
32.0357142857143 0.0053968253968254 0.646111111111111
32.1428571428571 0.00388888888888889 0.65
32.25 0.00373015873015873 0.653730158730159
32.3571428571429 0.00277777777777778 0.656507936507936
32.4642857142857 0.00349206349206349 0.66
32.5714285714286 0.00222222222222222 0.662222222222222
32.6785714285714 0.00333333333333333 0.665555555555556
32.7857142857143 0.00253968253968254 0.668095238095238
32.8928571428571 0.00444444444444444 0.672539682539683
33 0.00357142857142857 0.676111111111111
33.1071428571429 0.00626984126984127 0.682380952380952
33.2142857142857 0.00166666666666667 0.684047619047619
33.3214285714286 0.00174603174603175 0.685793650793651
33.4285714285714 0.00388888888888889 0.68968253968254
33.5357142857143 0.00126984126984127 0.690952380952381
33.6428571428571 0.00349206349206349 0.694444444444444
33.75 0.00634920634920635 0.700793650793651
33.8571428571429 0.00142857142857143 0.702222222222222
33.9642857142857 0.00547619047619048 0.707698412698413
34.0714285714286 0.00373015873015873 0.711428571428571
34.1785714285714 0.00126984126984127 0.712698412698413
34.2857142857143 0.00412698412698413 0.716825396825397
34.3928571428571 0.0019047619047619 0.718730158730159
34.5 0.0030952380952381 0.721825396825397
34.6071428571429 0.00428571428571429 0.726111111111111
34.7142857142857 0.0030952380952381 0.729206349206349
34.8214285714286 0.00341269841269841 0.732619047619048
34.9285714285714 0.00277777777777778 0.735396825396825
35.0357142857143 0.0023015873015873 0.737698412698413
35.1428571428571 0.00261904761904762 0.74031746031746
35.25 0.00261904761904762 0.742936507936508
35.3571428571429 0.00412698412698413 0.747063492063492
35.4642857142857 0.0023015873015873 0.749365079365079
35.5714285714286 0.0030952380952381 0.752460317460317
35.6785714285714 0.00373015873015873 0.756190476190476
35.7857142857143 0.00198412698412698 0.758174603174603
35.8928571428571 0.00182539682539683 0.76
36 0.00333333333333333 0.763333333333333
36.1071428571429 0.00119047619047619 0.764523809523809
36.2142857142857 0.00293650793650794 0.767460317460318
36.3214285714286 0.00492063492063492 0.772380952380952
36.4285714285714 0.00214285714285714 0.77452380952381
36.5357142857143 0.00380952380952381 0.778333333333333
36.6428571428571 0.00277777777777778 0.781111111111111
36.75 0.00150793650793651 0.782619047619048
36.8571428571429 0.00293650793650794 0.785555555555556
36.9642857142857 0.00253968253968254 0.788095238095238
37.0714285714286 0.00142857142857143 0.78952380952381
37.1785714285714 0.005 0.79452380952381
37.2857142857143 0.00214285714285714 0.796666666666667
37.3928571428571 0.00222222222222222 0.798888888888889
37.5 0.00293650793650794 0.801825396825397
37.6071428571429 0.0023015873015873 0.804126984126984
37.7142857142857 0.00238095238095238 0.806507936507937
37.8214285714286 0.00246031746031746 0.808968253968254
37.9285714285714 0.00261904761904762 0.811587301587302
38.0357142857143 0.0023015873015873 0.813888888888889
38.1428571428571 0.00134920634920635 0.815238095238095
38.25 0.00476190476190476 0.82
38.3571428571429 0.00198412698412698 0.821984126984127
38.4642857142857 0.00142857142857143 0.823412698412698
38.5714285714286 0.00277777777777778 0.826190476190476
38.6785714285714 0.000714285714285714 0.826904761904762
38.7857142857143 0.00341269841269841 0.83031746031746
38.8928571428571 0.00436507936507937 0.83468253968254
39 0.00126984126984127 0.835952380952381
39.1071428571429 0.00404761904761905 0.84
39.2142857142857 0.00174603174603175 0.841746031746032
39.3214285714286 0.000873015873015873 0.842619047619048
39.4285714285714 0.00277777777777778 0.845396825396825
39.5357142857143 0.00142857142857143 0.846825396825397
39.6428571428571 0.00158730158730159 0.848412698412698
39.75 0.0030952380952381 0.851507936507937
39.8571428571429 0.0026984126984127 0.854206349206349
39.9642857142857 0.0019047619047619 0.856111111111111
40.0714285714286 0.00214285714285714 0.858253968253968
40.1785714285714 0.00150793650793651 0.859761904761905
40.2857142857143 0.00126984126984127 0.861031746031746
40.3928571428571 0.00150793650793651 0.862539682539683
40.5 0.00325396825396825 0.865793650793651
40.6071428571428 0.0026984126984127 0.868492063492064
40.7142857142857 0.00174603174603175 0.870238095238095
40.8214285714286 0.00301587301587302 0.873253968253968
40.9285714285714 0.00126984126984127 0.874523809523809
41.0357142857143 0.00103174603174603 0.875555555555556
41.1428571428571 0.00253968253968254 0.878095238095238
41.25 0.000714285714285714 0.878809523809524
41.3571428571428 0.00134920634920635 0.88015873015873
41.4642857142857 0.00317460317460317 0.883333333333333
41.5714285714286 0.00126984126984127 0.884603174603175
41.6785714285714 0.00277777777777778 0.887380952380952
41.7857142857143 0.0023015873015873 0.88968253968254
41.8928571428571 0.000555555555555556 0.890238095238095
42 0.00222222222222222 0.892460317460318
42.1071428571428 0.00198412698412698 0.894444444444444
42.2142857142857 0.00174603174603175 0.896190476190476
42.3214285714286 0.0030952380952381 0.899285714285714
42.4285714285714 0.00111111111111111 0.900396825396825
42.5357142857143 0.00222222222222222 0.902619047619048
42.6428571428571 0.00134920634920635 0.903968253968254
42.75 0.00150793650793651 0.90547619047619
42.8571428571428 0.000952380952380952 0.906428571428571
42.9642857142857 0.00119047619047619 0.907619047619048
43.0714285714286 0.0019047619047619 0.90952380952381
43.1785714285714 0.00174603174603175 0.911269841269841
43.2857142857143 0.00119047619047619 0.912460317460318
43.3928571428571 0.00253968253968254 0.915
43.5 0.000793650793650794 0.915793650793651
43.6071428571428 0.00134920634920635 0.917142857142857
43.7142857142857 0.00222222222222222 0.919365079365079
43.8214285714286 0.000634920634920635 0.92
43.9285714285714 0.00246031746031746 0.922460317460317
44.0357142857143 0.0023015873015873 0.924761904761905
44.1428571428571 0.000634920634920635 0.925396825396825
44.25 0.0023015873015873 0.927698412698413
44.3571428571428 0.00206349206349206 0.929761904761905
44.4642857142857 0.000714285714285714 0.93047619047619
44.5714285714286 0.00111111111111111 0.931587301587302
44.6785714285714 0.000952380952380952 0.932539682539683
44.7857142857143 0.000396825396825397 0.932936507936508
44.8928571428571 0.00222222222222222 0.93515873015873
45 0.00166666666666667 0.936825396825397
45.1071428571428 0.00103174603174603 0.937857142857143
45.2142857142857 0.00182539682539683 0.93968253968254
45.3214285714286 0.000873015873015873 0.940555555555556
45.4285714285714 0.00119047619047619 0.941746031746032
45.5357142857143 0.00103174603174603 0.942777777777778
45.6428571428571 0.00111111111111111 0.943888888888889
45.75 0.00126984126984127 0.94515873015873
45.8571428571428 0.00111111111111111 0.946269841269841
45.9642857142857 0.00198412698412698 0.948253968253968
46.0714285714286 0.000396825396825397 0.948650793650794
46.1785714285714 0.000634920634920635 0.949285714285714
46.2857142857143 0.00158730158730159 0.950873015873016
46.3928571428571 0.000555555555555556 0.951428571428571
46.5 0.00150793650793651 0.952936507936508
46.6071428571428 0.00206349206349206 0.955
46.7142857142857 0.000238095238095238 0.955238095238095
46.8214285714286 0.00198412698412698 0.957222222222222
46.9285714285714 0.00103174603174603 0.958253968253968
47.0357142857143 0.000476190476190476 0.958730158730159
47.1428571428571 0.00119047619047619 0.959920634920635
47.25 0.000555555555555556 0.96047619047619
47.3571428571428 0.000952380952380952 0.961428571428571
47.4642857142857 0.00174603174603175 0.963174603174603
47.5714285714286 0.000793650793650794 0.963968253968254
47.6785714285714 0.00126984126984127 0.965238095238095
47.7857142857143 0.000634920634920635 0.965873015873016
47.8928571428571 0.000714285714285714 0.966587301587302
48 0.000952380952380952 0.967539682539683
48.1071428571428 0.000238095238095238 0.967777777777778
48.2142857142857 0.00150793650793651 0.969285714285714
48.3214285714286 0.000793650793650794 0.970079365079365
48.4285714285714 0.000793650793650794 0.970873015873016
48.5357142857143 0.00158730158730159 0.972460317460317
48.6428571428571 0.000396825396825397 0.972857142857143
48.75 0.000952380952380952 0.973809523809524
48.8571428571428 0.000952380952380952 0.974761904761905
48.9642857142857 0.000158730158730159 0.974920634920635
49.0714285714286 0.000873015873015873 0.975793650793651
49.1785714285714 0.00142857142857143 0.977222222222222
49.2857142857143 0.000396825396825397 0.977619047619048
49.3928571428571 0.000714285714285714 0.978333333333333
49.5 0.000793650793650794 0.979126984126984
49.6071428571428 0.000238095238095238 0.979365079365079
49.7142857142857 0.000952380952380952 0.98031746031746
49.8214285714286 0.000793650793650794 0.981111111111111
49.9285714285714 0.000238095238095238 0.981349206349206
50.0357142857143 0.00134920634920635 0.982698412698413
50.1428571428571 0.000873015873015873 0.983571428571429
50.25 0.000555555555555556 0.984126984126984
50.3571428571428 0.000952380952380952 0.985079365079365
50.4642857142857 0.000317460317460317 0.985396825396825
50.5714285714286 0.000634920634920635 0.986031746031746
50.6785714285714 0.000476190476190476 0.986507936507937
50.7857142857143 0.000396825396825397 0.986904761904762
50.8928571428571 0.000873015873015873 0.987777777777778
51 0.000238095238095238 0.988015873015873
51.1071428571428 0.00126984126984127 0.989285714285714
51.2142857142857 0.000158730158730159 0.989444444444444
51.3214285714286 0.000158730158730159 0.989603174603175
51.4285714285714 0.000873015873015873 0.990476190476191
51.6428571428571 0.000317460317460317 0.990793650793651
51.75 0.00142857142857143 0.992222222222222
51.8571428571428 0.000158730158730159 0.992380952380952
51.9642857142857 0.000873015873015873 0.993253968253968
52.0714285714286 0.000158730158730159 0.993412698412698
52.1785714285714 0.000396825396825397 0.993809523809524
52.2857142857143 0.000317460317460317 0.994126984126984
52.3928571428571 0.000476190476190476 0.994603174603175
52.5 0.000317460317460317 0.994920634920635
52.6071428571428 0.000396825396825397 0.99531746031746
52.7142857142857 0.000317460317460317 0.995634920634921
52.8214285714286 0.000238095238095238 0.995873015873016
52.9285714285714 7.93650793650794e-05 0.995952380952381
53.0357142857143 0.000317460317460317 0.996269841269841
53.1428571428571 0.000158730158730159 0.996428571428571
53.25 0.000476190476190476 0.996904761904762
53.3571428571428 0.000238095238095238 0.997142857142857
53.4642857142857 0.000158730158730159 0.997301587301587
53.5714285714286 0.000317460317460317 0.997619047619048
53.6785714285714 0.000238095238095238 0.997857142857143
54 0.000476190476190476 0.998333333333333
54.2142857142857 0.000158730158730159 0.998492063492064
54.3214285714286 0.000476190476190476 0.998968253968254
54.5357142857143 0.000476190476190476 0.999444444444444
54.8571428571428 0.000238095238095238 0.99968253968254
55.1785714285714 0.000317460317460317 1

Limiting Distribution

We know that the Kruskal-Wallis statistic asymptotically follows a chi square distribution with degrees of freedom k-1 when all the sample sizes are sufficiently large and none of the sizes are too small.

Now let us check whether this is valid or not.

\(\mathcal{N}(0,10)\;\;\;n_1=25,n_2=20,n_3=22,n_4=30\)

Here we shall check only for normal distribution as we have already shown that the statistic is distribution free under null hypothesis, hence showing for normal only will suffice.

H<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rnorm(n1,0,10)
  
  n2<-20
  s2<-rnorm(n2,0,10)
  
  n3<-22
  s3<-rnorm(n3,0,10)
  
  n4<-30
  s4<-rnorm(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}
 c1=0
c2=0
c3=0
c4=0
c5=0
c6=0
c7=0

for(k in 1:10000){
if(H[k]<=2) {
  c1=c1+1
}
else if (H[k]>2 && H[k]<=4){
  c2=c2+1
}
else if (H[k]>4 && H[k]<=6){
  c3=c3+1
}
  else if (H[k]>6 && H[k]<=8){
    c4=c4+1
  }
  else if (H[k]>8 && H[k]<=10){
    c5=c5+1
  }
  else if (H[k]>10 && H[k]<=12){
    c6=c6+1
  }
  else if (H[k]>12){
    c7=c7+1
  }
  
}
d1<-pchisq(2,3)*10000
d2<-(pchisq(4,3)-pchisq(2,3))*10000
d3<-(pchisq(6,3)-pchisq(4,3))*10000
d4<-(pchisq(8,3)-pchisq(6,3))*10000
d5<-(pchisq(10,3)-pchisq(8,3))*10000
d6<-(pchisq(12,3)-pchisq(10,3))*10000
d7<-(1-pchisq(12,3))*10000

chi<-(c1-d1)*(c1-d1)/d1+(c2-d2)*(c2-d2)/d2+(c3-d3)*(c3-d3)/d3+(c4-d4)*(c4-d4)/d4+(c5-d5)*(c5-d5)/d5+(c6-d6)*(c6-d6)/d6+(c7-d7)*(c7-d7)/d7
chi
[1] 2.975674
qchisq(0.01,6,lower.tail = FALSE)
[1] 16.81189

\(\mathcal{N}(0,10)\;\;\;n_1=2 ,n_2=3,n_3=2,n_4=3\)

H<-NULL
for(k in 1:10000)
{
  n1<-22
  s1<-rnorm(n1,0,10)
  
  n2<-3
  s2<-rnorm(n2,0,10)
  
  n3<-2
  s3<-rnorm(n3,0,10)
  
  n4<-3
  s4<-rnorm(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}
 c1=0
c2=0
c3=0
c4=0


for(k in 1:10000){
if(H[k]<=2) {
  c1=c1+1
}
else if (H[k]>2 && H[k]<=4){
  c2=c2+1
}
else if (H[k]>4 && H[k]<=6){
  c3=c3+1
}
  else if (H[k]>6){
    c4=c4+1
  }
 
  
}
d1<-pchisq(2,3)*10000
d2<-(pchisq(4,3)-pchisq(2,3))*10000
d3<-(pchisq(6,3)-pchisq(4,3))*10000
d4<-(1-pchisq(6,3))*10000

chi<-(c1-d1)*(c1-d1)/d1+(c2-d2)*(c2-d2)/d2+(c3-d3)*(c3-d3)/d3+(c4-d4)*(c4-d4)/d4
chi
[1] 125.1252
qchisq(0.01,3,lower.tail = FALSE)
[1] 11.34487

Checking Robustness

Now we are interested to check whether this test performs well evn in presence of outliers or not, i.e, the test is Robust or not.
Here what we shall we shall randomly mix some outliers with the original data points. Now if the distribution of the statistic under \(H_0\) still remians same we can say that the Test is Robust,…

\(n_1=25 ,n_2=20,n_3=22,n_4=30\)

set.seed(1234)
H<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rcauchy(n1,0,10) + runif(n1,min= 10, max = 20)
  
  n2<-20
  s2<-rcauchy(n2,0,10)+ runif(n2,min= 10, max = 20)
  
  n3<-22
  s3<-rcauchy(n3,0,10)+ runif(n3,min= 10, max = 20)
  
  n4<-30
  s4<-rcauchy(n4,0,10)+ runif(n4,min= 10, max = 20)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}


c1=0
c2=0
c3=0
c4=0
c5=0
c6=0
c7=0

for(k in 1:10000){
  if(H[k]<=2) {
    c1=c1+1
  }
  else if (H[k]>2 && H[k]<=4){
    c2=c2+1
  }
  else if (H[k]>4 && H[k]<=6){
    c3=c3+1
  }
  else if (H[k]>6 && H[k]<=8){
    c4=c4+1
  }
  else if (H[k]>8 && H[k]<=10){
    c5=c5+1
  }
  else if (H[k]>10 && H[k]<=12){
    c6=c6+1
  }
  else if (H[k]>12){
    c7=c7+1
  }
  
}


H<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rcauchy(n1,0,10)
  
  n2<-20
  s2<-rcauchy(n2,0,10)
  
  n3<-22
  s3<-rcauchy(n3,0,10)
  
  n4<-30
  s4<-rcauchy(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}

e1=0
e2=0
e3=0
e4=0
e5=0
e6=0
e7=0

for(k in 1:10000){
  if(H[k]<=2) {
    e1=e1+1
  }
  else if (H[k]>2 && H[k]<=4){
    e2=e2+1
  }
  else if (H[k]>4 && H[k]<=6){
    e3=e3+1
  }
  else if (H[k]>6 && H[k]<=8){
    e4=e4+1
  }
  else if (H[k]>8 && H[k]<=10){
    e5=e5+1
  }
  else if (H[k]>10 && H[k]<=12){
    e6=e6+1
  }
  else if (H[k]>12){
    e7=e7+1
  }
  
}


chi<-(c1-e1)*(c1-e1)/e1+(c2-e2)*(c2-e2)/e2+(c3-e3)*(c3-e3)/e3+(c4-e4)*(c4-e4)/e4+(c5-e5)*(c5-e5)/e5+(c6-e6)*(c6-e6)/e6+(c7-e7)*(c7-e7)/e7
chi
[1] 4.204111
qchisq(0.01,6,lower.tail = FALSE)
[1] 16.81189

\(n_1=2 ,n_2=3,n_3=2,n_4=3\)

H<-NULL
for(k in 1:10000)
{
  n1<-2
  s1<-rcauchy(n1,0,10) + runif(n1,min= 10, max = 20)
  
  n2<-3
  s2<-rcauchy(n2,0,10)+ runif(n2,min= 10, max = 20)
  
  n3<-2
  s3<-rcauchy(n3,0,10)+ runif(n3,min= 10, max = 20)
  
  n4<-3
  s4<-rcauchy(n4,0,10)+ runif(n4,min= 10, max = 20)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}


c1=0
c2=0
c3=0
c4=0

for(k in 1:10000){
  if(H[k]<=2) {
    c1=c1+1
  }
  else if (H[k]>2 && H[k]<=4){
    c2=c2+1
  }
  else if (H[k]>4 && H[k]<=6){
    c3=c3+1
  }
  else if (H[k]>6 ){
    c4=c4+1
  }
  
  
}


H<-NULL
for(k in 1:10000)
{
  n1<-2
  s1<-rcauchy(n1,0,10)
  
  n2<-3
  s2<-rcauchy(n2,0,10)
  
  n3<-2
  s3<-rcauchy(n3,0,10)
  
  n4<-3
  s4<-rcauchy(n4,0,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  H[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}

e1=0
e2=0
e3=0
e4=0


for(k in 1:10000){
  if(H[k]<=2) {
    e1=e1+1
  }
  else if (H[k]>2 && H[k]<=4){
    e2=e2+1
  }
  else if (H[k]>4 && H[k]<=6){
    e3=e3+1
  }
  else if (H[k]>6 ){
    e4=e4+1
  }
  
  
}


chi<-(c1-e1)*(c1-e1)/e1+(c2-e2)*(c2-e2)/e2+(c3-e3)*(c3-e3)/e3+(c4-e4)*(c4-e4)/e4+(c5-e5)*(c5-e5)/e5+(c6-e6)*(c6-e6)/e6+(c7-e7)*(c7-e7)/e7
chi
[1] 7.087954
qchisq(0.01,3,lower.tail = FALSE)
[1] 11.34487

Power Comparisons

We shall consider the cases of small sample and large sample separately.
- In case small sample tests we shall consider the cutoff point corresponding to the exact test.
- Whereas in case of large sample tests we shall consider the cut off point corresponding to the asymptotic test.

Sample Sample Power

\(n_1=2 ,n_2=3,n_3=2,n_4=3\)

When we are varying \(\theta_2\) in 0,1,2,3,4,5 We considered the small sample tests having level 0.05

library(kableExtra)
P=NULL
for (l in c(0,1,2,3,4,5)){
  

h<-NULL
for(k in 1:10000)
{
  n1<-1
  s1<-rnorm(n1,0,1)
  
  n2<-2
  s2<-rnorm(n2,l,1)
  
  n3<-2
  s3<-rnorm(n3,0,1)
  
  n4<-2
  s4<-rnorm(n4,0,1)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  h[k]<- 12*(n1*(r1-(N+1)/2)*(r1-(N+1)/2) + n2*(r2-(N+1)/2)*(r2-(N+1)/2) + n3*(r3-(N+1)/2)*(r3-(N+1)/2) + n4*(r4-(N+1)/2)*(r4-(N+1)/2))/(N*(N+1)) 
}

  I=0
  for(k in 1:10000){
    if (h[k]>5.35714285714285){
      I=I+1
    }
  }
  P[l+1]= I/10000
  
}


l = 0:5
Power = P
send = cbind(l,Power)
colnames(send)[1] = 'theta_2'

kbl(send,format = 'html',caption = 'small sample Table') %>%
  kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
small sample Table
theta_2 Power
0 0.0668
1 0.1191
2 0.2188
3 0.2960
4 0.3262
5 0.3365

Power Curve

plot(0:5,Power,type = 'l',col = 'navyblue',lwd=2,ylab = 'power',xlab='theta 2')
grid()

\(n_1=2 ,n_2=3,n_3=2,n_4=3\)

When we are varying \(\theta_2\) in 0,1,2,3,4,5 and \(\theta_3\) in 0,1,2,3,4 We considered the small sample tests having level 0.05

library(kableExtra)
P=NULL
for(m in c(0,1,2,3,4,5)){
for (l in c(0,1,2,3,4,5)){
  

h<-NULL
for(k in 1:10000)
{
  n1<-1
  s1<-rnorm(n1,0,1)
  
  n2<-2
  s2<-rnorm(n2,l,1)
  
  n3<-2
  s3<-rnorm(n3,m,1)
  
  n4<-2
  s4<-rnorm(n4,0,1)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  h[k]<- 12*(n1*(r1-(N+1)/2)*(r1-(N+1)/2) + n2*(r2-(N+1)/2)*(r2-(N+1)/2) + n3*(r3-(N+1)/2)*(r3-(N+1)/2) + n4*(r4-(N+1)/2)*(r4-(N+1)/2))/(N*(N+1)) 
}

  I=0
  for(k in 1:10000){
    if (h[k]>5.35714285714285){
      I=I+1
    }
  }
  P[m*6+(l+1)]= I/10000
  
}
}
Power = P
m = rep(0:5,each = 6)
l = rep(0:5, times = 6)

send = cbind(m, l, Power)

kbl(send,format = 'html',caption = 'Small sample power',
    col.names = c('$\\theta_3$','$\\theta_2$','Power'  )
    
    ) %>%
  kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
Small sample power
$\theta_3$ $\theta_2$ Power
0 0 0.0686
0 1 0.1130
0 2 0.2215
0 3 0.3056
0 4 0.3316
0 5 0.3252
1 0 0.1111
1 1 0.1089
1 2 0.2068
1 3 0.3647
1 4 0.4592
1 5 0.4910
2 0 0.2149
2 1 0.2098
2 2 0.2132
2 3 0.3709
2 4 0.5886
2 5 0.7269
3 0 0.3015
3 1 0.3606
3 2 0.3611
3 3 0.2966
3 4 0.4566
3 5 0.7291
4 0 0.3344
4 1 0.4483
4 2 0.5998
4 3 0.4619
4 4 0.3313
4 5 0.4848
5 0 0.3354
5 1 0.4965
5 2 0.7293
5 3 0.7394
5 4 0.4961
5 5 0.3299

Power Curve

plot(Power,type = 'l',col = 'navyblue',lwd=2,ylab = 'power')
grid()

\(n_1=2 ,n_2=3,n_3=2,n_4=3\)

When we are varying \(\theta_2\) in 0,1,2,3,4,5 and \(\theta_3\) in 0,1,2,3,4 & \(\theta_4\) in 0,1,2,3,4,5 We considered the small sample tests having level 0.05

P=NULL
for(o in c(0,1,2,3,4,5)){
for(m in c(0,1,2,3,4,5)){
for (l in c(0,1,2,3,4,5)){
  

h<-NULL
for(k in 1:10000)
{
  n1<-1
  s1<-rnorm(n1,0,1)
  
  n2<-2
  s2<-rnorm(n2,l,1)
  
  n3<-2
  s3<-rnorm(n3,m,1)
  
  n4<-2
  s4<-rnorm(n4,o,1)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  h[k]<- 12*(n1*(r1-(N+1)/2)*(r1-(N+1)/2) + n2*(r2-(N+1)/2)*(r2-(N+1)/2) + n3*(r3-(N+1)/2)*(r3-(N+1)/2) + n4*(r4-(N+1)/2)*(r4-(N+1)/2))/(N*(N+1)) 
}

  I=0
  for(k in 1:10000){
    if (h[k]>5.35714285714285){
      I=I+1
    }
  }
  P[o*36+m*6+(l+1)]= I/10000
  
}
}
}
Power = P
o = rep(0:5,each = 36)
m = rep(0:5,each = 6,times=6)
l = rep(0:5, times = 36)

send = cbind(o,m, l, Power)

kbl(send,format = 'html',caption = 'Small sample power',
    col.names = c('$\\theta_4$','$\\theta_3$','$\\theta_2$','Power'  )
    
    ) %>%
  kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
Small sample power
$\theta_4$ $\theta_3$ $\theta_2$ Power
0 0 0 0.0649
0 0 1 0.1146
0 0 2 0.2261
0 0 3 0.3004
0 0 4 0.3315
0 0 5 0.3333
0 1 0 0.1125
0 1 1 0.1112
0 1 2 0.2118
0 1 3 0.3607
0 1 4 0.4520
0 1 5 0.4876
0 2 0 0.2237
0 2 1 0.2067
0 2 2 0.2215
0 2 3 0.3673
0 2 4 0.5907
0 2 5 0.7303
0 3 0 0.2972
0 3 1 0.3600
0 3 2 0.3639
0 3 3 0.3048
0 3 4 0.4519
0 3 5 0.7320
0 4 0 0.3268
0 4 1 0.4562
0 4 2 0.5925
0 4 3 0.4603
0 4 4 0.3265
0 4 5 0.4894
0 5 0 0.3295
0 5 1 0.4824
0 5 2 0.7211
0 5 3 0.7255
0 5 4 0.4936
0 5 5 0.3289
1 0 0 0.1094
1 0 1 0.1122
1 0 2 0.2099
1 0 3 0.3583
1 0 4 0.4680
1 0 5 0.4858
1 1 0 0.1110
1 1 1 0.0637
1 1 2 0.1155
1 1 3 0.2204
1 1 4 0.2945
1 1 5 0.3218
1 2 0 0.2009
1 2 1 0.1115
1 2 2 0.1141
1 2 3 0.2125
1 2 4 0.3666
1 2 5 0.4561
1 3 0 0.3688
1 3 1 0.2224
1 3 2 0.2155
1 3 3 0.2173
1 3 4 0.3737
1 3 5 0.5846
1 4 0 0.4545
1 4 1 0.3012
1 4 2 0.3585
1 4 3 0.3640
1 4 4 0.2981
1 4 5 0.4653
1 5 0 0.4875
1 5 1 0.3240
1 5 2 0.4689
1 5 3 0.5893
1 5 4 0.4622
1 5 5 0.3311
2 0 0 0.2204
2 0 1 0.2177
2 0 2 0.2186
2 0 3 0.3715
2 0 4 0.5923
2 0 5 0.7328
2 1 0 0.2096
2 1 1 0.1063
2 1 2 0.1129
2 1 3 0.2130
2 1 4 0.3653
2 1 5 0.4577
2 2 0 0.2213
2 2 1 0.1113
2 2 2 0.0708
2 2 3 0.1179
2 2 4 0.2210
2 2 5 0.2997
2 3 0 0.3650
2 3 1 0.2033
2 3 2 0.1155
2 3 3 0.1130
2 3 4 0.2093
2 3 5 0.3612
2 4 0 0.5889
2 4 1 0.3571
2 4 2 0.2250
2 4 3 0.2065
2 4 4 0.2158
2 4 5 0.3532
2 5 0 0.7280
2 5 1 0.4626
2 5 2 0.2926
2 5 3 0.3668
2 5 4 0.3616
2 5 5 0.3061
3 0 0 0.3071
3 0 1 0.3545
3 0 2 0.3553
3 0 3 0.2975
3 0 4 0.4556
3 0 5 0.7326
3 1 0 0.3525
3 1 1 0.2154
3 1 2 0.2065
3 1 3 0.2189
3 1 4 0.3585
3 1 5 0.5932
3 2 0 0.3569
3 2 1 0.2117
3 2 2 0.1147
3 2 3 0.1138
3 2 4 0.2066
3 2 5 0.3616
3 3 0 0.2999
3 3 1 0.2201
3 3 2 0.1118
3 3 3 0.0699
3 3 4 0.1122
3 3 5 0.2132
3 4 0 0.4603
3 4 1 0.3614
3 4 2 0.2075
3 4 3 0.1129
3 4 4 0.1122
3 4 5 0.2009
3 5 0 0.7402
3 5 1 0.5912
3 5 2 0.3530
3 5 3 0.2185
3 5 4 0.2083
3 5 5 0.2190
4 0 0 0.3423
4 0 1 0.4639
4 0 2 0.6001
4 0 3 0.4627
4 0 4 0.3264
4 0 5 0.4844
4 1 0 0.4489
4 1 1 0.2926
4 1 2 0.3685
4 1 3 0.3581
4 1 4 0.2958
4 1 5 0.4593
4 2 0 0.5919
4 2 1 0.3596
4 2 2 0.2207
4 2 3 0.2099
4 2 4 0.2145
4 2 5 0.3640
4 3 0 0.4554
4 3 1 0.3643
4 3 2 0.2154
4 3 3 0.1124
4 3 4 0.1133
4 3 5 0.2058
4 4 0 0.3307
4 4 1 0.3001
4 4 2 0.2161
4 4 3 0.1128
4 4 4 0.0692
4 4 5 0.1129
4 5 0 0.4842
4 5 1 0.4640
4 5 2 0.3540
4 5 3 0.2076
4 5 4 0.1093
4 5 5 0.1215
5 0 0 0.3386
5 0 1 0.4804
5 0 2 0.7271
5 0 3 0.7268
5 0 4 0.4927
5 0 5 0.3339
5 1 0 0.4819
5 1 1 0.3211
5 1 2 0.4656
5 1 3 0.5941
5 1 4 0.4635
5 1 5 0.3257
5 2 0 0.7223
5 2 1 0.4685
5 2 2 0.2922
5 2 3 0.3603
5 2 4 0.3639
5 2 5 0.3013
5 3 0 0.7299
5 3 1 0.5931
5 3 2 0.3630
5 3 3 0.2246
5 3 4 0.2060
5 3 5 0.2160
5 4 0 0.4929
5 4 1 0.4548
5 4 2 0.3628
5 4 3 0.2098
5 4 4 0.1110
5 4 5 0.1125
5 5 0 0.3430
5 5 1 0.3379
5 5 2 0.3000
5 5 3 0.2176
5 5 4 0.1185
5 5 5 0.0662

Power Curve

plot(Power,type = 'l',col = 'navyblue',lwd=2,ylab = 'power')
grid()

Large Sample Powers

\(n_1=25 ,n_2=20,n_3=22,n_4=30\)

When we are varying \(\theta_2\) in 1,2,3,4,5 We considered the Large sample tests having level 0.05

P=NULL
for(l in c(1,2,3,4,5))
{
h<-NULL
for(k in 1:1000)
{
  n1<-25
  s1<-rexp(n1,1)
  
  n2<-20
  s2<-rexp(n2,l)
  
  n3<-22
  s3<-rexp(n3,1)
  
  n4<-30
  s4<-rexp(n4,1)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  h[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}
I=0
for(k in 1:1000){
  if (h[k]>qchisq(0.05,3,lower.tail = F)){
    I=I+1
  }
}
P[l]= I/1000

}
l = 1:5
Power = P
send = cbind(l,Power)
colnames(send)[1] = 'theta_2'

kbl(send,format = 'html',caption = 'Large sample power Table') %>%
  kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
Large sample power Table
theta_2 Power
1 0.053
2 0.455
3 0.891
4 0.988
5 1.000

Power Curve

plot(1:5,Power,type = 'l',col = 'navyblue',lwd=2,ylab = 'power',xlab='theta 2')
grid()

\(n_1=25 ,n_2=20,n_3=22,n_4=30\)

When we are varying \(\theta_2\) in 1,2,3,4,5 & \(\theta_3\) in 1,2,3,4,5 & \(\theta_4\) in 1,2,3,4,5 We considered the small sample tests having level 0.05

P=NULL
for(o in c(1,2,3,4,5)){
for(m in c(1,2,3,4,5)){
for(l in c(1,2,3,4,5))
{
h<-NULL
for(k in 1:1000)
{
  n1<-25
  s1<-rexp(n1,1)
  
  n2<-20
  s2<-rexp(n2,l)
  
  n3<-22
  s3<-rexp(n3,m)
  
  n4<-30
  s4<-rexp(n4,o)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  h[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}
I=0
for(k in 1:1000){
  if (h[k]>qchisq(0.05,3,lower.tail = F)){
    I=I+1
  }
}
P[(o-1)*25+(m-1)*5+l]= I/1000

}
}
}
Power = P
o = rep(1:5,each = 25)
m = rep(0:5,each = 5,times=5)
l = rep(0:5, times = 25)

send = cbind(o,m, l, Power)
Warning in cbind(o, m, l, Power): number of rows of result is not a multiple of
vector length (arg 1)
kbl(send,format = 'html',caption = 'Small sample power',
    col.names = c('$\\theta_4$','$\\theta_3$','$\\theta_2$','Power'  )
    
    ) %>%
  kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
Small sample power
$\theta_4$ $\theta_3$ $\theta_2$ Power
1 0 0 0.046
1 0 1 0.469
1 0 2 0.884
1 0 3 0.988
1 0 4 1.000
1 1 5 0.490
1 1 0 0.629
1 1 1 0.913
1 1 2 0.982
1 1 3 0.997
1 2 4 0.905
1 2 5 0.906
1 2 0 0.974
1 2 1 0.993
1 2 2 0.999
1 3 3 0.992
1 3 4 0.987
1 3 5 0.995
1 3 0 0.998
1 3 1 1.000
1 4 2 1.000
1 4 3 0.999
1 4 4 0.998
1 4 5 1.000
1 4 0 1.000
2 5 1 0.579
2 5 2 0.677
2 5 3 0.889
2 5 4 0.979
2 5 5 0.998
2 0 0 0.678
2 0 1 0.523
2 0 2 0.755
2 0 3 0.903
2 0 4 0.967
2 1 5 0.912
2 1 0 0.770
2 1 1 0.869
2 1 2 0.929
2 1 3 0.973
2 2 4 0.976
2 2 5 0.935
2 2 0 0.956
2 2 1 0.986
2 2 2 0.991
2 3 3 0.998
2 3 4 0.985
2 3 5 0.980
2 3 0 0.992
2 3 1 0.997
3 4 2 0.962
3 4 3 0.934
3 4 4 0.975
3 4 5 0.990
3 4 0 0.998
3 5 1 0.924
3 5 2 0.830
3 5 3 0.882
3 5 4 0.949
3 5 5 0.982
3 0 0 0.972
3 0 1 0.880
3 0 2 0.891
3 0 3 0.941
3 0 4 0.979
3 1 5 0.993
3 1 0 0.947
3 1 1 0.943
3 1 2 0.961
3 1 3 0.989
3 2 4 0.997
3 2 5 0.986
3 2 0 0.976
3 2 1 0.991
3 2 2 0.994
4 3 3 0.997
4 3 4 0.992
4 3 5 0.994
4 3 0 0.998
4 3 1 0.999
4 4 2 0.993
4 4 3 0.951
4 4 4 0.955
4 4 5 0.982
4 4 0 0.995
4 5 1 0.996
4 5 2 0.963
4 5 3 0.962
4 5 4 0.983
4 5 5 0.990
4 0 0 1.000
4 0 1 0.979
4 0 2 0.980
4 0 3 0.984
4 0 4 0.992
4 1 5 1.000
4 1 0 0.991
4 1 1 0.992
4 1 2 0.991
4 1 3 0.993
5 2 4 1.000
5 2 5 0.998
5 2 0 0.998
5 2 1 1.000
5 2 2 1.000
5 3 3 0.999
5 3 4 0.989
5 3 5 0.995
5 3 0 0.991
5 3 1 0.998
5 4 2 0.999
5 4 3 0.993
5 4 4 0.978
5 4 5 0.995
5 4 0 0.999
5 5 1 1.000
5 5 2 0.997
5 5 3 0.994
5 5 4 0.989
5 5 5 0.994
5 0 0 1.000
5 0 1 1.000
5 0 2 0.996
5 0 3 0.992
5 0 4 0.992
1 1 5 0.046
1 1 0 0.469
1 1 1 0.884
1 1 2 0.988
1 1 3 1.000
1 2 4 0.490
1 2 5 0.629
1 2 0 0.913
1 2 1 0.982
1 2 2 0.997
1 3 3 0.905
1 3 4 0.906
1 3 5 0.974
1 3 0 0.993
1 3 1 0.999
1 4 2 0.992
1 4 3 0.987
1 4 4 0.995
1 4 5 0.998
1 4 0 1.000
1 5 1 1.000
1 5 2 0.999
1 5 3 0.998
1 5 4 1.000
1 5 5 1.000

Power Curve

plot(Power,type = 'l',col = 'navyblue',lwd=2,ylab = 'power')
grid()

Checking Consistency

If we can successfully show that the power tends to one as the sample size increases then we can claim that the test is Consistent

Power Curve with respect to sample size

P=NULL
for(n in c(20,30,40,50,60,70,80,90,100))
{
h<-NULL
for(k in 1:1000)
{
  
  s1<-rexp(n,1)
  
  
  s2<-rexp(n,2)
  
  
  s3<-rexp(n,1)
  
  
  s4<-rexp(n,1)
  
  N<-4*n
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n])/n
  r2<-mean(r[(n+1):(2*n)])
  r3<-mean(r[(2*n+1):(3*n)])
  r4<-mean(r[(3*n+1):(4*n)])
  
  
  h[k]<- ((12/((N+1)*N))*(n*r1**2+n*r2*r2+n*r3*r3+n*r4*r4)) - 3*(N+1)
}
I=0
for(k in 1:1000){
  if (h[k]>qchisq(0.05,3,lower.tail = F)){
    I=I+1
  }
}
P[(n/10)-1]= I/1000

}
n = c(20,30,40,50,60,70,80,90,100)
Power = P
send = cbind(n,Power)
colnames(send)=c('n','Power')

kbl(send,format = 'html',caption = 'Power Table') %>%
  kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
Power Table
n Power
20 0.431
30 0.645
40 0.775
50 0.888
60 0.918
70 0.970
80 0.983
90 0.995
100 0.998

Power Curve

plot(c(20,30,40,50,60,70,80,90,100),Power,type = 'l',col = 'navyblue',lwd=2,ylab = 'power',xlab='sample size')
grid()

Similarity with nonparametric two sample tests when k=2

We can show that when k=2 we shall have comparable results using Kruskal Wallis test as we would have had in two sample location nonparametric tests such as: - Wilcoxon Rank Sum Test - Mann Whitney Test

Power Comparison whit Mann-whitney Test

Consider a two sample location problem with \(n_1=25\) and \(n_2=20\)
We shall test the location parameters and calculate the power for \(\theta_1=0\) and \(\theta_2\) varying in 0,1,2,3,4,5,6 for both the Kruskal Wallis test and the Mann Whitney test.

library(kableExtra)
P=NULL
for(l in c(0,1,2,3,4,5))
{
h<-NULL
for(k in 1:10000)
{
  n1<-25
  s1<-rnorm(n1,0,10)
  
  n2<-20
  s2<-rnorm(n2,l,10)
  
  
  N<-n1+n2
  s<-c(s1,s2)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  
  
  
  h[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2)) - 3*(N+1)
}
I=0
for(k in 1:10000){
  if (h[k]>qchisq(0.05,1,lower.tail = F)){
    I=I+1
  }
}
P[l+1]= I/10000

}

Q<-NULL
for(l in c(0,1,2,3,4,5))
{
  U<-NULL
  S<-NULL
  for(m in 1:10000){
    n1=25
    n2=20
    s1<-rnorm(n1,0,10)
    s2<-rnorm(n2,l,10)
    k=1
    store<-NULL
    for(i in 1:n1){
      for(j in 1:n2){
        store[k]=s1[i]<s2[j]
        k=k+1
      }
    }
    U[m]<-sum(store)
    S[m]<- (U[m] - (n1*n2)/2)/sqrt(n1*n2*(n1+n2+1)/12)
  }
  I=0
  for(m in 1:10000){
    if (abs(S[m])>qnorm(0.025,lower.tail = F)){
      I=I+1
    }
  }
  Q[l+1]=I/10000
  
}

KW<- P
MW<- Q
l<-0:5
send<-cbind(l,KW,MW)
kbl(send,format = 'html',caption = 'Power Comparison',
    col.names = c('$\\theta_2$','KW power','MW Power'  ) )%>%
  kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
Power Comparison
$\theta_2$ KW power MW Power
0 0.0505 0.0514
1 0.0630 0.0642
2 0.0962 0.0998
3 0.1603 0.1609
4 0.2502 0.2425
5 0.3505 0.3637

Power Curve

plot(0:5,KW,type = 'l',col = 'navyblue',lwd=2,ylab = 'power',xlab='theta 2')
grid()
lines(0:5,MW,type='l',col='red',lwd=2)
legend('topleft',legend = c('Mann Whitney','Kruskal Wallis'),col = c('navyblue','red'),lty = c(1,1))

Comparison Of Kruskal Wallis Test with the parametric Test

In this section we shall consider comparing Kruskal Wallis test with its Parametric Counter Part.
Note that if it was known to us that the Samples are drawn from the populations which follows Normal distribution we would have opted for the Parametric test of Anova
We shall judge the validity of this decision by comparing the powers of the two tests

Power Comparison of Kruskal Wallis Test with Anova Test

library(kableExtra)
P=NULL
for(o in c(1,2,3,4,5)){
for(m in c(1,2,3,4,5)){
for(l in c(1,2,3,4,5))
{
h<-NULL
for(k in 1:1000)
{
  n1<-25
  s1<-rnorm(n1,0,10)
  
  n2<-20
  s2<-rnorm(n2,l,10)
  
  n3<-22
  s3<-rnorm(n3,m,10)
  
  n4<-30
  s4<-rnorm(n4,o,10)
  
  N<-n1+n2+n3+n4
  s<-c(s1,s2,s3,s4)
  r<-rank(s)
  r1<-sum(r[1:n1])/n1
  r2<-mean(r[(n1+1):(n1+n2)])
  r3<-mean(r[(n1+n2+1):(n1+n2+n3)])
  r4<-mean(r[(n1+n2+n3+1):(n1+n2+n3+n4)])
  
  
  h[k]<- ((12/((N+1)*N))*(n1*r1**2+n2*r2*r2+n3*r3*r3+n4*r4*r4)) - 3*(N+1)
}
I=0
for(k in 1:1000){
  if (h[k]>qchisq(0.05,3,lower.tail = F)){
    I=I+1
  }
}
P[(o-1)*25+(m-1)*5+l]= I/1000

}
}
}

Q<-NULL
for(o in c(0,1,2,3,4,5)){
  for(m in c(0,1,2,3,4,5)){
    for(l in c(0,1,2,3,4,5))
    {
      
f<-NULL
for(k in 1:1000)
{
  n1<-25
  s1<-rnorm(n1,0,10)
  
  n2<-20
  s2<-rnorm(n2,l,10)
  
  n3<-22
  s3<-rnorm(n3,m,10)
  
  n4<-30
  s4<-rnorm(n4,o,10)
  
  s<-c(s1,s2,s3,s4)
  N<-n1+n2+n3+n4
  SSW<-(n1-1)*var(s1)+(n2-1)*var(s2)+(n3-1)*var(s3)+(n4-1)*var(s4)
  TSS<-(N-1)*var(s)
  SSB<-TSS-SSW
  f[k]<-(SSB/3)/(SSW/N-4)
  
}
I=0
for(k in 1:1000){
  if (f[k]>qf(0.05,3,N-4,lower.tail = F)){
    I=I+1
  }
}
Q[(o-1)*25+(m-1)*5+l]= I/1000

    }
    
  }
}
KW= P
An=Q
o=rep(1:5,each=25)
m = rep(1:5,each =5,times=5)
l = rep(1:5, times = 25)

send = cbind(o,m, l,KW,An)

kbl(send,format = 'html',caption = 'Power Comparison',
    col.names = c('$\\theta_4$','$\\theta_3$','$\\theta_2$',' KW Power','Anova Power'  )
    
    ) %>%
  kable_styling(bootstrap_options = c('striped','hover','condensed','responsive'))
Power Comparison
$\theta_4$ $\theta_3$ $\theta_2$ KW Power Anova Power
1 1 1 0.052 0.083
1 1 2 0.077 0.107
1 1 3 0.103 0.140
1 1 4 0.166 0.219
1 1 5 0.278 0.104
1 2 1 0.074 0.105
1 2 2 0.098 0.105
1 2 3 0.111 0.152
1 2 4 0.153 0.229
1 2 5 0.253 0.157
1 3 1 0.118 0.142
1 3 2 0.135 0.148
1 3 3 0.149 0.206
1 3 4 0.202 0.227
1 3 5 0.294 0.251
1 4 1 0.187 0.237
1 4 2 0.191 0.211
1 4 3 0.196 0.264
1 4 4 0.262 0.317
1 4 5 0.300 0.134
1 5 1 0.299 0.115
1 5 2 0.257 0.152
1 5 3 0.276 0.200
1 5 4 0.347 0.285
1 5 5 0.371 0.106
2 1 1 0.070 0.122
2 1 2 0.083 0.127
2 1 3 0.100 0.142
2 1 4 0.154 0.219
2 1 5 0.264 0.135
2 2 1 0.089 0.131
2 2 2 0.082 0.124
2 2 3 0.087 0.135
2 2 4 0.179 0.200
2 2 5 0.233 0.190
2 3 1 0.106 0.173
2 3 2 0.133 0.153
2 3 3 0.118 0.175
2 3 4 0.167 0.227
2 3 5 0.248 0.262
2 4 1 0.170 0.215
2 4 2 0.159 0.238
2 4 3 0.169 0.230
2 4 4 0.208 0.280
2 4 5 0.269 0.231
2 5 1 0.269 0.208
2 5 2 0.238 0.205
2 5 3 0.230 0.225
2 5 4 0.282 0.285
2 5 5 0.336 0.191
3 1 1 0.126 0.139
3 1 2 0.140 0.163
3 1 3 0.155 0.186
3 1 4 0.221 0.233
3 1 5 0.255 0.218
3 2 1 0.115 0.182
3 2 2 0.129 0.147
3 2 3 0.143 0.158
3 2 4 0.173 0.210
3 2 5 0.228 0.231
3 3 1 0.143 0.210
3 3 2 0.160 0.204
3 3 3 0.147 0.200
3 3 4 0.173 0.229
3 3 5 0.237 0.318
3 4 1 0.185 0.234
3 4 2 0.166 0.244
3 4 3 0.191 0.235
3 4 4 0.234 0.283
3 4 5 0.248 0.361
3 5 1 0.298 0.330
3 5 2 0.268 0.320
3 5 3 0.248 0.332
3 5 4 0.291 0.359
3 5 5 0.308 0.299
4 1 1 0.216 0.237
4 1 2 0.218 0.268
4 1 3 0.249 0.238
4 1 4 0.245 0.310
4 1 5 0.328 0.310
4 2 1 0.202 0.252
4 2 2 0.195 0.248
4 2 3 0.192 0.227
4 2 4 0.229 0.297
4 2 5 0.297 0.315
4 3 1 0.232 0.270
4 3 2 0.194 0.245
4 3 3 0.191 0.233
4 3 4 0.214 0.253
4 3 5 0.285 0.397
4 4 1 0.242 0.327
4 4 2 0.238 0.271
4 4 3 0.228 0.260
4 4 4 0.240 0.323
4 4 5 0.302 0.475
4 5 1 0.331 0.441
4 5 2 0.306 0.449
4 5 3 0.295 0.477
4 5 4 0.272 0.504
4 5 5 0.327 0.449
5 1 1 0.333 0.401
5 1 2 0.317 0.368
5 1 3 0.326 0.373
5 1 4 0.351 0.420
5 1 5 0.413 0.427
5 2 1 0.314 0.379
5 2 2 0.282 0.359
5 2 3 0.289 0.331
5 2 4 0.306 0.403
5 2 5 0.350 0.466
5 3 1 0.328 0.377
5 3 2 0.270 0.336
5 3 3 0.283 0.337
5 3 4 0.308 0.361
5 3 5 0.329 0.476
5 4 1 0.352 0.443
5 4 2 0.309 0.374
5 4 3 0.284 0.363
5 4 4 0.304 0.405
5 4 5 0.338 0.526
5 5 1 0.424 0.461
5 5 2 0.355 0.429
5 5 3 0.354 0.408
5 5 4 0.342 0.436
5 5 5 0.358 0.447

Limitations Of Kruskal Wallis Test

Note that Kruskal Wallis test can only be used against the alternative hypothesis of type\[ H_1 :\,\,\,\,\theta_i\neq\,\theta_j\,\,\,\,\,\,\,for\,\,some\,\,\,i\,\neq\,j\] But this test is not applicable for testing the ordered alternatives i.e., \[H_1:\theta_1\le \theta_2\le...\le\theta_k\] or
\[H_1:\theta_1\ge \theta_2\ge...\ge \theta_k\] With strict inequality atleast once To solve this kind of problems we look upon to tests like JONCKHEERE TEST.

JOCNKHEERE TEST

consider the Testing Problem \[H_0:\theta_1= \theta_2=...= \theta_k\] vs \[H_1:\theta_1\le\theta_2\le...\le\theta_k\] with strict inequality atleast once. Jonckheere(1954) proposed a test which is based on the pairwise Mann Whitney Statistic.
This approach has the unique feature that the comparisons of samples i and j does not depends on the rest of the combined data.

Test Statistic

Define, \[W_{ij}\;= no.of(X_{vj}> X{ui})\;\;\;\;\;\;\;\;v=1(1)n_j;\,\,u=1(1)n_i\] Therefore \[ \mathcal{J} = \sum_{i=1}^{k}\sum_{j>i}W_{ij}\] We reject for large Values of \(\mathcal{J}\)

Expectation and Variance of Jonckheere Statistic under \(H_0\)

It can be shown that under \(H_0\) \[E(J) = N^{2} - \sum_{j=1}^{k}n_j^2/4\] \[ Var(J)= \dfrac{1}{72}[ N^{2}(2N+3)\,-\,\sum_{j=1}^{k}n_j^2(2n+3)]\]

Null Distribution

When the Sample sizes are sufficiently large \[(J-E(J))/(Var(J))\sim N(0,1)\;\;\;\;\;for\;\;large\;\;n \]