From 1fc32e7cc32ce16ac2a021ee05b783ff7addebf8 Mon Sep 17 00:00:00 2001 From: Klaus Schliep Date: Mon, 8 Feb 2021 16:02:53 +0100 Subject: [PATCH] remove phangorn dependency and add upgma --- DESCRIPTION | 1 - NAMESPACE | 6 ++++-- R/bruvo.r | 5 ++--- R/upgma.R | 22 ++++++++++++++++++++++ man/bruvo.boot.Rd | 4 ++-- man/upgma.Rd | 33 +++++++++++++++++++++++++++++++++ vignettes/mlg.Rmd | 1 - 7 files changed, 63 insertions(+), 9 deletions(-) create mode 100644 R/upgma.R create mode 100644 man/upgma.Rd diff --git a/DESCRIPTION b/DESCRIPTION index 5d7813ef..ddebe1cf 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -46,7 +46,6 @@ Imports: utils, vegan, ggplot2, - phangorn, ape (>= 3.1-1), igraph (>= 1.0.0), methods, diff --git a/NAMESPACE b/NAMESPACE index 8705129e..0b5cd621 100755 --- a/NAMESPACE +++ b/NAMESPACE @@ -99,6 +99,7 @@ export(rrmlg) export(samp.ia) export(shufflepop) export(test_replen) +export(upgma) export(visible) export(win.ia) exportClasses(MLG) @@ -123,6 +124,8 @@ importFrom(ade4,is.euclid) importFrom(ade4,lingoes) importFrom(ade4,quasieuclid) importFrom(ape,add.scale.bar) +importFrom(ape,as.phylo) +importFrom(ape,as.phylo.hclust) importFrom(ape,axisPhylo) importFrom(ape,boot.phylo) importFrom(ape,is.ultrametric) @@ -172,8 +175,6 @@ importFrom(igraph,print.igraph) importFrom(magrittr,"%>%") importFrom(pegas,as.loci) importFrom(pegas,loci2genind) -importFrom(phangorn,midpoint) -importFrom(phangorn,upgma) importFrom(polysat,Genotypes) importFrom(polysat,Missing) importFrom(polysat,Ploidies) @@ -185,6 +186,7 @@ importFrom(stats,as.formula) importFrom(stats,dbinom) importFrom(stats,df) importFrom(stats,dist) +importFrom(stats,hclust) importFrom(stats,median) importFrom(stats,printCoefmat) importFrom(stats,quantile) diff --git a/R/bruvo.r b/R/bruvo.r index 6e1a3521..fa3defb1 100755 --- a/R/bruvo.r +++ b/R/bruvo.r @@ -297,7 +297,7 @@ bruvo.between <- function(query, ref, replen = 1, add = TRUE, loss = TRUE, by_lo #' desired. #' #' @param tree any function that can generate a tree from a distance matrix. -#' Default is \code{\link[phangorn]{upgma}}. +#' Default is \code{\link{upgma}}. #' #' @param showtree \code{logical} if \code{TRUE}, a tree will be plotted with #' nodelabels. @@ -331,7 +331,7 @@ bruvo.between <- function(query, ref, replen = 1, add = TRUE, loss = TRUE, by_lo #' replacement, recalculate the tree, and tally up the bootstrap support #' (measured in percent success). While this function can take any tree #' function, it has native support for two algorithms: \code{\link[ape]{nj}} -#' and \code{\link[phangorn]{upgma}}. If you want to use any other functions, +#' and \code{\link{upgma}}. If you want to use any other functions, #' you must load the package before you use them (see examples). #' #' @note \strong{Please refer to the documentation for bruvo.dist for details on @@ -383,7 +383,6 @@ bruvo.between <- function(query, ref, replen = 1, add = TRUE, loss = TRUE, by_lo #' } #' #==============================================================================# -#' @importFrom phangorn upgma midpoint #' @importFrom ape nodelabels nj boot.phylo plot.phylo axisPhylo ladderize #' @importFrom ape add.scale.bar nodelabels tiplabels is.ultrametric # / \ diff --git a/R/upgma.R b/R/upgma.R new file mode 100644 index 00000000..336e00f6 --- /dev/null +++ b/R/upgma.R @@ -0,0 +1,22 @@ +#' UPGMA +#' +#' UPGMA clustering. Just a wrapper function around \code{\link[stats]{hclust}}. +#' +#' @param d A distance matrix. +#' @return A phylogenetic tree of class \code{phylo}. +#' @author Klaus Schliep \email{klaus.schliep@@gmail.com} +#' @seealso \code{\link{hclust}}, \code{\link{as.phylo}} +#' @importFrom ape as.phylo.hclust as.phylo +#' @importFrom stats hclust as.dist +#' @keywords cluster +#' @examples +#' +#' library(ape) +#' data(woodmouse) +#' dm <- dist.dna(woodmouse) +#' tree <- upgma(dm) +#' plot(tree) +#' +#' @rdname upgma +#' @export +"upgma" <- function(d) as.phylo(hclust(as.dist(d), method = "average")) diff --git a/man/bruvo.boot.Rd b/man/bruvo.boot.Rd index f277baf1..fb3126f4 100755 --- a/man/bruvo.boot.Rd +++ b/man/bruvo.boot.Rd @@ -34,7 +34,7 @@ the genome loss model presented in Bruvo et al. 2004.} desired.} \item{tree}{any function that can generate a tree from a distance matrix. -Default is \code{\link[phangorn]{upgma}}.} +Default is \code{\link{upgma}}.} \item{showtree}{\code{logical} if \code{TRUE}, a tree will be plotted with nodelabels.} @@ -67,7 +67,7 @@ This function will calculate a tree based off of Bruvo's distance replacement, recalculate the tree, and tally up the bootstrap support (measured in percent success). While this function can take any tree function, it has native support for two algorithms: \code{\link[ape]{nj}} - and \code{\link[phangorn]{upgma}}. If you want to use any other functions, + and \code{\link{upgma}}. If you want to use any other functions, you must load the package before you use them (see examples). } \note{ diff --git a/man/upgma.Rd b/man/upgma.Rd new file mode 100644 index 00000000..5cda987d --- /dev/null +++ b/man/upgma.Rd @@ -0,0 +1,33 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/upgma.R +\name{upgma} +\alias{upgma} +\title{UPGMA} +\usage{ +upgma(d) +} +\arguments{ +\item{d}{A distance matrix.} +} +\value{ +A phylogenetic tree of class \code{phylo}. +} +\description{ +UPGMA clustering. Just a wrapper function around \code{\link[stats]{hclust}}. +} +\examples{ + +library(ape) +data(woodmouse) +dm <- dist.dna(woodmouse) +tree <- upgma(dm) +plot(tree) + +} +\seealso{ +\code{\link{hclust}}, \code{\link{as.phylo}} +} +\author{ +Klaus Schliep \email{klaus.schliep@gmail.com} +} +\keyword{cluster} diff --git a/vignettes/mlg.Rmd b/vignettes/mlg.Rmd index ed655308..ab123d04 100644 --- a/vignettes/mlg.Rmd +++ b/vignettes/mlg.Rmd @@ -197,7 +197,6 @@ them as different. If we calculate the pairwise euclidean distances between samples, we see that "new", "mut" and, "C" are very similar to each other: ```{r, fig.width = 5, fig.height = 5} -library("phangorn") library("ape") raw_dist <- function(x){ dist(genind2df(x, usepop = FALSE))