R : Copyright 2005, The R Foundation for Statistical Computing Version 2.1.0 Patched (2005-05-12), 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. > ### *
> ### > 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("pvclust-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('pvclust') > > assign(".oldSearch", search(), env = .CheckExEnv) > assign(".oldNS", loadedNamespaces(), env = .CheckExEnv) > cleanEx(); ..nameEx <- "lung" > > ### * lung > > flush(stderr()); flush(stdout()) > > ### Name: lung > ### Title: DNA Microarray Data of Lung Tumors > ### Aliases: lung > ### Keywords: datasets > > ### ** Examples > > ## Reading the data > data(lung) > > ## Multiscale Bootstrap Resampling > lung.pv <- pvclust(lung, nboot=100) Bootstrap (r = 0.5)... Done. Bootstrap (r = 0.6)... Done. Bootstrap (r = 0.7)... Done. Bootstrap (r = 0.8)... Done. Bootstrap (r = 0.9)... Done. Bootstrap (r = 1.0)... Done. Bootstrap (r = 1.1)... Done. Bootstrap (r = 1.2)... Done. Bootstrap (r = 1.3)... Done. Bootstrap (r = 1.4)... Done. > > ## CAUTION: nboot=100 may be too small for actual use. > ## We suggest nboot=1000 or larger. > ## plot/print functions will be useful for diagnostics. > > ## Plot the result > plot(lung.pv, cex=0.8, cex.pv=0.7) > > ask.bak <- par()$ask > par(ask=TRUE) > > pvrect(lung.pv, alpha=0.9) > msplot(lung.pv, edges=c(51,62,68,71)) Hit to see next plot: > par(ask=ask.bak) > > ## Print a cluster with high p-value > lung.pp <- pvpick(lung.pv, alpha=0.9) > lung.pp$clusters[[2]] [1] "306-99_Adeno" "306-99_node" > > ## Print its edge number > lung.pp$edges[2] [1] 4 > > ## We recommend parallel computing for large dataset as this one > ## Not run: > ##D library(snow) > ##D cl <- makeCluster(10, type="MPI") > ##D lung.pv <- parPvclust(cl, lung, nboot=1000) > ## End(Not run) > > > > graphics::par(get("par.postscript", env = .CheckExEnv)) > cleanEx(); ..nameEx <- "pvclust" > > ### * pvclust > > flush(stderr()); flush(stdout()) > > ### Name: pvclust > ### Title: Calculating P-values for Hierchical Clustering > ### Aliases: pvclust parPvclust > ### Keywords: cluster > > ### ** Examples > > ## using Boston data in package MASS > library(MASS) > data(Boston) > > ## multiscale bootstrap resampling > boston.pv <- pvclust(Boston, nboot=100) Bootstrap (r = 0.5)... Done. Bootstrap (r = 0.6)... Done. Bootstrap (r = 0.7)... Done. Bootstrap (r = 0.8)... Done. Bootstrap (r = 0.9)... Done. Bootstrap (r = 1.0)... Done. Bootstrap (r = 1.1)... Done. Bootstrap (r = 1.2)... Done. Bootstrap (r = 1.3)... Done. Bootstrap (r = 1.4)... Done. > > ## CAUTION: nboot=100 may be too small for actual use. > ## We suggest nboot=1000 or larger. > ## plot/print functions will be useful for diagnostics. > > ## plot dendrogram with p-values > plot(boston.pv) > > ask.bak <- par()$ask > par(ask=TRUE) > > ## highlight clusters with high au p-values > pvrect(boston.pv) > > ## print the result of multiscale bootstrap resampling > print(boston.pv, digits=3) Cluster method: average Distance : correlation Estimates on edges: au bp se.au se.bp v c pchi 1 1.000 1.000 0.000 0.000 0.000 0.000 0.000 2 0.966 0.942 0.027 0.008 -1.702 0.128 0.410 3 1.000 1.000 0.000 0.000 0.000 0.000 0.000 4 0.793 0.924 0.103 0.009 -1.126 -0.310 0.703 5 1.000 1.000 0.000 0.000 0.000 0.000 0.000 6 0.760 0.682 0.079 0.016 -0.589 0.116 0.126 7 0.755 0.735 0.082 0.015 -0.659 0.030 0.220 8 1.000 0.998 0.005 0.005 -3.086 0.213 1.000 9 1.000 1.000 0.000 0.000 0.000 0.000 0.000 10 0.829 0.871 0.078 0.011 -1.042 -0.091 0.964 11 0.940 0.986 0.094 0.006 -1.873 -0.321 0.947 12 0.774 0.803 0.083 0.013 -0.801 -0.050 0.914 13 1.000 1.000 0.000 0.000 0.000 0.000 0.000 > > ## plot diagnostic for curve fitting > msplot(boston.pv, edges=c(2,4,6,7)) Hit to see next plot: > par(ask=ask.bak) > > ## Print clusters with high p-values > boston.pp <- pvpick(boston.pv) > boston.pp $clusters $clusters[[1]] [1] "rm" "medv" $clusters[[2]] [1] "zn" "dis" $clusters[[3]] [1] "crim" "indus" "nox" "age" "rad" "tax" "ptratio" [8] "lstat" $edges [1] 3 5 9 > > ## Not run: > ##D ## parallel computation via snow package > ##D library(snow) > ##D cl <- makeCluster(10, type="MPI") > ##D > ##D ## parallel version of pvclust > ##D boston.pv <- parPvclust(cl, Boston, nboot=1000) > ## End(Not run) > > > > graphics::par(get("par.postscript", env = .CheckExEnv)) > ### *