mulwmwcrit:

Usage Arguments Examples

Usage

1
mulwmwcrit(mm1, mm2, plotit = TRUE, cop = 3, iter = 1000, alpha = 0.05, SEED = NA, pr = FALSE)

Arguments

mm1
mm2
plotit
cop
iter
alpha
SEED
pr

Examples

 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
##---- 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 (mm1, mm2, plotit = TRUE, cop = 3, iter = 1000, alpha = 0.05, 
    SEED = NA, pr = FALSE) 
{
    if (!is.matrix(mm1)) 
        stop("Data are assumed to be stored in a matrix having two or more columns. For univariate data, use the function outbox or out")
    if (is.na(SEED)) 
        set.seed(2)
    if (!is.na(SEED)) 
        set.seed(SEED)
    val <- NA
    n1 <- nrow(mm1)
    n2 <- nrow(mm2)
    for (it in 1:iter) {
        ivec1 <- sample(c(1:n1), replace = TRUE)
        ivec2 <- sample(c(1:n2), replace = TRUE)
        m1 <- mm1[ivec1, ]
        m2 <- mm2[ivec2, ]
        if (cop == 1) {
            if (ncol(m1) > 2) {
                center1 <- dmean(m1, tr = 0.5)
                center2 <- dmean(m2, tr = 0.5)
            }
            if (ncol(m1) == 2) {
                tempd <- NA
                for (i in 1:nrow(m1)) tempd[i] <- depth(m1[i, 
                  1], m1[i, 2], m1)
                mdep <- max(tempd)
                flag <- (tempd == mdep)
                if (sum(flag) == 1) 
                  center1 <- m1[flag, ]
                if (sum(flag) > 1) 
                  center1 <- apply(m1[flag, ], 2, mean)
                for (i in 1:nrow(m2)) tempd[i] <- depth(m2[i, 
                  1], m2[i, 2], m2)
                mdep <- max(tempd)
                flag <- (tempd == mdep)
                if (sum(flag) == 1) 
                  center2 <- m2[flag, ]
                if (sum(flag) > 1) 
                  center2 <- apply(m2[flag, ], 2, mean)
            }
        }
        if (cop == 2) {
            center1 <- cov.mcd(m1)$center
            center2 <- cov.mcd(m2)$center
        }
        if (cop == 3) {
            center1 <- apply(m1, 2, median)
            center2 <- apply(m2, 2, median)
        }
        center <- (center1 + center2)/2
        B <- center1 - center2
        if (sum(center1^2) > sum(center2^2)) 
            B <- (0 - 1) * B
        BB <- B^2
        bot <- sum(BB)
        disx <- NA
        disy <- NA
        if (bot != 0) {
            for (j in 1:nrow(m1)) {
                AX <- m1[j, ] - center
                tempx <- sum(AX * B) * B/bot
                disx[j] <- sign(sum(AX * B)) * sqrt(sum(tempx^2))
            }
            for (j in 1:nrow(m2)) {
                AY <- m2[j, ] - center
                tempy <- sum(AY * B) * B/bot
                disy[j] <- sign(sum(AY * B)) * sqrt(sum(tempy^2))
            }
        }
        m <- outer(disx, disy, FUN = "-")
        m <- sign(m)
        val[it] <- (1 - mean(m))/2
        if (bot == 0) 
            val[it] <- 0.5
        if (pr) 
            print(paste("Iteration ", it, " of ", iter, " complete"))
    }
    val <- sort(val)
    low <- round(alpha * iter/2) + 1
    up <- iter - low
    crit <- NA
    crit[1] <- val[low]
    crit[2] <- val[up]
    crit
  }

musto101/wilcox_R documentation built on May 23, 2019, 10:52 a.m.