mcskew:

In musto101/wilcox_R: What the package does (short line)

Usage

mcskew(z)

Arguments

z

Examples

##---- 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 (z) 
{
    n = length(z)
    y1 = 0
    y2 = 0
    left = 0
    right = 0
    q = 0
    p = 0
    eps = 1e-13
    z = -z
    xmed = pull(z, n, floor(n/2) + 1)
    if (n%%2 == 0) {
        xmed = (xmed + pull(z, n, floor(n/2)))/2
    }
    z = z - xmed
    y = -sort(z)
    y1 = y[y > -eps]
    y2 = y[y <= eps]
    h1 = length(y1)
    h2 = length(y2)
    left[1:h2] = 1
    right[1:h2] = h1
    nl = 0
    nr = h1 * h2
    knew = floor(nr/2) + 1
    IsFound = 0
    while ((nr - nl > n) & (IsFound == 0)) {
        weight = 0
        work = 0
        j = 1
        for (i in 1:h2) {
            if (left[i] <= right[i]) {
                weight[j] = right[i] - left[i] + 1
                k = left[i] + floor(weight[j]/2)
                work[j] = calwork(y1[k], y2[i], k, i, h1 + 1, 
                  eps)
                j = j + 1
            }
        }
        trial = whimed(work, weight, j - 1)
        j = 1
        for (i in h2:1) {
            while ((j <= h1) & (calwork(y1[min(j, h1)], y2[i], 
                j, i, h1 + 1, eps) > trial)) {
                j = j + 1
            }
            p[i] = j - 1
        }
        j = h1
        for (i in 1:h2) {
            while ((j >= 1) & (calwork(y1[max(j, 1)], y2[i], 
                j, i, h1 + 1, eps) < trial)) {
                j = j - 1
            }
            q[i] = j + 1
        }
        sump = sum(p[1:h2])
        sumq = sum(q[1:h2]) - h2
        if (knew <= sump) {
            right[1:h2] = p[1:h2]
            nr = sump
        }
        else {
            if (knew > sumq) {
                left[1:h2] = q[1:h2]
                nl = sumq
            }
            else {
                medc = trial
                IsFound = 1
            }
        }
    }
    if (IsFound == 0) {
        work = 0
        j = 1
        for (i in 1:h2) {
            if (left[i] <= right[i]) {
                for (jj in left[i]:right[i]) {
                  work[j] = 0 - calwork(y1[jj], y2[i], jj, i, 
                    h1 + 1, eps)
                  j = j + 1
                }
            }
        }
        medc = 0 - pull(work, j - 1, knew - nl)
    }
    medc
  }

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