title: R code for generalized (extreme Studentized deviate) ESD outlier test
body: <p>Here is the test I'm interested in:
  <a href="http://www.itl.nist.gov/div898/handbook/eda/section3/eda35h3.htm" rel="nofollow">http://www.itl.nist.gov/div898/handbook/eda/section3/eda35h3.htm</a></p>
  
  <p>How can I adapt this code into a function that accepts a vector of numeric values and returns a logical vector specifying which data points to remove?</p>
  
  <p>I have attempted to do this below, but I'm getting stuck because when I sort the vector to return, it doesn't line up with the input vector data.</p>
  
  <pre><code># input data
  y = c(-0.25, 0.68, 0.94, 1.15, 1.20, 1.26, 1.26,
        1.34, 1.38, 1.43, 1.49, 1.49, 1.55, 1.56,
        1.58, 1.65, 1.69, 1.70, 1.76, 1.77, 1.81,
        1.91, 1.94, 1.96, 1.99, 2.06, 2.09, 2.10,
        2.14, 2.15, 2.23, 2.24, 2.26, 2.35, 2.37,
        2.40, 2.47, 2.54, 2.62, 2.64, 2.90, 2.92,
        2.92, 2.93, 3.21, 3.26, 3.30, 3.59, 3.68,
        4.30, 4.64, 5.34, 5.42, 6.01)
  
  ## Generate normal probability plot.
  qqnorm(y)
  
  removeoutliers = function(dfinputcol) {
  
    y = as.vector(dfinputcol)
  
    ## Create function to compute the test statistic.
    rval = function(y){
      ares = abs(y - mean(y))/sd(y)
      df = data.frame(y, ares)
      r = max(df$ares)
      list(r, df)}
  
    ## Define values and vectors.
    n = length(y)
    alpha = 0.05
    lam = c(1:10)
    R = c(1:10)
  
    ## Compute test statistic until r=10 values have been
    ## removed from the sample.
    for (i in 1:10){
  
      if(i==1){
        rt = rval(y)
        R[i] = unlist(rt[1])
        df = data.frame(rt[2])
        newdf = df[df$ares!=max(df$ares),]}
  
      else if(i!=1){
        rt = rval(newdf$y)
        R[i] = unlist(rt[1])
        df = data.frame(rt[2])
        newdf = df[df$ares!=max(df$ares),]}
  
      ## Compute critical value.
      p = 1 - alpha/(2*(n-i+1))
      t = qt(p,(n-i-1))
      lam[i] = t*(n-i) / sqrt((n-i-1+t**2)*(n-i+1))
  
    }
    ## Print results.
    newdf = data.frame(c(1:10),R,lam)
    names(newdf)=c("Outliers","TestStat.", "CriticalVal.")
  
    # determine how many outliers to remove
    toremove = max(newdf$Outliers[newdf$TestStat. &gt; newdf$CriticalVal.])
  
    # create vector of same size as input vector
    logicalvectorTifshouldremove = logical(length=length(y))
  
    # but how to determine which outliers to remove?
    # set largest data points as outliers to remove.. but could be the smallest in some data sets..
    logicalvectorTifshouldremove = replace(logicalvectorTifshouldremove, tail(sort(y), toremove), TRUE)
  
    return (logicalvectorTifshouldremove)
  }
  
  # this should have 3 data points set to TRUE .. but it has 2 and they aren't the correct ones
  output = removeoutliers(y)
  length(output[output==T])
  </code></pre>