bruvo.boot: Create a tree using Bruvo's Distance with non-parametric...

View source: R/bruvo.r

bruvo.bootR Documentation

Create a tree using Bruvo's Distance with non-parametric bootstrapping.

Description

Create a tree using Bruvo's Distance with non-parametric bootstrapping.

Usage

bruvo.boot(
  pop,
  replen = 1,
  add = TRUE,
  loss = TRUE,
  sample = 100,
  tree = "upgma",
  showtree = TRUE,
  cutoff = NULL,
  quiet = FALSE,
  root = NULL,
  ...
)

Arguments

pop

a genind or genclone object

replen

a vector of integers indicating the length of the nucleotide repeats for each microsatellite locus.

add

if TRUE, genotypes with zero values will be treated under the genome addition model presented in Bruvo et al. 2004.

loss

if TRUE, genotypes with zero values will be treated under the genome loss model presented in Bruvo et al. 2004.

sample

an integer indicated the number of bootstrap replicates desired.

tree

any function that can generate a tree from a distance matrix. Default is upgma.

showtree

logical if TRUE, a tree will be plotted with nodelabels.

cutoff

integer the cutoff value for bootstrap node label values (between 0 and 100).

quiet

logical defaults to FALSE. If TRUE, a progress bar and messages will be suppressed.

root

logical This is a parameter passed on to boot.phylo. If the tree argument produces a rooted tree (e.g. "upgma"), then this value should be TRUE. If it produces an unrooted tree (e.g. "nj"), then the value should be FALSE. By default, it is set to NULL, which will assume an unrooted phylogeny unless the function name contains "upgma".

...

any argument to be passed on to boot.phylo. eg. quiet = TRUE.

Details

This function will calculate a tree based off of Bruvo's distance and then utilize boot.phylo to randomly sample loci with 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: nj and upgma. If you want to use any other functions, you must load the package before you use them (see examples).

Value

a tree of class phylo with nodelables

Note

Please refer to the documentation for bruvo.dist for details on the algorithm. If the user does not provide a vector of appropriate length for replen , it will be estimated by taking the minimum difference among represented alleles at each locus. IT IS NOT RECOMMENDED TO RELY ON THIS ESTIMATION.

Author(s)

Zhian N. Kamvar, Javier F. Tabima

References

Ruzica Bruvo, Nicolaas K. Michiels, Thomas G. D'Souza, and Hinrich Schulenburg. A simple method for the calculation of microsatellite genotype distances irrespective of ploidy level. Molecular Ecology, 13(7):2101-2106, 2004.

See Also

bruvo.dist, nancycats, upgma, nj, boot.phylo, nodelabels, tab, missingno.

Examples

# Please note that the data presented is assuming that the nancycat dataset 
# contains all dinucleotide repeats, it most likely is not an accurate
# representation of the data.

# Load the nancycats dataset and construct the repeat vector.
data(nancycats)
ssr <- rep(2, 9)

# Analyze the 1st population in nancycats

bruvo.boot(popsub(nancycats, 1), replen = ssr)

## Not run: 

# Always load the library before you specify the function.
library("ape")

# Estimate the tree based off of the BIONJ algorithm.

bruvo.boot(popsub(nancycats, 9), replen = ssr, tree = bionj)

# Utilizing  balanced FastME
bruvo.boot(popsub(nancycats, 9), replen = ssr, tree = fastme.bal)

# To change parameters for the tree, wrap it in a function.
# For example, let's build the tree without utilizing subtree-prune-regraft

myFastME <- function(x) fastme.bal(x, nni = TRUE, spr = FALSE, tbr = TRUE)
bruvo.boot(popsub(nancycats, 9), replen = ssr, tree = myFastME)


## End(Not run)


poppr documentation built on May 29, 2024, 5:54 a.m.