R/smmcrit.R

smmcrit<-function(nuhat,C){
#
#  Determine the .95 quantile of the C-variate Studentized maximum
#  modulus distribution using linear interpolation on inverse
#  degrees of freedom
#  If C=1, this function returns the .975 quantile of Student's t
#  distribution.
#
if(C-round(C)!=0)stop("The number of contrasts, C, must be an  integer")
if(C>=29)stop("C must be less than or equal to 28")
if(C<=0)stop("C must be greater than or equal to 1")
if(nuhat<2)stop("The degrees of freedom must be greater than or equal to 2")
if(C==1)smmcrit<-qt(.975,nuhat)
if(C>=2){
C<-C-1
m1<-matrix(0,20,27)
m1[1,]<-c(5.57,6.34,6.89,7.31,7.65,7.93,8.17,8.83,8.57,
8.74,8.89,9.03,9.16,9.28,9.39,9.49,9.59, 9.68,
9.77,9.85,9.92,10.00,10.07,10.13,10.20,10.26,10.32)
m1[2,]<-c(3.96,4.43,4.76,5.02,5.23,5.41,5.56,5.69,5.81,
5.92,6.01,6.10,6.18,6.26,6.33,6.39,6.45,6.51,
6.57,6.62,6.67,6.71,6.76,6.80,6.84,6.88, 6.92)
m1[3,]<-c(3.38,3.74,4.01,4.20,4.37,4.50,4.62,4.72,4.82,
4.89,4.97,5.04,5.11,5.17,5.22,5.27,5.32, 5.37,
5.41,5.45,5.49,5.52,5.56,5.59,5.63,5.66,5.69)
m1[4,]<-c(3.09,3.39,3.62,3.79,3.93,4.04,4.14,4.23,4.31,
4.38,4.45,4.51,4.56,4.61,4.66,4.70,4.74,4.78,
4.82,4.85,4.89,4.92,4.95,4.98,5.00,5.03,5.06)
m1[5,]<-c(2.92,3.19,3.39,3.54,3.66,3.77,3.86,3.94,4.01,
4.07,4.13,4.18,4.23,4.28,4.32,4.36,4.39,4.43,
4.46,4.49,4.52,4.55,4.58,4.60,4.63,4.65,4.68)
m1[6,]<-c(2.80,3.06,3.24,3.38,3.49,3.59,3.67,3.74,3.80,
3.86,3.92,3.96,4.01,4.05,4.09,4.13,4.16,4.19,
4.22,4.25,4.28,4.31,4.33,4.35,4.38,4.39,4.42)
m1[7,]<-c(2.72,2.96,3.13,3.26,3.36,3.45,3.53,3.60,3.66,
3.71,3.76,3.81,3.85,3.89,3.93,3.96,3.99, 4.02,
4.05,4.08,4.10,4.13,4.15,4.18,4.19,4.22,4.24)
m1[8,]<-c(2.66,2.89,3.05,3.17,3.27,3.36,3.43,3.49,3.55,
3.60,3.65,3.69,3.73,3.77,3.80,3.84,3.87,3.89,
3.92,3.95,3.97,3.99,4.02,4.04,4.06,4.08,4.09)
m1[9,]<-c(2.61,2.83,2.98,3.10,3.19,3.28,3.35,3.41,3.47,
3.52,3.56,3.60,3.64,3.68,3.71,3.74,3.77,3.79,
3.82,3.85,3.87,3.89,3.91,3.94,3.95, 3.97,3.99)
m1[10,]<-c(2.57,2.78,2.93,3.05,3.14,3.22,3.29,3.35,3.40,
3.45,3.49,3.53,3.57,3.60,3.63,3.66,3.69,3.72,
3.74,3.77,3.79,3.81,3.83,3.85,3.87,3.89,3.91)
m1[11,]<-c(2.54,2.75,2.89,3.01,3.09,3.17,3.24,3.29,3.35,
3.39,3.43,3.47,3.51,3.54,3.57,3.60,3.63,3.65,
3.68,3.70,3.72,3.74,3.76,3.78,3.80,3.82,3.83)
m1[12,]<-c(2.49,2.69,2.83,2.94,3.02,3.09,3.16,3.21,3.26,
3.30,3.34,3.38,3.41,3.45,3.48,3.50,3.53,3.55,
3.58,3.59,3.62,3.64,3.66,3.68,3.69,3.71,3.73)
m1[13,]<-c(2.46,2.65,2.78,2.89,2.97,3.04,3.09,3.15,3.19,
3.24,3.28,3.31,3.35,3.38,3.40,3.43,3.46,3.48,
3.50,3.52,3.54,3.56,3.58,3.59,3.61,3.63,3.64)
m1[14,]<-c(2.43,2.62,2.75,2.85,2.93,2.99,3.05,3.11,3.15,
3.19,3.23,3.26,3.29,3.32,3.35,3.38,3.40,3.42,
3.44,3.46,3.48,3.50,3.52,3.54,3.55,3.57,3.58)
m1[15,]<-c(2.41,2.59,2.72,2.82,2.89,2.96,3.02,3.07,3.11,
3.15,3.19,3.22,3.25,3.28,3.31,3.33,3.36,3.38,
3.39,3.42,3.44,3.46,3.47,3.49,3.50,3.52,3.53)
m1[16,]<-c(2.38,2.56,2.68,2.77,2.85,2.91,2.97,3.02,3.06,
3.09,3.13,3.16,3.19,3.22,3.25,3.27,3.29,3.31,
3.33,3.35,3.37,3.39,3.40,3.42,3.43,3.45,3.46)
m1[17,]<-c(2.35,2.52,2.64,2.73,2.80,2.87,2.92,2.96,3.01,
3.04,3.07,3.11,3.13,3.16,3.18,3.21,3.23,3.25,
3.27,3.29,3.30,3.32,3.33,3.35,3.36,3.37,3.39)
m1[18,]<-c(2.32,2.49,2.60,2.69,2.76,2.82,2.87,2.91,2.95,
2.99,3.02,3.05,3.08,3.09,3.12,3.14,3.17, 3.18,
3.20,3.22,3.24,3.25,3.27,3.28,3.29,3.31,3.32)
m1[19,]<-c(2.29,2.45,2.56,2.65,2.72,2.77,2.82,2.86,2.90,
2.93,2.96,2.99,3.02,3.04,3.06,3.08,3.10, 3.12,
3.14,3.16,3.17,3.19,3.20,3.21,3.23,3.24,3.25)
m1[20,]<-c(2.24,2.39,2.49,2.57,2.63,2.68,2.73,2.77,2.79,
2.83,2.86,2.88,2.91,2.93,2.95,2.97,2.98, 3.01,
3.02,3.03,3.04,3.06,3.07,3.08,3.09,3.11,3.12)
if(nuhat>=200)smmcrit<-m1[20,C]
if(nuhat<200){
nu<-c(2,3,4,5,6,7,8,9,10,11,12,14,16,18,20,24,30,40,60,200)
temp<-abs(nu-nuhat)
find<-order(temp)
if(temp[find[1]]==0)smmcrit<-m1[find[1],C]
if(temp[find[1]]!=0){
if(nuhat>nu[find[1]]){
smmcrit<-m1[find[1],C]-
(1/nu[find[1]]-1/nuhat)*(m1[find[1],C]-m1[find[1]+1,C])/
(1/nu[find[1]]-1/nu[find[1]+1])
}
if(nuhat<nu[find[1]]){
smmcrit<-m1[find[1]-1,C]-
(1/nu[find[1]-1]-1/nuhat)*(m1[find[1]-1,C]-m1[find[1],C])/
(1/nu[find[1]-1]-1/nu[find[1]])
}
}}
}
smmcrit
}

Try the WRS2 package in your browser

Any scripts or data that you put into this service are public.

WRS2 documentation built on May 2, 2019, 4:46 p.m.