num <- 50
res <- array(NA,dim=c(num,3))
q1 <- function(data,n,p,i,f) {
np <- n*p;
i <<- floor(np)
f <<- np - i
qvalue <- (1-f)*data[i] + f*data[i+1]
}
q2 <- function(data,n,p,i,f) {
np <- (n+1)*p
i <<- floor(np)
f <<- np - i
qvalue <- (1-f)*data[i] + f*data[i+1]
}
q3 <- function(data,n,p,i,f) {
np <- n*p
i <<- floor(np)
f <<- np - i
if (f==0) {
qvalue <- data[i]
} else {
qvalue <- data[i+1]
}
}
q4 <- function(data,n,p,i,f) {
np <- n*p
i <<- floor(np)
f <<- np - i
if (f==0) {
qvalue <- (data[i]+data[i+1])/2
} else {
qvalue <- data[i+1]
}
}
q5 <- function(data,n,p,i,f) {
np <- (n-1)*p
i <<- floor(np)
f <<- np - i
if (f==0) {
qvalue <- data[i+1]
} else {
qvalue <- data[i+1] + f*(data[i+2]-data[i+1])
}
}
q6 <- function(data,n,p,i,f) {
np <- n*p+0.5
i <<- floor(np)
f <<- np - i
qvalue <- data[i]
}
q7 <- function(data,n,p,i,f) {
np <- (n+1)*p
i <<- floor(np)
f <<- np - i
if (f==0) {
qvalue <- data[i]
} else {
qvalue <- f*data[i] + (1-f)*data[i+1]
}
}
q8 <- function(data,n,p,i,f) {
np <- (n+1)*p
i <<- floor(np)
f <<- np - i
if (f==0) {
qvalue <- data[i]
} else {
if (f == 0.5) {
qvalue <- (data[i]+data[i+1])/2
} else {
if (f < 0.5) {
qvalue <- data[i]
} else {
qvalue <- data[i+1]
}
}
}
}
iqd <- function(x,def) {
x <-sort(x[!is.na(x)])
n<-length(x)
if (def==1) {
qvalue1 <- q1(x,n,0.25,i,f)
qvalue3 <- q1(x,n,0.75,i,f)
}
if (def==2) {
qvalue1 <- q2(x,n,0.25,i,f)
qvalue3 <- q2(x,n,0.75,i,f)
}
if (def==3) {
qvalue1 <- q3(x,n,0.25,i,f)
qvalue3 <- q3(x,n,0.75,i,f)
}
if (def==4) {
qvalue1 <- q4(x,n,0.25,i,f)
qvalue3 <- q4(x,n,0.75,i,f)
}
if (def==5) {
qvalue1 <- q5(x,n,0.25,i,f)
qvalue3 <- q5(x,n,0.75,i,f)
}
if (def==6) {
qvalue1 <- q6(x,n,0.25,i,f)
qvalue3 <- q6(x,n,0.75,i,f)
}
if (def==7) {
qvalue1 <- q7(x,n,0.25,i,f)
qvalue3 <- q7(x,n,0.75,i,f)
}
if (def==8) {
qvalue1 <- q8(x,n,0.25,i,f)
qvalue3 <- q8(x,n,0.75,i,f)
}
iqdiff <- qvalue3 - qvalue1
return(c(iqdiff,iqdiff/2,iqdiff/(qvalue3 + qvalue1)))
}
range <- max(x) - min(x)
lx <- length(x)
biasf <- (lx-1)/lx
varx <- var(x)
bvarx <- varx*biasf
sdx <- sqrt(varx)
mx <- mean(x)
bsdx <- sqrt(bvarx)
x2 <- x*x
mse0 <- sum(x2)/lx
xmm <- x-mx
xmm2 <- xmm*xmm
msem <- sum(xmm2)/lx
axmm <- abs(x - mx)
medx <- median(x)
axmmed <- abs(x - medx)
xmmed <- x - medx
xmmed2 <- xmmed*xmmed
msemed <- sum(xmmed2)/lx
qarr <- array(NA,dim=c(8,3))
for (j in 1:8) {
qarr[j,] <- iqd(x,j)
}
sdpo <- 0
adpo <- 0
for (i in 1:(lx-1)) {
for (j in (i+1):lx) {
ldi <- x[i]-x[j]
aldi <- abs(ldi)
sdpo = sdpo + ldi * ldi
adpo = adpo + aldi
}
}
denom <- (lx*(lx-1)/2)
sdpo = sdpo / denom
adpo = adpo / denom
gmd <- 0
for (i in 1:lx) {
for (j in 1:lx) {
ldi <- abs(x[i]-x[j])
gmd = gmd + ldi
}
}
gmd <- gmd / (lx*(lx-1))
sumx <- sum(x)
pk <- x / sumx
ck <- cumsum(pk)
dk <- array(NA,dim=lx)
for (i in 1:lx) {
if (ck[i] <= 0.5) dk[i] <- ck[i] else dk[i] <- 1 - ck[i]
}
bigd <- sum(dk) * 2 / (lx-1)
iod <- 1 - sum(pk*pk)
res[1,] <- c('Absolute range','absolute.htm', range)
res[2,] <- c('Relative range (unbiased)','relative.htm', range/sd(x))
res[3,] <- c('Relative range (biased)','relative.htm', range/sqrt(varx*biasf))
res[4,] <- c('Variance (unbiased)','unbiased.htm', varx)
res[5,] <- c('Variance (biased)','biased.htm', bvarx)
res[6,] <- c('Standard Deviation (unbiased)','unbiased1.htm', sdx)
res[7,] <- c('Standard Deviation (biased)','biased1.htm', bsdx)
res[8,] <- c('Coefficient of Variation (unbiased)','variation.htm', sdx/mx)
res[9,] <- c('Coefficient of Variation (biased)','variation.htm', bsdx/mx)
res[10,] <- c('Mean Squared Error (MSE versus 0)','mse.htm', mse0)
res[11,] <- c('Mean Squared Error (MSE versus Mean)','mse.htm', msem)
res[12,] <- c('Mean Absolute Deviation from Mean (MAD Mean)', 'mean2.htm', sum(axmm)/lx)
res[13,] <- c('Mean Absolute Deviation from Median (MAD Median)', 'median1.htm', sum(axmmed)/lx)
res[14,] <- c('Median Absolute Deviation from Mean', 'mean3.htm', median(axmm))
res[15,] <- c('Median Absolute Deviation from Median', 'median2.htm', median(axmmed))
res[16,] <- c('Mean Squared Deviation from Mean', 'mean1.htm', msem)
res[17,] <- c('Mean Squared Deviation from Median', 'median.htm', msemed)
load(file='createtable')
mylink1 <- 'Interquartile Difference'
mylink2 <- paste(mylink1,'(Weighted Average at Xnp)',sep=' ')
res[18,] <- c('', mylink2, qarr[1,1])
mylink2 <- paste(mylink1,'(Weighted Average at X(n+1)p)',sep=' ')
res[19,] <- c('', mylink2, qarr[2,1])
mylink2 <- paste(mylink1,'(Empirical Distribution Function)',sep=' ')
res[20,] <- c('', mylink2, qarr[3,1])
mylink2 <- paste(mylink1,'(Empirical Distribution Function - Averaging)',sep=' ')
res[21,] <- c('', mylink2, qarr[4,1])
mylink2 <- paste(mylink1,'(Empirical Distribution Function - Interpolation)',sep=' ')
res[22,] <- c('', mylink2, qarr[5,1])
mylink2 <- paste(mylink1,'(Closest Observation)',sep=' ')
res[23,] <- c('', mylink2, qarr[6,1])
mylink2 <- paste(mylink1,'(True Basic - Statistics Graphics Toolkit)',sep=' ')
res[24,] <- c('', mylink2, qarr[7,1])
mylink2 <- paste(mylink1,'(MS Excel (old versions))',sep=' ')
res[25,] <- c('', mylink2, qarr[8,1])
mylink1 <- 'Semi Interquartile Difference'
mylink2 <- paste(mylink1,'(Weighted Average at Xnp)',sep=' ')
res[26,] <- c('', mylink2, qarr[1,2])
mylink2 <- paste(mylink1,'(Weighted Average at X(n+1)p)',sep=' ')
res[27,] <- c('', mylink2, qarr[2,2])
mylink2 <- paste(mylink1,'(Empirical Distribution Function)',sep=' ')
res[28,] <- c('', mylink2, qarr[3,2])
mylink2 <- paste(mylink1,'(Empirical Distribution Function - Averaging)',sep=' ')
res[29,] <- c('', mylink2, qarr[4,2])
mylink2 <- paste(mylink1,'(Empirical Distribution Function - Interpolation)',sep=' ')
res[30,] <- c('', mylink2, qarr[5,2])
mylink2 <- paste(mylink1,'(Closest Observation)',sep=' ')
res[31,] <- c('', mylink2, qarr[6,2])
mylink2 <- paste(mylink1,'(True Basic - Statistics Graphics Toolkit)',sep=' ')
res[32,] <- c('', mylink2, qarr[7,2])
mylink2 <- paste(mylink1,'(MS Excel (old versions))',sep=' ')
res[33,] <- c('', mylink2, qarr[8,2])
mylink1 <- 'Coefficient of Quartile Variation'
mylink2 <- paste(mylink1,'(Weighted Average at Xnp)',sep=' ')
res[34,] <- c('', mylink2, qarr[1,3])
mylink2 <- paste(mylink1,'(Weighted Average at X(n+1)p)',sep=' ')
res[35,] <- c('', mylink2, qarr[2,3])
mylink2 <- paste(mylink1,'(Empirical Distribution Function)',sep=' ')
res[36,] <- c('', mylink2, qarr[3,3])
mylink2 <- paste(mylink1,'(Empirical Distribution Function - Averaging)',sep=' ')
res[37,] <- c('', mylink2, qarr[4,3])
mylink2 <- paste(mylink1,'(Empirical Distribution Function - Interpolation)',sep=' ')
res[38,] <- c('', mylink2, qarr[5,3])
mylink2 <- paste(mylink1,'(Closest Observation)',sep=' ')
res[39,] <- c('', mylink2, qarr[6,3])
mylink2 <- paste(mylink1,'(True Basic - Statistics Graphics Toolkit)',sep=' ')
res[40,] <- c('', mylink2, qarr[7,3])
mylink2 <- paste(mylink1,'(MS Excel (old versions))',sep=' ')
res[41,] <- c('', mylink2, qarr[8,3])
res[42,] <- c('Number of all Pairs of Observations', 'pair_numbers.htm', lx*(lx-1)/2)
res[43,] <- c('Squared Differences between all Pairs of Observations', 'squared_differences.htm', sdpo)
res[44,] <- c('Mean Absolute Differences between all Pairs of Observations', 'mean_abs_differences.htm', adpo)
res[45,] <- c('Gini Mean Difference', 'gini_mean_difference.htm', gmd)
res[46,] <- c('Leik Measure of Dispersion', 'leiks_d.htm', bigd)
res[47,] <- c('Index of Diversity', 'diversity.htm', iod)
res[48,] <- c('Index of Qualitative Variation', 'qualitative_variation.htm', iod*lx/(lx-1))
res[49,] <- c('Coefficient of Dispersion', 'dispersion.htm', sum(axmm)/lx/medx)
res[50,] <- c('Observations', '', lx)
print(res)
a<-table.start()
a<-table.row.start(a)
a<-table.element(a,'Variability - Ungrouped Data',2,TRUE)
a<-table.row.end(a)
for (i in 1:num) {
a<-table.row.start(a)
if (res[i,1] != '') {
a<-table.element(a,res[i,1],header=TRUE)
} else {
a<-table.element(a,res[i,2],header=TRUE)
}
a<-table.element(a,signif(as.numeric(res[i,3],6)))
a<-table.row.end(a)
}
a<-table.end(a)
table.save(a,file='mytable.tab')
|