diff --git a/DESCRIPTION b/DESCRIPTION index c0ccd5e..17d8f05 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -54,9 +54,9 @@ Imports: future.apply Depends: R (>= 2.10) -URL: https://github.com/maribraga/evolnets, - https://maribraga.github.io/evolnets/ -BugReports: https://github.com/maribraga/evolnets/issues +URL: https://github.com/evonetslab/evolnets, + https://evonetslab.github.io/evolnets/ +BugReports: https://github.com/evonetslab/evolnets/issues VignetteBuilder: knitr Config/testthat/edition: 3 Config/Needs/website: rmarkdown diff --git a/README.Rmd b/README.Rmd index 59dc7ae..a923d4b 100644 --- a/README.Rmd +++ b/README.Rmd @@ -35,7 +35,7 @@ if(!require("devtools", quietly = TRUE)) { library(devtools) } -devtools::install_github("maribraga/evolnets") +devtools::install_github("evonetslab/evolnets") ``` ## About @@ -55,4 +55,4 @@ rates, ancestral states and samples. `NODF_posterior_at_ages( )`. -See the full documentation at the [evolnets' website](https://maribraga.github.io/evolnets/). +See the full documentation at the [evolnets' website](https://evonetslab.github.io/evolnets/). diff --git a/README.md b/README.md index 37cc2bb..0afd1bb 100644 --- a/README.md +++ b/README.md @@ -24,7 +24,7 @@ if(!require("devtools", quietly = TRUE)) { library(devtools) } -devtools::install_github("maribraga/evolnets") +devtools::install_github("evonetslab/evolnets") ``` ## About @@ -46,4 +46,4 @@ rates, ancestral states and samples. `Q_posterior_at_ages( )`, `NODF_posterior_at_ages( )`. See the full documentation at the [evolnets’ -website](https://maribraga.github.io/evolnets/). +website](https://evonetslab.github.io/evolnets/). diff --git a/_pkgdown.yml b/_pkgdown.yml index 05f2a97..ad820c0 100644 --- a/_pkgdown.yml +++ b/_pkgdown.yml @@ -1,4 +1,4 @@ -url: https://maribraga.github.io/evolnets/ +url: https://evonetslab.github.io/evolnets/ template: bootstrap: 5 diff --git a/man/evolnets-package.Rd b/man/evolnets-package.Rd index 2baf307..c358aca 100644 --- a/man/evolnets-package.Rd +++ b/man/evolnets-package.Rd @@ -11,9 +11,9 @@ RevBayes offers models to infer host-repertoire evolution, but no tools to parse \seealso{ Useful links: \itemize{ - \item \url{https://github.com/maribraga/evolnets} - \item \url{https://maribraga.github.io/evolnets/} - \item Report bugs at \url{https://github.com/maribraga/evolnets/issues} + \item \url{https://github.com/evonetslab/evolnets} + \item \url{https://evonetslab.github.io/evolnets/} + \item Report bugs at \url{https://github.com/evonetslab/evolnets/issues} } } diff --git a/vignettes/articles/evolnets.Rmd b/vignettes/articles/evolnets.Rmd index a8c8530..00e0983 100644 --- a/vignettes/articles/evolnets.Rmd +++ b/vignettes/articles/evolnets.Rmd @@ -40,7 +40,6 @@ We will need a few packages: ```{r setup, message = FALSE} library(evolnets) library(ape) -library(treeio) library(dplyr) library(ggplot2) ``` @@ -52,45 +51,35 @@ and four pieces of data: - `history.txt` output from RevBayes, which contains the sampled histories of host-repertoire evolution during MCMC - Matrix of extant interactions for plotting -`evolnets` includes example data which we'll use in this vignette. But when you use your own data, here is how you read them into R. The evolnets function `read_history` is designed specifically to read .txt files produced by RevBayes. As more inference methods become available, we will expand the scope of this function. - -```{r, eval=FALSE} -history <- read_history('history.txt', burnin = 0.1) -``` +`evolnets` includes example data which we'll use in this vignette. But when you use your own data, simply replace the `system.file()` call with the path for your files. Your phylogenetic trees need to be `phylo` objects, and many tree files can be read into R using the `ape` or `treeio` packages. One important step when reading in the symbiont tree is to retrieve node labels from RevBayes. R and RevBayes label the internal nodes of a tree in different ways, so it's important to export the symbiont tree from RevBayes to be able to use them in R (using the `write()` function in RevBayes). Here's how you do it: -```{r, eval=FALSE} -# we don't need node labels for the host tree, so you can read it using ape::read.tree() -host_tree <- read.tree('host_tree.tre') - -# read the symbiont tree exported from RevBayes with treeio::read.beast.newick() -treeRev <- read.beast.newick("symbiont_tree.tre") - -# treeRev contains the tree -tree <- treeRev@phylo +```{r} +# read the symbiont tree exported from RevBayes +tree <- read_tree_from_revbayes(system.file("extdata", "tree_pieridae.tre", package = "evolnets")) -# and more information, including -# a data table with the node order from RevBayes ($index) and R ($node) -index_node <- treeRev@data %>% - mutate(node = as.numeric(node)) %>% - arrange(node) +# we don't need node labels for the host tree, so you can read it using ape::read.tree() +host_tree <- ape::read.tree(system.file("extdata", "host_tree_pieridae.phy", package = "evolnets")) +``` -# use node indices from RevBayes as node labels -indices <- index_node %>% - filter(node > Ntip(tree)) %>% - pull(index) -names(indices) <- NULL +The evolnets function `read_history` is designed specifically to read .txt files produced by RevBayes. As more inference methods become available, we will expand the scope of this function. -tree$node.label <- paste0("Index_",indices) +```{r} +history <- read_history(system.file("extdata", "history_thin_pieridae.txt", package = "evolnets")) ``` -As I said before, in this vignette we'll use the data that comes with evolnets. This data comes from a [paper](https://onlinelibrary.wiley.com/doi/10.1111/ele.13842) where we studied the evolution of interactions between Pieridae butterflies and their Angiosperm host plants. The results are not identical to the ones in the paper because we will use a subset of the MCMC samples in this example in order to reduce file size and speed up the analysis. +Lastly, we need a matrix of extant interations. -```{r data} -data(tree, host_tree, history, extant_net) +```{r} +matrix <- read.csv( + system.file("extdata", "interaction_matrix_pieridae.csv", package = "evolnets"), + row.names = 1) |> + as.matrix() ``` +As I said before, in this vignette we'll use data that comes with evolnets. This data comes from a [paper](https://onlinelibrary.wiley.com/doi/10.1111/ele.13842) where we studied the evolution of interactions between Pieridae butterflies and their Angiosperm host plants. The results are not identical to the ones in the paper because we will use a subset of the MCMC samples in this example in order to reduce file size and speed up the analysis. + Now we can start extracting information from the inferred history. ## Rate of evolution @@ -102,13 +91,25 @@ Let's calculate the average number of events (host gains and host losses) across (gl_events <- count_gl(history)) ``` -We estimated that 151 events happened across the diversification of Pieridae, being 75 host gains and 76 host losses. Similarly, we can calculate the rate of host-repertoire evolution across the branches of the symbiont tree, which is the number of events divided by the sum of branch lengths of the symbiont tree. In this case, we inferred that the rate of evolution is around 6 events every 100 million years, along each branch of the Pieridae tree. +```{r, echo = FALSE} +tot_events <- round(n_events$mean, digits = 1) +gains <- gl_events[1] +losses <- gl_events[2] +``` + +We estimated that, on average, `r tot_events` events happened across the diversification of Pieridae, being `r gains` host gains and `r losses` host losses. Similarly, we can calculate the rate of host-repertoire evolution across the branches of the symbiont tree, which is the number of events divided by the sum of branch lengths of the symbiont tree. ```{r} (rate <- effective_rate(history,tree)) (gl_rates <- rate_gl(history, tree)) ``` +```{r, echo = FALSE} +rate_m <- round(rate$mean, digits = 2) * 100 +``` + +In this case, we inferred that the rate of evolution is around `r rate_m` events every 100 million years, along each branch of the Pieridae tree. + ## Ancestral states at internal nodes In line with a more traditional approach of ancestral state reconstruction, we can extract the posterior probability for each interaction between ancestral symbiont species and all hosts. First we need to choose which internal nodes we want to include (default is all nodes). For that, we need to look at the labels at the internal nodes of the tree file exported by RevBayes. @@ -142,7 +143,7 @@ at_nodes <- posterior_at_nodes(history, tree, host_tree) ``` ```{r, fig.width=7, fig.height=7} -p <- plot_module_matrix2(extant_net, at_nodes, tree, host_tree) +p <- plot_matrix_phylo(matrix, at_nodes, tree, host_tree) # adjust text size in the second panel p[[2]] <- p[[2]] + theme(axis.text = element_text(size = 4)) @@ -164,31 +165,37 @@ The main goal of `evolnets` is to reconstruct the evolution of networks. We can ```{r} ages <- c(60, 50, 40, 0) at_ages <- posterior_at_ages(history, ages, tree, host_tree) -pp_at_ages <- at_ages[[2]] ``` -Then, we can make different interaction matrices based on two things: the minimum posterior probability for an interaction to be included in the network and whether we want a binary or a weighted network. +Then, we can make different interaction matrices based on two things: the minimum posterior probability for an interaction to be included in the network, and whether we want a binary or a weighted network. ```{r} -weighted_nets_50 <- get_summary_network(pp_at_ages, ages, pt = 0.5, weighted = TRUE) -binary_nets_90 <- get_summary_network(pp_at_ages, ages, pt = 0.9, weighted = FALSE) +weighted_nets_50 <- get_summary_networks(at_ages, threshold = 0.5, weighted = TRUE) +binary_nets_90 <- get_summary_networks(at_ages, threshold = 0.9, weighted = FALSE) + +#plot + ``` + + + + ### Sampled networks ```{r} -samples_at_ages <- at_ages$samples +samples_at_ages <- get_sampled_networks(at_ages) ``` Then, we can calculate network properties for each sampled network. ```{r} # calculate posterior distribution of nestedness -Nz <- index_at_ages(samples_at_ages, index = "NODF", nnull = 10) -``` +Nz <- index_at_ages_samples(samples_at_ages, index = "NODF", nnull = 10) + +#plot + -```{r, eval=FALSE} -# calculate posterior distribution of modularity (extremely slow) -Qz <- index_at_ages(samples_at_ages, index = "Q", nnull = 10) ``` +You can do the same for modularity, but the algorithm is much slower.