version 1.0

DESCRIPTION

@@ -0,0 +1,18 @@
+Package: ExcessMass
+Type: Package
+Title: Excess Mass Calculation and Plots
+Version: 1.0
+Date: 2017-05-10
+Author: Marc-Daniel Mildenberger
+Maintainer: Marc-Daniel Mildenberger <mildenberger.stat@web.de>
+Description: Implementation of a function which calculates the empirical excess mass 
+	for given \eqn{\lambda} and given maximal number of modes (excessm()). Offering 
+	powerful plot features to visualize empirical excess mass (exmplot()). This 
+	includes the possibility of drawing several plots (with different maximal 
+	number of modes / cut off values) in a single graph.
+License: LGPL
+Suggests: MASS
+Packaged: 2017-05-16 15:06:53 UTC; d91137
+NeedsCompilation: no
+Repository: CRAN
+Date/Publication: 2017-05-16 15:16:45 UTC

MD5

@@ -0,0 +1,17 @@
+aa0d7408fd8947d9f2a2900a056ff937 *DESCRIPTION
+422c3e0fa5253f35a30974603b237a05 *NAMESPACE
+a4d418804f9d4ecfeaf82c69e6547885 *R/colorizing.R
+9aee9e186c6cb52cbab95659708a8a21 *R/excessm.R
+939ffefc019202a4cb03f00b05003545 *R/exmplot.R
+374b60dc8ae4e6b93f8266d6fce4b1e6 *R/exmsilhouette.R
+d9a92773c344a157ce28ea08aa53f1ee *R/lambdaweight.R
+90014b446bf071e0dca5e814e334607a *R/localmax.R
+ddbc9fde5ba3feb40cec14141bb45fb1 *R/mexmsilhouette.R
+ce090839325e8cdf28fde021fd1d9db4 *R/searchMaxLambda.R
+516b0262616bc860bd3aab469c52cd94 *man/ExcessMass-package.Rd
+bf21a73daf6831b3da85752292b71569 *man/excessm.Rd
+ea3b13dbb040747ee3d389e1b621276f *man/exmplot.Rd
+29ff7deed87b38d2bfb0f721de65c99f *man/exmsilhouette.Rd
+6bd9b3bb6851afd00461b7f82ea52e59 *man/mexmsilhouette.Rd
+4bc2b458a21838b564688592d1e8caf6 *man/searchMaxLambda.Rd
+3452bf4960e44b1410698822c0a200c1 *tests/ex1.R

NAMESPACE

@@ -0,0 +1,2 @@
+importFrom("graphics", "lines", "mtext", "par", "plot", "title")
+exportPattern("^[[:alpha:]]+")

R/colorizing.R

@@ -0,0 +1,16 @@
+colorizing <- function(excm){
+  color="black"
+    if(excm > .75)
+      {color="mediumpurple4"
+      }else{if(excm > .5)
+        {color="green3"
+        }else{if(excm > .25)
+          {color="turquoise2"
+          }else{if(excm > .05)
+            {color="blue"
+            }
+          }
+        }
+      }   
+  return(color)
+}

R/excessm.R

@@ -0,0 +1,149 @@
+excessm <-function(x, lambda, M=1, UpToM=FALSE){
+if(all(is.na(x))){stop('Input vector is empty.')}else{
+  x <- sort(x)
+  L <- lambdaweight(x,lambda)
+  localmaxres <- localmax(L)
+  if(localmaxres[[1]]!=1 & localmaxres[[2]]!=length(x))
+  {
+    m=matrix(c(1,localmaxres[[1]],1,localmaxres[[1]],localmaxres[[2]],0,localmaxres[[2]],
+               length(x),1),3,3, byrow = TRUE)
+  }else{if((localmaxres[[1]]==1) & (localmaxres[[2]]==length(x)))
+  {
+    m=matrix(c(1,length(x),0),1,3, byrow = TRUE)
+  }else{
+    if(localmaxres[[1]]==localmaxres[[2]])
+    {stop('Interval which maximizes excess mass is empty set')
+    }else{if(localmaxres[[1]]==1){
+      m=matrix(c(1,localmaxres[[2]],0,localmaxres[[2]],length(x),1),2,3, byrow = TRUE)
+    }else{ 
+      m=matrix(c(1,localmaxres[[1]],1,localmaxres[[1]],length(x),0),2,3, byrow = TRUE)
+    }     
+    }
+  }         
+  }
+  if(UpToM==TRUE)
+  {
+    u <- matrix(nrow=(dim(m)[1]-sum(m[,3])), ncol=4,dimnames=list(1:(dim(m)[1]-sum(m[,3])),
+                                                                  c("Start Index","End Index","Start Point","End Point")))
+    for(i in 1:(dim(m)[1]))
+    {
+      if(m[i,3]==0)
+      {
+        u[1,1] <- m[i,1]
+        u[1,2] <- m[i,2]
+        u[1,3] <- x[m[i,1]]
+        u[1,4] <- x[m[i,2]]
+      }
+    }
+    r <- list(u)
+  }
+  variation <- localmaxres[[3]]
+  if(M>1)
+  {
+    for(i in 2:M)
+    {
+      n=matrix(nrow=(dim(m)[1]),ncol=3)
+      for(j in (1:(dim(m)[1]))) 
+      {
+        localmaxres <- localmax(L,m[[j,1]],m[[j,2]],m[[j,3]])
+        n[j,1] <- localmaxres[[1]]
+        n[j,2] <- localmaxres[[2]]
+        n[j,3] <- localmaxres[[3]]
+      }
+      maxemi <- which.max(abs(n[,3])) 
+
+      if(n[maxemi,3] == 0){warning('Number of intervals is smaller than M')
+                           break}
+
+      if((n[maxemi,1] != m[maxemi,1]) & (n[maxemi,2] != m[maxemi,2])) 
+      {
+        helpm <- matrix(nrow = ((dim(m)[1])+2), ncol=3)
+        helpm[maxemi,1] <- m[[maxemi,1]]
+        helpm[maxemi,2] <- n[maxemi,1]
+        helpm[maxemi+1,1] <- n[maxemi,1]
+        helpm[maxemi+1,2] <- n[maxemi,2]
+        helpm[maxemi+2,1] <- n[maxemi,2]
+        helpm[maxemi+2,2] <- m[[maxemi,2]]
+        if((n[maxemi,3])>0)
+        {
+          helpm[maxemi,3] <- TRUE
+          helpm[maxemi+1,3] <- FALSE
+          helpm[maxemi+2,3] <- TRUE
+        }else{  
+          helpm[maxemi,3] <- FALSE
+          helpm[maxemi+1,3] <- TRUE
+          helpm[maxemi+2,3] <- FALSE
+        }
+        if(maxemi < (dim(m)[1])){
+            helpm[(maxemi+1):(dim(m)[1])+2,]=m[(maxemi+1):(dim(m)[1]),]
+        }
+      }else{      
+        helpm <- matrix(nrow = ((dim(m)[1])+1), ncol=3)
+        if((n[maxemi,1]) == 1) 
+        {
+          helpm[maxemi,1] <- 1 
+          helpm[maxemi,2] <- n[maxemi,2]
+          helpm[maxemi+1,1] <- n[maxemi,2]
+          helpm[maxemi+1,2] <- m[[maxemi,2]]        
+          helpm[maxemi,3] <- FALSE
+          helpm[maxemi+1,3] <- TRUE
+        }else{                        
+          helpm[maxemi,1] <- m[[maxemi,1]]
+          helpm[maxemi,2] <- n[maxemi,1]
+          helpm[maxemi+1,1] <- n[maxemi,1]
+          helpm[maxemi+1,2] <- length(x)
+          helpm[maxemi,3] <- TRUE
+          helpm[maxemi+1,3] <- FALSE
+        }
+        if(maxemi < (dim(m)[1])){    
+            helpm[((maxemi+1):(dim(m)[1])+1),]=m[(maxemi+1):(dim(m)[1]),]
+        }
+      }
+      if(maxemi > 1){ 
+        helpm[1:(maxemi-1),]=m[1:(maxemi-1),]
+      }  
+      m <- helpm
+      variation[i] <- variation[i-1] + abs(n[maxemi,3])
+      if(UpToM==TRUE)
+      {
+        j <- 1
+        u <- matrix(nrow=(dim(m)[1]-sum(m[,3])), ncol=4,dimnames=list(1:(dim(m)[1]-sum(m[,3])),
+                                                                      c("Start Index","End Index","Start Point","End Point")))
+        for(i in 1:(dim(m)[1]))
+        {
+          if(m[i,3]==0)
+          {
+            u[j,1] <- m[i,1]
+            u[j,2] <- m[i,2]
+            u[j,3] <- x[m[i,1]]
+            u[j,4] <- x[m[i,2]]
+            j <- j+1
+          }
+        }
+        r[[length(r)+1]] <- u
+      }
+    }
+  }
+  if(UpToM==TRUE)
+  {
+    r <- list("intervals"=r, "excess_mass"=variation)
+  }else{
+    j <- 1
+    u <- matrix(nrow=(dim(m)[1]-sum(m[,3])), ncol=4,dimnames=list(1:(dim(m)[1]-sum(m[,3])),
+                                                                  c("Start Index","End Index","Start Point","End Point")))
+    for(i in 1:(dim(m)[1]))
+    {
+      if(m[i,3]==0)
+      {
+        u[j,1] <- m[i,1]
+        u[j,2] <- m[i,2]
+        u[j,3] <- x[m[i,1]]
+        u[j,4] <- x[m[i,2]]
+        j <- j+1
+      }
+    }  
+    r <- list("intervals"=u,"excess_mass"=variation)
+  }
+  return(r)
+  }
+}

R/exmplot.R

@@ -0,0 +1,103 @@
+exmplot <- function(xdata, M=1, CutOff=1, steps=50, Lambda=NULL){
+if(all(is.na(xdata))){stop('Input vector is empty.')}else{
+  options(warn=-1)
+  if(is.null(Lambda) == TRUE){
+    mLambda <- searchMaxLambda(xdata,4*CutOff)
+    r <- array(list(), c(steps,2))
+    for(i in 1:steps)
+    {
+      r[i,] <- excessm(xdata, (i* (mLambda/steps)), M=M)
+    }
+    #because of monotonicity smallest excess mass is achied under M=1, largest Lambda
+    #because of monotonicity largest excess mass is achied under M=M, smallest Lambda
+    plot(1, type="n", xlim=c(r[[steps,2]][[1]],1.025*r[[1,2]][[length(r[1,2])]]), ylim=c(mLambda/steps,mLambda), 
+         ylab="Lambda", xlab="Excess Mass", main='Excess Mass Lambda Plot', cex.lab=1.1, cex.main=1.65)
+    for(i in 1:M){
+      for(j in 1:(steps-1)){
+        l1 = length(r[[j,2]])
+        l2 = length(r[[j+1,2]])
+        if(l1 >= i ){
+          if(l2 >= i){
+            lines(x=c(r[[j,2]][[i]], r[[j+1,2]][[i]]),y=c(j* mLambda/steps, (j+1)* mLambda/steps), lwd=1.7)
+          }else{
+            lines(x=c(r[[j,2]][[i]], r[[j+1,2]][[l2]]),y=c(j* mLambda/steps, (j+1)* mLambda/steps), lwd=1.7)
+          }          
+        }else{
+          if(l2 >= i){
+            lines(x=c(r[[j,2]][[l1]], r[[j+1,2]][[i]]),y=c(j* mLambda/steps, (j+1)* mLambda/steps), lwd=1.7)
+          }else{
+            lines(x=c(r[[j,2]][[l1]], r[[j+1,2]][[l2]]),y=c(j* mLambda/steps, (j+1)* mLambda/steps), lwd=1.7)
+          }
+        }
+      }
+    }
+    max_dist <- vector(length=M-1)
+    max_dist_Lambda <- vector(length=M-1)
+    for(i in 2:M){
+      max_val <- 0
+      max_ind <- 0
+      for(j in 1:steps){
+        if(length(r[[j,2]]) >= i){ #otherwise the excess mass for i and i+1 are equal
+          curval <-  r[[j,2]][[i]] - r[[j,2]][[i-1]]
+          if(curval > max_val){
+            max_val <- curval
+            max_ind <- j
+          }
+        }
+      }
+      max_dist[i-1] <- max_val
+      max_dist_Lambda[i-1] <- max_ind* mLambda/steps #Lambda for which the maximal distance is achieved 
+    }
+
+  }else{
+    Lambda <- sort(Lambda)
+    r <- array(list(), c(length(Lambda),2))
+    for(i in 1:length(Lambda))
+    {
+      r[i,] <- excessm(xdata, Lambda[i], M=M)
+    }
+    plot(1, type="n", xlim=c(r[[length(Lambda),2]][[1]],1.025*r[[1,2]][[length(r[1,2])]]), ylim=c(Lambda[1],Lambda[length(Lambda)]), 
+         ylab="Lambda", xlab="Excess Mass", main='Excess Mass Lambda Plot', cex.lab=1.1, cex.main=1.65)
+    for(i in 1:M){
+      for(j in 1:(length(Lambda)-1)){
+        l1 = length(r[[j,2]])
+        l2 = length(r[[j+1,2]])
+        if(l1 >= i ){
+          if(l2 >=i){
+            lines(x=c(r[[j,2]][[i]], r[[j+1,2]][[i]]),y=c(Lambda[j], Lambda[j+1]), lwd=1.7)
+          }else{
+            lines(x=c(r[[j,2]][[i]], r[[j+1,2]][[l2]]),y=c(Lambda[j], Lambda[j+1]), lwd=1.7)
+          }          
+        }else{
+          if(l2 >= i){
+            lines(x=c(r[[j,2]][[l1]], r[[j+1,2]][[i]]),y=c(Lambda[j], Lambda[j+1]), lwd=1.7)
+          }else{
+            lines(x=c(r[[j,2]][[l1]], r[[j+1,2]][[l2]]),y=c(Lambda[j], Lambda[j+1]), lwd=1.7)
+          }
+        }
+      }
+    }
+    max_dist <- vector(length=M-1)
+    max_dist_Lambda <- vector(length=M-1)
+    for(i in 2:M){
+      max_val <- 0
+      max_ind <- 0
+      for(j in 1:length(Lambda)){
+        if(length(r[[j,2]]) >= i){ 
+          curval <- r[[j,2]][[i]] - r[[j,2]][[i-1]]
+          if(curval > max_val){
+            max_val <- curval
+            max_ind <- j
+          }
+        }
+      }
+      max_dist[i-1] <- max_val
+      max_dist_Lambda[i-1] <- Lambda[max_ind] #Lambda for which the maximal distance is achieved 
+    }
+
+  }  
+
+  u <- list("max_dist"=max_dist, "max_dist_Lambda"=max_dist_Lambda)
+  return(u)
+}
+}