R : Copyright 2005, The R Foundation for Statistical Computing
Version 2.1.1  (2005-06-20), ISBN 3-900051-07-0

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for a HTML browser interface to help.
Type 'q()' to quit R.

> ### * <HEADER>
> ###
> attach(NULL, name = "CheckExEnv")
> assign(".CheckExEnv", as.environment(2), pos = length(search())) # base
> ## add some hooks to label plot pages for base and grid graphics
> setHook("plot.new", ".newplot.hook")
> setHook("persp", ".newplot.hook")
> setHook("grid.newpage", ".gridplot.hook")
> 
> assign("cleanEx",
+        function(env = .GlobalEnv) {
+ 	   rm(list = ls(envir = env, all.names = TRUE), envir = env)
+            RNGkind("default", "default")
+ 	   set.seed(1)
+    	   options(warn = 1)
+ 	   delayedAssign("T", stop("T used instead of TRUE"),
+ 		  assign.env = .CheckExEnv)
+ 	   delayedAssign("F", stop("F used instead of FALSE"),
+ 		  assign.env = .CheckExEnv)
+ 	   sch <- search()
+ 	   newitems <- sch[! sch %in% .oldSearch]
+ 	   for(item in rev(newitems))
+                eval(substitute(detach(item), list(item=item)))
+ 	   missitems <- .oldSearch[! .oldSearch %in% sch]
+ 	   if(length(missitems))
+ 	       warning("items ", paste(missitems, collapse=", "),
+ 		       " have been removed from the search path")
+        },
+        env = .CheckExEnv)
> assign("..nameEx", "__{must remake R-ex/*.R}__", env = .CheckExEnv) # for now
> assign("ptime", proc.time(), env = .CheckExEnv)
> grDevices::postscript("ks-Examples.ps")
> assign("par.postscript", graphics::par(no.readonly = TRUE), env = .CheckExEnv)
> options(contrasts = c(unordered = "contr.treatment", ordered = "contr.poly"))
> options(warn = 1)    
> library('ks')
Loading required package: mvtnorm
ks version 1.3.2 (2005) 
> 
> assign(".oldSearch", search(), env = .CheckExEnv)
> assign(".oldNS", loadedNamespaces(), env = .CheckExEnv)
> cleanEx(); ..nameEx <- "Hbcv"
> 
> ### * Hbcv
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: Hbcv, Hbcv.diag
> ### Title: Biased cross-validation (BCV) bandwidth matrix selector for
> ###   bivariate data
> ### Aliases: Hbcv Hbcv.diag
> ### Keywords: smooth
> 
> ### ** Examples
> 
> data(unicef)
> Hbcv(unicef)
          [,1]      [,2]
[1,] 1087.0682 135.33067
[2,]  135.3307  23.58613
> Hbcv.diag(unicef)
         [,1]     [,2]
[1,] 1072.781 0.000000
[2,]    0.000 9.298426
> 
> 
> 
> cleanEx(); ..nameEx <- "Hkda"
> 
> ### * Hkda
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: Hkda, Hkda.diag
> ### Title: Bandwidth matrix selectors for kernel discriminant analysis for
> ###   multivariate data
> ### Aliases: Hkda Hkda.diag
> ### Keywords: smooth
> 
> ### ** Examples
> 
> library(MASS)
> data(iris)
> iris.mat <- rbind(iris[,,1], iris[,,2], iris[,,3])
> ir <- iris.mat[,c(1,2)]
> ir.gr <- iris.mat[,5]
> 
> Hkda(ir, ir.gr, bw="scv", pre="scale")
            [,1]        [,2]
[1,] 0.008621007 0.001417174
[2,] 0.001417174 0.011827099
[3,] 0.035043481 0.007630606
[4,] 0.007630606 0.013803803
[5,] 0.047315514 0.001424001
[6,] 0.001424001 0.012836026
> Hkda.diag(ir, ir.gr, bw="plugin", pre="scale")
            [,1]        [,2]
[1,] 0.005176874 0.000000000
[2,] 0.000000000 0.006327773
[3,] 0.017391151 0.000000000
[4,] 0.000000000 0.007579993
[5,] 0.031348803 0.000000000
[6,] 0.000000000 0.008827732
> 
> 
> 
> cleanEx(); ..nameEx <- "Hlscv"
> 
> ### * Hlscv
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: Hlscv, Hlscv.diag
> ### Title: Least-squares cross-validation (LSCV) bandwidth matrix selector
> ###   for multivariate data
> ### Aliases: Hlscv Hlscv.diag
> ### Keywords: smooth
> 
> ### ** Examples
> 
> data(unicef)
> Hlscv(unicef)
          [,1]       [,2]
[1,] 1140.9209 -142.03486
[2,] -142.0349   24.75457
> Hlscv.diag(unicef)
         [,1]    [,2]
[1,] 1125.925 0.00000
[2,]    0.000 9.75909
> 
> 
> 
> cleanEx(); ..nameEx <- "Hmise"
> 
> ### * Hmise
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: Hmise.mixt, Hamise.mixt
> ### Title: MISE- and AMISE-optimal bandwidth matrix selectors for normal
> ###   mixture densities
> ### Aliases: Hmise.mixt Hamise.mixt
> ### Keywords: smooth
> 
> ### ** Examples
> 
> mus <- rbind(c(0,0,0), c(2,2,2))
> Sigma <- matrix(c(1, 0.7, 0.7, 0.7, 1, 0.7, 0.7, 0.7, 1), nr=3, nc=3) 
> Sigmas <- rbind(Sigma, Sigma)
> props <- c(1/2, 1/2)
> samp <- 1000
> Hmise.mixt(mus, Sigmas, props, samp)
          [,1]      [,2]      [,3]
[1,] 0.2800151 0.2355075 0.2355013
[2,] 0.2355075 0.2798914 0.2354426
[3,] 0.2355013 0.2354426 0.2798487
> Hamise.mixt(mus, Sigmas, props, samp)
          [,1]      [,2]      [,3]
[1,] 0.2213553 0.1798598 0.1798577
[2,] 0.1798598 0.2213521 0.1798569
[3,] 0.1798577 0.1798569 0.2213505
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "Hpi"
> 
> ### * Hpi
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: Hpi, Hpi.diag
> ### Title: Plug-in bandwidth matrix selector for multivariate data
> ### Aliases: Hpi Hpi.diag
> ### Keywords: smooth
> 
> ### ** Examples
> 
> ###bivariate example
> 
> data(unicef)
> Hpi(unicef, nstage=1, pilot="amse", pre="scale")
          [,1]       [,2]
[1,] 391.03446 -34.733738
[2,] -34.73374   9.898164
> Hpi(unicef, nstage=2, pilot="samse", pre="sphere")
          [,1]       [,2]
[1,]  797.5755 -106.63338
[2,] -106.6334   19.56761
> Hpi.diag(unicef, nstage=2, pilot="amse", pre="scale") 
         [,1]     [,2]
[1,] 201.5118 0.000000
[2,]   0.0000 6.242821
> 
> ### trivariate example
> 
> mus <- rbind(c(0,0,0), c(2,2,2))
> Sigma <- matrix(c(1, 0.7, 0.7, 0.7, 1, 0.7, 0.7, 0.7, 1), nr=3, nc=3) 
> Sigmas <- rbind(Sigma, Sigma)
> x <- rmvnorm.mixt(n=100, mus=mus, Sigmas=Sigmas, props=c(1/2, 1/2))
> Hpi(x)
          [,1]      [,2]      [,3]
[1,] 0.4429662 0.3764609 0.3932967
[2,] 0.3764609 0.4268842 0.3897707
[3,] 0.3932967 0.3897707 0.4817339
> Hscv(x)
          [,1]      [,2]      [,3]
[1,] 0.8161399 0.7085630 0.7618059
[2,] 0.7085630 0.8094901 0.7444281
[3,] 0.7618059 0.7444281 0.9001240
> 
> 
> 
> cleanEx(); ..nameEx <- "Hscv"
> 
> ### * Hscv
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: Hscv
> ### Title: Smoothed cross-validation (SCV) bandwidth matrix selector for
> ###   bivariate data
> ### Aliases: Hscv
> ### Keywords: smooth
> 
> ### ** Examples
> 
> data(unicef)
> Hscv(unicef)
          [,1]       [,2]
[1,] 1322.3429 -191.76583
[2,] -191.7658   34.98503
> 
> 
> 
> cleanEx(); ..nameEx <- "kda"
> 
> ### * kda
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: kda, pda, compare
> ### Title: Kernel and parametric discriminant analysis
> ### Aliases: kda pda compare
> ### Keywords: smooth
> 
> ### ** Examples
> 
> 
> ### bivariate example - restricted iris dataset  
> library(MASS)
> data(iris)
> iris.mat <- rbind(iris[,,1], iris[,,2], iris[,,3])
> ir <- iris.mat[,c(1,2)]
> ir.gr <- iris.mat[,5]
> 
> H <- Hkda(ir, ir.gr, bw="plugin", pre="scale")
> kda.gr <- kda(ir, ir.gr, H, ir)
> lda.gr <- pda(ir, ir.gr, ir, type="line")
> qda.gr <- pda(ir, ir.gr, ir, type="quad")
> 
> compare(kda.gr, ir.gr)
$cross
           setosa versicolor virginica    
setosa        150          0         0 150
versicolor      0        117        30 147
virginica       0         33       120 153
              150        150       150 450

$error
[1] 0.14

> compare(qda.gr, ir.gr)
$cross
           setosa versicolor virginica    
setosa        147          0         0 147
versicolor      3        111        48 162
virginica       0         39       102 141
              150        150       150 450

$error
[1] 0.2

> compare(lda.gr, ir.gr)
$cross
           setosa versicolor virginica    
setosa        147          0         0 147
versicolor      3        108        45 156
virginica       0         42       105 147
              150        150       150 450

$error
[1] 0.2

> 
> ### multivariate example - full iris dataset
> ir <- iris[,1:4]
> ir.gr <- iris[,5]
> 
> H <- Hkda(ir, ir.gr, bw="plugin", pre="scale")
> kda.gr <- kda(ir, ir.gr, H, ir)
> lda.gr <- pda(ir, ir.gr, ir, type="line")
> qda.gr <- pda(ir, ir.gr, ir, type="quad")
> 
> compare(kda.gr, ir.gr)
$cross
           setosa versicolor virginica    
setosa         50          0         0  50
versicolor      0         50         0  50
virginica       0          0        50  50
               50         50        50 150

$error
[1] 0

> compare(qda.gr, ir.gr)
$cross
           setosa versicolor virginica    
setosa         50          0         0  50
versicolor      0         48         1  49
virginica       0          2        49  51
               50         50        50 150

$error
[1] 0.02

> compare(lda.gr, ir.gr)
$cross
           setosa versicolor virginica    
setosa         50          0         0  50
versicolor      0         48         1  49
virginica       0          2        49  51
               50         50        50 150

$error
[1] 0.02

> 
> 
> 
> cleanEx(); ..nameEx <- "kda.kde"
> 
> ### * kda.kde
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: kda.kde, pda.pde
> ### Title: Density estimates for kernel and parametric discriminant
> ###   analysis for multivariate data
> ### Aliases: kda.kde pda.pde
> ### Keywords: smooth
> 
> ### ** Examples
> 
> library(MASS)
> data(iris)
> iris.mat <- rbind(iris[,,1], iris[,,2], iris[,,3])
> ir <- iris.mat[,c(1,2)]
> ir.gr <- iris.mat[,5]
> 
> H <- Hkda(ir, ir.gr, bw="plugin", pre="scale")
> kda.gr <- kda(ir, ir.gr, H, ir)
> fhat <- kda.kde(ir, ir.gr, H, gridsize=c(250,250))
> 
> qda.gr <- pda(ir, ir.gr, ir, type="quad")
> qda.fhat <- pda.pde(ir, ir.gr, gridsize=c(250,250))
> 
> 
> 
> 
> cleanEx(); ..nameEx <- "kde"
> 
> ### * kde
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: kde
> ### Title: Kernel density estimate for multivariate data
> ### Aliases: kde
> ### Keywords: smooth
> 
> ### ** Examples
> 
> ### bivariate example
> data(unicef)
> H.pi <- Hpi(unicef, nstage=1)
> fhat <- kde(unicef, H.pi)
> 
> ### trivariate example
> mus <- rbind(c(0,0,0), c(2,2,2))
> Sigma <- matrix(c(1, 0.7, 0.7, 0.7, 1, 0.7, 0.7, 0.7, 1), nr=3, nc=3) 
> Sigmas <- rbind(Sigma, Sigma)
> props <- c(1/2, 1/2)
> x <- rmvnorm.mixt(n=100, mus=mus, Sigmas=Sigmas, props=props)
> H.pi <- Hpi(x)
> fhat <- kde(x, H.pi, eval.levels=seq(-3,3, length=9)) 
> 
> ### 4-variate example
> library(MASS)
> data(iris)
> ir <- iris[,1:4][iris[,5]=="setosa",]
> H.scv <- Hscv(ir)
> fhat <- kde(ir, H.scv, eval.points=ir)  
> 
> 
> 
> cleanEx(); ..nameEx <- "mise"
> 
> ### * mise
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: ise, mise, amise
> ### Title: ISE, MISE and AMISE of kernel density estimates for normal
> ###   mixture densities
> ### Aliases: ise.mixt mise.mixt amise.mixt
> ### Keywords: smooth
> 
> ### ** Examples
> 
> samp <- 50
> d <- 3
> mus <- rbind(rep(0,d), rep(1,d))
> Sigmas <- 0.25*rbind(diag(d), diag(d))
> props <- c(2/3, 1/3)
> x <- rmvnorm.mixt(samp, mus, Sigmas, props)
> H <- Hpi(x)
> ise.mixt(x, H, mus, Sigmas, props)
[1] 0.03099587
> mise.mixt(H, mus, Sigmas, props, samp)
[1] 0.03117458
> 
> 
> 
> cleanEx(); ..nameEx <- "mvnorm"
> 
> ### * mvnorm
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: rmvnorm.mixt, dmvnorm.mixt
> ### Title: Multivariate normal mixture distribution
> ### Aliases: rmvnorm.mixt dmvnorm.mixt
> ### Keywords: distribution
> 
> ### ** Examples
> 
> mus <- rbind(c(-3/2,0), c(3/2,0))
> Sigmas <- rbind(diag(c(1/16, 1)), rbind(c(1/16, 1/18), c(1/18, 1/16)))
> props <- c(2/3, 1/3)
> x <- rmvnorm.mixt(1000, mus, Sigmas, props)
> dens <- dmvnorm.mixt(x, mus, Sigmas, props)
> 
> 
> 
> cleanEx(); ..nameEx <- "mvt"
> 
> ### * mvt
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: rmvt.mixt, dmvt.mixt
> ### Title: Multivariate t mixture distribution
> ### Aliases: rmvt.mixt dmvt.mixt
> ### Keywords: distribution
> 
> ### ** Examples
> 
> mus <- rbind(c(-3/2,0), c(3/2,0))
> Sigmas <- rbind(diag(c(1/16, 1)), rbind(c(1/16, 1/18), c(1/18, 1/16)))
> props <- c(2/3, 1/3)
> dfs <- c(7,3)
> x <- rmvt.mixt(1000, mus, Sigmas, dfs, props)
> dens <- dmvt.mixt(x, mus, Sigmas, dfs, props)
> 
> 
> 
> cleanEx(); ..nameEx <- "plot.dade"
> 
> ### * plot.dade
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: plot.dade
> ### Title: Density estimates and partition plot for discriminant analysis
> ###   for bivariate data
> ### Aliases: plot.dade
> ### Keywords: hplot
> 
> ### ** Examples
> 
> 
> library(MASS)
> data(iris)
> iris.mat <- rbind(iris[,,1], iris[,,2], iris[,,3])
> ir <- iris.mat[,c(1,2)]
> ir.gr <- iris.mat[,5]
> xlab <- "Sepal length (mm)"
> ylab <- "Sepal width (mm)"
> xlim <- c(4,8)
> ylim <- c(2,4.5)
> 
> H <- Hkda(ir, ir.gr, bw="plugin", pre="scale")
> fhat <- kda.kde(ir, ir.gr, H, gridsize=c(250,250))
> lda.fhat <- pda.pde(ir, ir.gr, type="line")
> qda.fhat <- pda.pde(ir, ir.gr, type="quad")
> 
> layout(rbind(c(1,2), c(3,4)))
> plot(fhat, cont=0, xlab=xlab, ylab=ylab, xlim=xlim, ylim=ylim, 
+     pch=c(1,5,10))
> plot(fhat, ncont=6, xlab=xlab, ylab=ylab, xlim=xlim, ylim=ylim,
+      col=c("transparent", "grey", "#8f8f8f"), drawlabels=FALSE)
> plot(lda.fhat, ncont=6, xlim=xlim, ylim=ylim, xlab=xlab, ylab=ylab, 
+      disp="")
> plot(qda.fhat, ncont=6, xlim=xlim, ylim=ylim,  xlab=xlab, ylab=ylab, 
+      lty=c(2,5,3))
> layout(1)
> 
> 
> 
> cleanEx(); ..nameEx <- "plot.kde"
> 
> ### * plot.kde
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: plot.kde
> ### Title: Kernel density estimate plot for 2- and 3-dimensional data
> ### Aliases: plot.kde
> ### Keywords: hplot
> 
> ### ** Examples
> 
> ### bivariate example
> data(unicef)
> H.scv <- Hscv(unicef)
> fhat <- kde(unicef, H.scv)
> 
> layout(rbind(c(1,2), c(3,4)))
> plot(fhat, display="slice", cont=seq(10,90, by=20), cex=0.3)
> plot(fhat, display="slice", ncont=5, cex=0.3, drawlabels=FALSE)
> plot(fhat, display="persp")
> plot(fhat, display="image", col=rev(heat.colors(15)))
> layout(1)
> 
> ### 3-variate example
> 
> mus <- rbind(c(0,0,0), c(2,2,2))
> Sigma <- matrix(c(1, 0.7, 0.7, 0.7, 1, 0.7, 0.7, 0.7, 1), nr=3, nc=3) 
> Sigmas <- rbind(Sigma, Sigma)
> props <- c(1/2, 1/2)
> x <- rmvnorm.mixt(n=100, mus=mus, Sigmas=Sigmas, props=props)
> H.pi <- Hpi(x)
> fhat <- kde(x, H.pi, eval.levels=seq(-3,3, length=9))  
> plot(fhat, disp="slice", ncont=6, cex=0.3)
> 
> 
> 
> cleanEx(); ..nameEx <- "pre.scale"
> 
> ### * pre.scale
> 
> flush(stderr()); flush(stdout())
> 
> ### Name: pre.sphere, pre.scale
> ### Title: Pre-sphering and pre-scaling
> ### Aliases: pre.sphere pre.scale
> ### Keywords: algebra
> 
> ### ** Examples
> 
> x <- rmvnorm.mixt(1000, mus=c(0,0), Sigmas=rbind(c(1,0.2), c(0.2, 0.5))) 
> x.sp <- pre.sphere(x)
> x.sc <- pre.scale(x)
> var(x.sp)
             [,1]         [,2]
[1,]  1.00000e+00 -7.63607e-16
[2,] -7.63607e-16  1.00000e+00
> var(x.sc)
          [,1]      [,2]
[1,] 1.0000000 0.2884959
[2,] 0.2884959 1.0000000
> 
> 
> 
> ### * <FOOTER>
> ###
> cat("Time elapsed: ", proc.time() - get("ptime", env = .CheckExEnv),"\n")
Time elapsed:  187.61 6.41 209.24 0 0 
> grDevices::dev.off()
null device 
          1 
> ###
> ### Local variables: ***
> ### mode: outline-minor ***
> ### outline-regexp: "\\(> \\)?### [*]+" ***
> ### End: ***
> quit('no')