1 |
x1 |
|
y1 |
|
x2 |
|
y2 |
|
fr1 |
|
fr2 |
|
tr |
|
nboot |
|
pts |
|
plotit |
|
xout |
|
outfun |
|
... |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 | ##---- Should be DIRECTLY executable !! ----
##-- ==> Define data, use random,
##-- or do help(data=index) for the standard data sets.
## The function is currently defined as
function (x1, y1, x2, y2, fr1 = 1, fr2 = 1, tr = 0.2, nboot = 599,
pts = NA, plotit = TRUE, xout = FALSE, outfun = outpro, ...)
{
if (is.na(pts[1])) {
isub <- c(1:5)
test <- c(1:5)
m1 = elimna(cbind(x1, y1))
x1 = m1[, 1]
y1 = m1[, 2]
m1 = elimna(cbind(x2, y2))
x2 = m1[, 1]
y2 = m1[, 2]
xorder <- order(x1)
y1 <- y1[xorder]
x1 <- x1[xorder]
xorder <- order(x2)
y2 <- y2[xorder]
x2 <- x2[xorder]
n1 <- 1
n2 <- 1
vecn <- 1
for (i in 1:length(x1)) n1[i] <- length(y1[near(x1, x1[i],
fr1)])
for (i in 1:length(x1)) n2[i] <- length(y2[near(x2, x1[i],
fr2)])
for (i in 1:length(x1)) vecn[i] <- min(n1[i], n2[i])
sub <- c(1:length(x1))
isub[1] <- min(sub[vecn >= 12])
isub[5] <- max(sub[vecn >= 12])
isub[3] <- floor((isub[1] + isub[5])/2)
isub[2] <- floor((isub[1] + isub[3])/2)
isub[4] <- floor((isub[3] + isub[5])/2)
mat <- matrix(NA, 5, 8)
dimnames(mat) <- list(NULL, c("X", "n1", "n2", "DIF",
"TEST", "ci.low", "ci.hi", "p.value"))
gv1 <- vector("list")
for (i in 1:5) {
j <- i + 5
temp1 <- y1[near(x1, x1[isub[i]], fr1)]
temp2 <- y2[near(x2, x1[isub[i]], fr2)]
temp1 <- temp1[!is.na(temp1)]
temp2 <- temp2[!is.na(temp2)]
mat[i, 2] <- length(temp1)
mat[i, 3] <- length(temp2)
gv1[[i]] <- temp1
gv1[[j]] <- temp2
}
I1 <- diag(5)
I2 <- 0 - I1
con <- rbind(I1, I2)
test <- linconb(gv1, con = con, tr = tr, nboot = nboot)
for (i in 1:5) {
mat[i, 1] <- x1[isub[i]]
}
mat[, 4] <- test$psihat[, 2]
mat[, 5] <- test$test[, 2]
mat[, 6] <- test$psihat[, 3]
mat[, 7] <- test$psihat[, 4]
mat[, 8] <- test$test[, 4]
}
if (!is.na(pts[1])) {
n1 <- 1
n2 <- 1
vecn <- 1
for (i in 1:length(pts)) {
n1[i] <- length(y1[near(x1, pts[i], fr1)])
n2[i] <- length(y2[near(x2, pts[i], fr2)])
if (n1[i] <= 5)
paste("Warning, there are", n1[i], " points corresponding to the design point X=",
pts[i])
if (n2[i] <= 5)
paste("Warning, there are", n2[i], " points corresponding to the design point X=",
pts[i])
}
mat <- matrix(NA, length(pts), 9)
dimnames(mat) <- list(NULL, c("X", "n1", "n2", "DIF",
"TEST", "se", "ci.low", "ci.hi", "p.value"))
gv <- vector("list", 2 * length(pts))
for (i in 1:length(pts)) {
g1 <- y1[near(x1, pts[i], fr1)]
g2 <- y2[near(x2, pts[i], fr2)]
g1 <- g1[!is.na(g1)]
g2 <- g2[!is.na(g2)]
j <- i + length(pts)
gv[[i]] <- g1
gv[[j]] <- g2
}
I1 <- diag(length(pts))
I2 <- 0 - I1
con <- rbind(I1, I2)
test <- linconb(gv, con = con, tr = tr, nboot = nboot)
mat[, 1] <- pts
mat[, 2] <- n1
mat[, 3] <- n2
mat[, 4] <- test$psihat[, 2]
mat[, 5] <- test$test[, 2]
mat[, 6] <- test$test[, 3]
mat[, 7] <- test$psihat[, 3]
mat[, 8] <- test$psihat[, 4]
mat[, 9] <- test$test[, 4]
}
if (plotit) {
if (xout) {
flag <- outfun(x1, ...)$keep
x1 <- x1[flag]
y1 <- y1[flag]
flag <- outfun(x2, ...)$keep
x2 <- x2[flag]
y2 <- y2[flag]
}
runmean2g(x1, y1, x2, y2, fr = fr1, est = mean, tr = tr)
}
list(output = mat, crit = test$crit)
}
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.