Skip to content

Commit

Permalink
version 1.0-0
Browse files Browse the repository at this point in the history
  • Loading branch information
@gaborcsardi
Jean-Francois Plante authored and gaborcsardi committed Sep 1, 2008
0 parents commit 788115f
Show file tree
Hide file tree
Showing 16 changed files with 6,479 additions and 0 deletions.
12 changes: 12 additions & 0 deletions DESCRIPTION
@@ -0,0 +1,12 @@
Package: depth
Type: Package
Title: Depth functions tools for multivariate analysis
Version: 1.0-0
Date: 2008-09-01
Author: Jean-Claude Masse <jcmasse@mat.ulaval.ca), Jean-Francois Plante <plante@utstat.toronto.edu>.
Maintainer: Jean-Francois Plante <plante@utstat.toronto.edu>
Description: Depth functions methodology applied to multivariate analysis. Besides allowing calculation of depth values and depth-based location estimators, the package includes functions for drawing contour plots and perspective plots of depth functions.
Depends: R (>= 2.4.0), grDevices, rgl
Suggests: robustbase, MASS
License: GPL-2
Packaged: Mon Oct 27 14:05:46 2008; pascalegiroux-denis
88 changes: 88 additions & 0 deletions R/ctrmean.R
@@ -0,0 +1,88 @@

ctrmean=function(x,alpha,eps=1e-8,mustdith=FALSE,maxdith=50,dithfactor=10,factor=.8){
UseMethod("ctrmean")
}

ctrmean.data.frame=function(x,alpha,eps=1e-8,mustdith=FALSE,maxdith=50,dithfactor=10,factor=.8){
x=as.matrix(x)
NextMethod("ctrmean",x)
}

ctrmean.list=function(x,alpha,eps=1e-8,mustdith=FALSE,maxdith=50,dithfactor=10,factor=.8){
m=length(x)
n=length(x[[1]])
y=matrix(0,n,m)
for(i in 1:m){
y[,i]=x[[i]]
if(length(x[[i]])!=n){ stop("When using a list, each element must be a vector of the same length.") }
}
x=y
NextMethod("ctrmean",x)
}

ctrmean.default=function(x,alpha,eps=1e-8,mustdith=FALSE,maxdith=50,dithfactor=10,factor=.8){

p=length(x[1,])
n=length(x[,1])
depth=ceiling(n*alpha)

if(p>n) { warning(message=paste("Is your data ",n," points in ",p," dimensions.\nIf not, you should transpose your data matrix.\n")) }
if(p!=2){ stop("Depth contours can only be calculated on bivaraite data.") }
if(length(depth)!=1 | round(depth)!=depth){ stop("The argument depth must be a single integer") }
if(depth<1 | depth>ceiling(n/2)){ stop(message=paste("Depth must be an integer between 1 and ",ceiling(n/2))) }

ndpth=1
y=x[,2]
x=x[,1]
maxnum=floor(4*n*sqrt(n)+1)


zz=.Fortran("halfmed",
as.numeric(x),
as.numeric(y),
as.integer(n),
integer(1),
numeric(2),
xcont=numeric(n*(n-1)/2),
ycont=numeric(n*(n-1)/2),
ncont=integer(ndpth),
as.integer(depth),
as.integer(ndpth),
as.integer(1),
as.integer(maxnum),
err=integer(1),
as.numeric(eps),
as.numeric(dithfactor),
as.integer(maxdith),
as.integer(mustdith),
missing=integer(ndpth),
as.numeric(factor),
PACKAGE="depth")


if (zz$err==-5) { warning("Ventilation was used on the data.") }
if (zz$err==-1) {
if(mustdith){
warning("Data are not in general position. Dithering was used.")
} else {
stop("Data are not in general position.")
}
}
if (zz$err==-2) { stop(message=paste(maxdith," ventilation steps and the data still fail to be in general position.")) }
if (zz$err==-4) { warning("A non critical numerical error occurred during the calculation of the depth contour.") }
if (zz$err>0) { stop(message=paste("CRITICAL numerical error in calculating contour")) }
if(zz$ncont==0){ warning("The requested contour does not exist.") }
n=zz$ncont
x=zz$xcont[c(1:n,1)]
y=zz$ycont[c(1:n,1)]


A=sum(x[-(n+1)]*y[-1]-x[-1]*y[-(n+1)])/2
xb=sum((x[-(n+1)]+x[-1])*(x[-(n+1)]*y[-1]-x[-1]*y[-(n+1)]))/6/A
yb=sum((y[-(n+1)]+y[-1])*(x[-(n+1)]*y[-1]-x[-1]*y[-(n+1)]))/6/A

c(xb,yb)

}


205 changes: 205 additions & 0 deletions R/depth.R
@@ -0,0 +1,205 @@
depth=function(u,x,method="Tukey",approx=FALSE,eps=1e-8,ndir=1000){
UseMethod("depth",x)
}

depth.matrix=function(u,x,method="Tukey",approx=FALSE,eps=1e-8,ndir=1000){
NextMethod("depth",x)
}

depth.data.frame=function(u,x,method="Tukey",approx=FALSE,eps=1e-8,ndir=1000){
x=as.matrix(x)
NextMethod("depth",x)
}

depth.list=function(u,x,method="Tukey",approx=FALSE,eps=1e-8,ndir=1000){
m=length(x)
n=length(x[[1]])
y=matrix(0,n,m)
for(i in 1:m){
y[,i]=x[[i]]
if(length(x[[i]])!=n){ stop("When using a list, each element must be a vector of the same length.") }
}
x=y
NextMethod("depth",x)
}

depth.default=function(u,x,method="Tukey",approx=FALSE,eps=1e-8,ndir=1000){

# Suppose n data points in p dimensions.
# Calculates depth of u in sample x.
# u= vector of p real numbers
# x= matrix n by p
# method= which definition of depth to use


match.arg(method,c("Tukey","Liu","Oja"))
if(!is.matrix(x)){
n=length(x)
if(length(u)>1){stop("Data is univariate, but u is mulativariate.")}
if(method=="Tukey"){
return(min(sum(x<=u),sum(x>=u))/n)
}
if(method=="Liu"){
return(sum(outer(u-x,u-x,function(a,b) a*b)<=0)/(n*(n-1)))
}
if(method=="Oja"){
return(.5*(1+mean(abs(u-x)))^(-1))
}
}
p=length(x[1,])
n=length(x[,1])

if(p>n) { warning(message=paste("Is your data ",n," points in ",p," dimensions.\nIf not, you should transpose your data matrix.\n")) }
if(length(u)!=p){ stop("Dimension mismatch between the data and the point u.") }
if(method!="Tukey"&approx==TRUE){ warning("An approximate algorithm is available only for Tukey's depth. Argument approx=TRUE is ignored.") }
if(method=="Liu"&p>2){ stop("Liu's depth can be calculated on bivariate datasets only.") }
if(method=="Tukey"&approx==FALSE&p>3){ warning("Tukey's depth can only be approximated in more than 3 dimensions. Argument approx=F is ignored.") }

# Calculation of Tukey's Median

if(method=="Tukey"){
if(p==2){

ans=.Fortran("depth",
as.numeric(u[1]),
as.numeric(u[2]),
as.integer(n),
as.numeric(x[,1]),
as.numeric(x[,2]),
numeric(n),
integer(n),
sdep=as.numeric(-1),
hdep=as.numeric(-1),
PACKAGE="depth")

rep=ans$hdep
}

if(p==3&&approx==FALSE) {

zz=.Fortran("stand",
as.integer(n),
x=as.numeric(x[,1]),
y=as.numeric(x[,2]),
z=as.numeric(x[,3]),
u=as.numeric(u[1]),
v=as.numeric(u[2]),
w=as.numeric(u[3]),
xn=numeric(n),
as.numeric(eps),
err=as.integer(0),
PACKAGE="depth")

if(zz$err!=0) {
if(zz$err>0&&zz$err<10) { stop(message=paste("At least half of the data set share the same value for variable ",zz$err," ."))}
if(zz$err>10){ stop(message=paste("Variable ",zz$err-10," has null covariance."))}
}

ans=.Fortran("depth3",
as.integer(n),
as.numeric(zz$u),
as.numeric(zz$v),
as.numeric(zz$w),
as.numeric(zz$x),
as.numeric(zz$y),
as.numeric(zz$z),
numeric(n),
integer(n),
numeric(n),
numeric(n),
as.numeric(eps),
ndim=as.integer(0),
ndep=as.integer(0),
PACKAGE="depth")

rep=ans$ndep/n
}

if(approx==TRUE||p>3){

zz=.Fortran("standpd",
as.integer(n),
as.integer(p),
n=as.integer(n),
p=as.integer(p),
x=as.numeric(x),
u=as.numeric(u),
as.numeric(eps),
err=integer(p),
ndep=as.numeric(n),
PACKAGE="depth")

if(zz$ndep==0||zz$p==0) {
cat("\nSample is singular.\n")
stop()
}


z=.Fortran("hdepth",
as.integer(n),
as.integer(p),
nnp=as.integer(p),
as.integer(ndir),
as.integer(n),
as.integer(p),
as.numeric(zz$x),
integer(p),
as.numeric(zz$u),
numeric(p),
numeric(p*p),
numeric(p),
numeric(p*p),
numeric(p),
as.numeric(eps),
err=as.numeric(0),
err2=numeric(ndir),
ndep=as.integer(0),
as.integer(0),
PACKAGE="depth")

if (z$err==-1) { warning(message=paste("\nThe dataset is singular.\nIts dimension was reduced to ",z$nnp, " dimensions.")) }
if (z$err==-2) { stop("Error in dimension reduction: eigenvectors are not independent.") }
for (i in 1:ndir){
if (z$err2[i]==-1) { stop(message=paste("\nNo eigenvalue for datum ",z$err2[i])) }
if (z$err2[i]==-2) { stop(message=paste("\nNull eigenvector for datum ",z$err2[i])) }
if (z$err2[i]>0) { stop(message=paste("\nError ",z$err2[i], " in the calculation of eigenvalues.")) }
}

rep=z$ndep/n
}
}

# Calculation of Liu's median

if(method=="Liu"){
ans=.Fortran("depth",
as.numeric(u[1]),
as.numeric(u[2]),
as.integer(n),
as.numeric(x[,1]),
as.numeric(x[,2]),
numeric(n),
integer(n),
sdep=as.numeric(-1),
hdep=as.numeric(-1),
PACKAGE="depth")

rep=ans$sdep

}

if(method=="Oja"){
ans=.Fortran("ojadepth",
as.numeric(x),
as.numeric(u),
as.integer(p),
as.integer(n),
depth=as.numeric(-1),
PACKAGE="depth")

rep=ans$depth

}

rep
}
2 changes: 2 additions & 0 deletions R/firstlib.R
@@ -0,0 +1,2 @@
.First.lib <- function(libpath, pkgname)
library.dynam("depth", pkgname, libpath)

0 comments on commit 788115f

Please sign in to comment.