Either a bammdata object or a credibleshiftset
object.
expectedNumberOfShifts
The expected number of shifts under the
prior.
threshold
The marginal posterior-to-prior odds ratio used as a
cutoff for distinguishing between "core" and "non-core" rate shifts.
Details
This function estimates the rate shift configuration with the
highest maximum a posteriori (MAP) probability. It returns a
bammdata object with a single sample. This can be plotted with
plot.bammdata, and individual rate shifts can then
be added with addBAMMshifts.
The parameters of this object are averaged over all samples in the
posterior that were assignable to the MAP shift configuration. All
non-core shifts have been excluded, such that the only shift
information contained in the object is from the "significant" rate
shifts, as determined by the relevant marginal posterior-to-prior odds
ratio threshold.
You can extract the same information from the credible set of shift
configurations. See credibleShiftSet for more
information.
Value
A class bammdata object with a single sample, corresponding
to the diversification rate shift configuration with the maximum a
posteriori probability. See getEventData for details.
data(whales, events.whales)
ed <- getEventData(whales, events.whales, burnin=0.1, nsamples=500)
# Get prior distribution on shifts-per-branch:
bp <- getBranchShiftPriors(whales, expectedNumberOfShifts = 1)
# Pass the event data object in to the function:
best <- getBestShiftConfiguration(ed, expectedNumberOfShifts = 1,
threshold = 5)
plot(best, lwd=2)
addBAMMshifts(best, cex=2)
# Now we can also work with the credible shift set:
css <- credibleShiftSet(ed, expectedNumberOfShifts = 1, threshold = 5)
summary(css)
# examine model-averaged shifts from MAP configuration-
# This gives us parameters, times, and associated nodes
# of each evolutionary rate regime (note that one of
# them corresponds to the root)
css$eventData[[1]];
# Get bammdata representation of MAP configuration:
best <- getBestShiftConfiguration(css, expectedNumberOfShifts = 1,
threshold = 5)
plot(best)
addBAMMshifts(best)
Results
R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)
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 an HTML browser interface to help.
Type 'q()' to quit R.
> library(BAMMtools)
Loading required package: ape
> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/BAMMtools/getBestShiftConfiguration.Rd_%03d_medium.png", width=480, height=480)
> ### Name: getBestShiftConfiguration
> ### Title: Get the best (sampled) rate shift configuration from a 'BAMM'
> ### analysis
> ### Aliases: getBestShiftConfiguration
> ### Keywords: models
>
> ### ** Examples
>
> data(whales, events.whales)
> ed <- getEventData(whales, events.whales, burnin=0.1, nsamples=500)
Processing event data from data.frame
Discarded as burnin: GENERATIONS < 995000
Analyzing 500 samples from posterior
Setting recursive sequence on tree...
Done with recursive sequence
>
> # Get prior distribution on shifts-per-branch:
> bp <- getBranchShiftPriors(whales, expectedNumberOfShifts = 1)
>
> # Pass the event data object in to the function:
> best <- getBestShiftConfiguration(ed, expectedNumberOfShifts = 1,
+ threshold = 5)
Processing event data from data.frame
Discarded as burnin: GENERATIONS < 0
Analyzing 1 samples from posterior
Setting recursive sequence on tree...
Done with recursive sequence
> plot(best, lwd=2)
> addBAMMshifts(best, cex=2)
>
> # Now we can also work with the credible shift set:
> css <- credibleShiftSet(ed, expectedNumberOfShifts = 1, threshold = 5)
>
> summary(css)
95 % credible set of rate shift configurations sampled with BAMM
Distinct shift configurations in credible set: 4
Frequency of 4 shift configurations with highest posterior probability:
rank probability cumulative Core_shifts
1 0.462 0.462 1
2 0.286 0.748 1
3 0.138 0.886 1
4 0.102 0.988 0
>
> # examine model-averaged shifts from MAP configuration-
> # This gives us parameters, times, and associated nodes
> # of each evolutionary rate regime (note that one of
> # them corresponds to the root)
> css$eventData[[1]];
node time lam1 lam2 mu1 mu2 index
1 88 0.0000 0.114234 -0.00794353 0.0455479 0 1
2 141 26.0462 0.315679 -0.00467391 0.0434990 0 2
>
> # Get bammdata representation of MAP configuration:
> best <- getBestShiftConfiguration(css, expectedNumberOfShifts = 1,
+ threshold = 5)
Processing event data from data.frame
Discarded as burnin: GENERATIONS < 0
Analyzing 1 samples from posterior
Setting recursive sequence on tree...
Done with recursive sequence
>
> plot(best)
> addBAMMshifts(best)
>
>
>
>
>
> dev.off()
null device
1
>