search.trend: Searching for evolutionary trends in phenotypes and rates
In pasraia/RRphylo: Phylogenetic Ridge Regression Methods for Comparative Studies

Description Usage Arguments Details Value Author(s) References See Also Examples

This function searches for evolutionary trends in the phenotypic mean and the evolutionary rates for the entire tree and individual clades.

1 2	search.trend(RR,y,x1=NULL,x1.residuals = FALSE, node=NULL,cov=NULL,nsim=100,clus=0.5,ConfInt=FALSE,foldername=NULL,filename)

`RR`	an object produced by `RRphylo`.
`y`	the named vector (or matrix if multivariate) of phenotypes.
`x1`	the additional predictor to be specified if the RR object has been created using an additional predictor (i.e. multiple version of `RRphylo`). `'x1'` vector must be as long as the number of nodes plus the number of tips of the tree, which can be obtained by running `RRphylo` on the predictor as well, and taking the vector of ancestral states and tip values to form the `x1`. Note: only one predictor at once can be specified.
`x1.residuals`	logical specifying whether the residuals of regression between `y` and `x1` should be inspected for a phenotypic trend (see details and examples below). Default is `FALSE`.
`node`	the node number of individual clades to be specifically tested and contrasted to each other. It is `NULL` by default. Notice the node number must refer to the dichotomic version of the original tree, as produced by `RRphylo`.
`cov`	the covariate values to be specified if the RR object has been created using a covariate for rates calculation. As for `RRphylo`, `'cov'` must be as long as the number of nodes plus the number of tips of the tree, which can be obtained by running `RRphylo` on the covariate as well, and taking the vector of ancestral states and tip values to form the covariate (see the example below).
`nsim`	number of simulations to be performed. It is set at 100 by default.
`clus`	the proportion of clusters to be used in parallel computing. To run the single-threaded version of `search.trend` set `clus` = 0.
`ConfInt`	if `TRUE`, the function returns 95% confidence intervals around phenotypes and rates produced according to the Brownian motion model of evolution. It is `FALSE` by default.
`foldername`	has been deprecated; please see the argument `filename` instead.
`filename`	a character indicating the name of the pdf file and the path where it is to be saved. If no path is indicated the file is stored in the working directory

The function simultaneously returns the regression of phenotypes and phenotypic evolutionary rates against age tested against Brownian motion simulations to assess significance. To this aim rates are rescaled in the 0-1 range and then logged. The function stores the rates (both rescaled and unscaled absolute values) versus age regression and the phenotype versus age regression plots as .pdf files. In the plots, the 95% confidence intervals of phenotypes and rates simulated under the Brownian motion for each node are plotted as shaded areas. Regression lines are printed for all regressions. To assess significance, slopes are compared to a family of simulated slopes (BMslopes, where the number of simulations is equal to nsim), generated under the Brownian motion, using the fastBM function in the package phytools. Individual nodes are compared to the rest of the tree in different ways depending on whether phenotypes or rates (always unscaled in this case) versus age regressions are tested. With the former, the regression slopes for individual clades and the slope difference between clades is contrasted to slopes obtained through Brownian motion simulations. For the latter, regression models are tested and contrasted to each other referring to estimated marginal means, by using the emmeans function in the package emmeans.

The multiple regression version of RRphylo allows to incorporate the effect of an additional predictor in the computation of evolutionary rates without altering the ancestral character estimation. Thus, when a multiple RRphylo output is fed to search.trend, the predictor effect is accounted for on the absolute evolutionary rates, but not on the phenotype. However, in some situations the user might want to factor out the predictor effect on phenotypes as well. Under the latter circumstance, by setting the argument x1.residuals = TRUE, the residuals of the response to predictor regression are used as to represent the phenotype.

The function returns a list object containing:

$trends.data a 'RRphyloList' object including:

$phenotypeVStime: a data frame of phenotypic values (or y versus x1 regression residuals if x1.residuals=TRUE) and their distance from the tree root for each node (i.e. ancestral states) and tip of the tree.
$absrateVStime: a data frame of RRphylo rates and the distance from the tree root (age). If y is multivariate, it also includes the multiple rates for each y vector. If node is specified, each branch is classified as belonging to an indicated clade.
$rescaledrateVStime: a data frame of rescaled RRphylo rates and the distance from the tree root (age). If y is multivariate, it also includes the multiple rates for each y vector. If node is specified, each branch is classified as belonging to an indicated clade. NAs correspond either to very small values or to outliers which are excluded from the analysis.

$phenotypic.regression results of phenotype (y versus x1 regression residuals) versus age regression. It reports a p-value for the regression slope between the variables (p.real), a p-value computed contrasting the real slope to Brownian motion simulations (p.random), and a parameter indicating the deviation of the phenotypic mean from the root value in terms of the number of standard deviations of the trait distribution (dev). dev is 0 under Brownian Motion. Only p.random should be inspected to assess significance.

$rate.regression results of the rates (rescaled absolute values) versus age regression. It reports a p-value for the regression between the variables (p.real), a p-value computed contrasting the real slope to Brownian motion simulations (p.random), and a parameter indicating the ratio between the range of phenotypic values and the range of such values halfway along the tree height, divided to the same figure under Brownian motion (spread). spread is 1 under Brownian Motion. Only p.random should be inspected to assess significance.

$ConfInts a 'RRphyloList' object including the 95% confidence intervals around phenotypes and rates (both rescaled and unscaled absolute rates) produced according to the Brownian motion model of evolution.

If specified, individual nodes are tested as the whole tree, the results are summarized in the objects:

$node.phenotypic.regression results of phenotype (or y versus x1 regression residuals) versus age regression through node. It reports the slope for the regression between the variables at node (slope), a p-value computed contrasting the real slope to Brownian motion simulations (p.random), the difference between estimated marginal means predictions for the group and for the rest of the tree (emm.difference), and a p-value for the emm.difference (p.emm).

$node.rate.regression results of the rates (absolute values) versus age regression through node. It reports the difference between estimated marginal means predictions for the group and for the rest of the tree (emm.difference), a p-value for the emm.difference (p.emm), the difference between regression slopes for the group and for the rest of the tree (slope.difference), and a p-value for the slope.difference (p.slope).

If more than one node is specified, the object $group.comparison reports the same results as $node.phenotypic.regression and $node.rate.regression obtained by comparing individual clades to each other.

Silvia Castiglione, Carmela Serio, Pasquale Raia, Alessandro Mondanaro, Marina Melchionna, Mirko Di Febbraro, Antonio Profico, Francesco Carotenuto

Castiglione, S., Serio, C., Mondanaro, A., Di Febbraro, M., Profico, A., Girardi, G., & Raia, P. (2019) Simultaneous detection of macroevolutionary patterns in phenotypic means and rate of change with and within phylogenetic trees including extinct species. PLoS ONE, 14: e0210101. https://doi.org/10.1371/journal.pone.0210101

search.trend vignette

 ## Not run: 
data("DataOrnithodirans")
DataOrnithodirans$treedino->treedino
DataOrnithodirans$massdino->massdino
cc<- 2/parallel::detectCores()

# Extract Pterosaurs tree and data
library(ape)
extract.clade(treedino,746)->treeptero
massdino[match(treeptero$tip.label,names(massdino))]->massptero
massptero[match(treeptero$tip.label,names(massptero))]->massptero

# Case 1. "RRphylo" whitout accounting for the effect of a covariate
RRphylo(tree=treeptero,y=log(massptero))->RRptero

# Case 1.1. "search.trend" whitout indicating nodes to be tested for trends
search.trend(RR=RRptero, y=log(massptero), nsim=100, clus=cc,
             filename=paste(tempdir(), "ST", sep="/"),cov=NULL,ConfInt=FALSE,node=NULL)

# Case 1.2. "search.trend" indicating nodes to be specifically tested for trends
search.trend(RR=RRptero, y=log(massptero), nsim=100, node=143, clus=cc,
             filename=paste(tempdir(), "STnode", sep="/"),cov=NULL,ConfInt=FALSE)


# Case 2. "RRphylo" accounting for the effect of a covariate
# "RRphylo" on the covariate in order to retrieve ancestral state values
RRphylo(tree=treeptero,y=log(massptero))->RRptero
c(RRptero$aces,log(massptero))->cov.values
names(cov.values)<-c(rownames(RRptero$aces),names(massptero))
RRphylo(tree=treeptero,y=log(massptero),cov=cov.values)->RRpteroCov

# Case 2.1. "search.trend" whitout indicating nodes to be tested for trends
search.trend(RR=RRpteroCov, y=log(massptero), nsim=100, clus=cc,
             filename=paste(tempdir(), "ST_cov", sep="/"),ConfInt=FALSE,cov=cov.values)

# Case 2.2. "search.trend" indicating nodes to be specifically tested for trends
search.trend(RR=RRpteroCov, y=log(massptero), nsim=100, node=143, clus=cc,
             filename=paste(tempdir(), "STnode_cov", sep="/"),ConfInt=FALSE,cov=cov.values)


# Case 3. "search.trend" on multiple "RRphylo"
data("DataCetaceans")
DataCetaceans$treecet->treecet
DataCetaceans$masscet->masscet
DataCetaceans$brainmasscet->brainmasscet
DataCetaceans$aceMyst->aceMyst

drop.tip(treecet,treecet$tip.label[-match(names(brainmasscet),treecet$tip.label)])->treecet.multi
masscet[match(treecet.multi$tip.label,names(masscet))]->masscet.multi

RRphylo(tree=treecet.multi,y=masscet.multi)->RRmass.multi
RRmass.multi$aces[,1]->acemass.multi
c(acemass.multi,masscet.multi)->x1.mass

RRphylo(tree=treecet.multi,y=brainmasscet,x1=x1.mass)->RRmulti

# incorporating the effect of body size at inspecting trends in absolute evolutionary rates
search.trend(RR=RRmulti, y=brainmasscet,x1=x1.mass,clus=cc,
             filename=paste(tempdir(), "STmulti_rate", sep="/"))

# incorporating the effect of body size at inspecting trends in both absolute evolutionary
# rates and phenotypic values (by using brain versus body mass regression residuals)
search.trend(RR=RRmulti, y=brainmasscet,x1=x1.mass,x1.residuals=TRUE,clus=cc,
             filename=paste(tempdir(), "STmulti", sep="/"))
   
## End(Not run)