Data collection and phylogeny

Toxin expression data were collected from 39 publications (Supplementary table 1). The selection criteria for toxins and taxa are the same as in (Barua and Mikheyev 2019). We scaled gene expression levels by the average within-species variance allowing us to standardize the measurements and carry out comparisons across toxins. This scaled dataset was used in all subsequent analysis. We used a previously described time-calibrated phylogeny of squamates based on two large datasets comprising of 44 nuclear genes for 161 squamates, and a dataset of 12 genes from 4161 squamate species, both these datasets represented families and subfamilies (Wiens et al. 2012; Pyron et al. 2013; Zheng and Wiens 2016; Barua and Mikheyev 2019).

library(tidyverse)
library(ggplot2)
library(reshape2)
library(geiger)
library(spptest)
library(BAMMtools)
library(ape)
library(phytools)
library(coda)
library(kableExtra)

Check tree

tree <- read.newick("./Final_tree.tre")
plotTree(tree)

Check if tree is ultrametric

is.ultrametric(tree)
## [1] TRUE


Load data

Toxin expression data was collected from 39 publications (online supp). Out of the 25 reported toxin families we selected only 10 as they were the most ubiquitous toxins amongst all snakes. We restricted our dataset to include components that are found in at least 50% of snakes and eliminated low-abundance toxin families (Supplementary Fig 1 & 2). Incidentally, this cutoff also eliminated many low-abundance toxin families (on average less than 1% of the venom). The abundance of these toxins would be more difficult to estimate, as they are closer to the signal to noise to threshold of gene expression experiments. Overall 10 out of 25 toxin families we retained. Toxins levels were recorded as per publication. Toxin values reported as absolute FPKM values, were converted to a percentages of total toxin transcript expression.

dat<-read.csv("./data/scaleddata/vPhen_data.csv", header = T)

Supplementary Table 1: Gene expression values and references

species TFTx BPP.CNP CRISP CTL GF KSPI LAAO ePLA2 vPLA2 SVMP SVSP Total cm Family Author Tech Ref
Agkistrodon_piscivorus 0.00000 0.04570 0.00000 0.06080 0.00000 0.00000 0.00000 0.00000 0.38960 0.34900 0.08300 0.92810 127.0 Viperidae Aronow (2018) Illumina https://search.proquest.com/openview/2997dede40d17e51e7dabfd03d9762df/1?pq-origsite=gscholar&cbl=18750&diss=y
Atropoides_picadoi 0.00400 0.16400 0.02600 0.04800 0.00600 0.00010 0.03100 0.00000 0.02300 0.56100 0.10900 0.97210 120.2 Viperidae Durban et al. (2011) Roche GS FLX https://doi.org/10.1186/1471-2164-12-259
Bitis_gabonica 0.00000 0.01813 0.00000 0.07591 0.03626 0.09064 0.05665 0.00000 0.05665 0.37729 0.14729 0.85882 120.0 Viperidae Francischetti et al. (2004) Beckman Coulter https://doi.org/10.1016/j.gene.2004.03.024
Boiga_irregularis 0.67475 0.01400 0.03835 0.02271 0.00025 0.00100 0.00000 0.00040 0.00000 0.24854 0.00000 1.00000 270.0 Colubridae McGivern et al. (2014) Illumina https://doi.org/10.1186/1471-2164-15-1061
Bothriechis_lateralis 0.01600 0.23600 0.03200 0.03100 0.01800 0.00020 0.03400 0.00000 0.03000 0.44900 0.10400 0.95020 90.0 Viperidae Durban et al. (2011) Roche GS FLX https://doi.org/10.1186/1471-2164-12-259
Bothriechis_schlegelii 0.03600 0.15300 0.03600 0.00200 0.04900 0.00120 0.03200 0.00000 0.09200 0.26300 0.20500 0.86920 55.0 Viperidae Durban et al. (2011) Roche GS FLX https://doi.org/10.1186/1471-2164-12-259
Bothrops_alternatus 0.00000 0.11800 0.02200 0.16500 0.02500 0.00000 0.01800 0.00000 0.01400 0.58500 0.05100 0.99800 100.0 Viperidae Pereira de Paula et al. (2014) MegaBACE (GE) https://doi.org/10.1016/j.cbd.2014.09.001
Bothrops_asper 0.00800 0.15000 0.00800 0.01000 0.00000 0.00030 0.04200 0.00000 0.17800 0.41600 0.13400 0.94630 150.0 Viperidae Durban et al. (2011) Roche GS FLX https://doi.org/10.1186/1471-2164-12-259
Bothrops_atrox 0.00000 0.07100 0.00500 0.04800 0.01000 0.00000 0.01400 0.00000 0.14600 0.61600 0.08200 0.99200 100.0 Viperidae Neiva et al. (2009) MegaBACE (GE) https://doi.org/10.1016/j.toxicon.2009.01.006
Bothrops_colombiensis 0.00000 0.01000 0.02000 0.04000 0.02000 0.00000 0.02000 0.00000 0.30000 0.37000 0.12000 0.90000 100.0 Viperidae Suntravat et al. (2016) ABI 3730 https://doi.org/10.1186/s12867-016-0059-7
Bothrops_insularis 0.00000 0.15800 0.01500 0.14200 0.05100 0.00000 0.03500 0.00000 0.05400 0.43200 0.11200 0.99900 150.0 Viperidae Valente et al. (2009) ABI 3100 https://doi.org/10.1016/j.jprot.2009.01.001
Bothrops_jararaca 0.00000 0.04250 0.00400 0.13650 0.02900 0.00200 0.04050 0.00000 0.03600 0.43600 0.26100 0.98750 160.0 Viperidae Goncalves-Machado et al. (2016) Roche GS FLX https://doi.org/10.1016/j.jprot.2015.04.029
Bothrops_jararacussu 0.00000 0.01695 0.03390 0.08475 0.00000 0.00000 0.01695 0.00000 0.56290 0.25820 0.01695 0.99060 70.0 Viperidae Kashima et al. (2004) ABI 3100 https://doi.org/10.1016/j.biochi.2004.02.002
Bothrops_moojeni 0.00060 0.00950 0.00570 0.02670 0.00000 0.00440 0.01270 0.00000 0.09020 0.46630 0.26180 0.87790 130.0 Viperidae Amorim et al. (2017) Illumina https://doi.org/10.1016/j.toxicon.2017.10.025
Bothrops_neuwiedi 0.00000 0.17400 0.00700 0.08600 0.00360 0.00000 0.02900 0.00000 0.26400 0.36000 0.05300 0.97660 65.0 Viperidae Rodigues et al. (2012) MegaBACE (GE) https://doi.org/10.1016/j.jprot.2012.03.028
Bungarus_flaviceps 0.39290 0.00000 0.00440 0.00440 0.00000 0.24060 0.00000 0.35320 0.00000 0.00000 0.00000 0.99550 185.0 Elapidae Siang et al. (2010) ABI PRISM https://dx.doi.org/10.1186/1471-2199-11-24
Bungarus_multicinctus 0.64700 0.00000 0.00700 0.03800 0.00200 0.10600 0.00100 0.17400 0.00000 0.00000 0.00100 0.97600 125.0 Elapidae Jiang et al. (2011) ABI 3100 https://doi.org/10.1186/1471-2164-12-1
Cerrophidion_godmani 0.00600 0.08300 0.02600 0.00200 0.00600 0.00000 0.05500 0.00000 0.07900 0.42200 0.28100 0.96000 55.0 Viperidae Durban et al. (2011) Roche GS FLX https://doi.org/10.1186/1471-2164-12-259
Crotalus_adamanteus 0.00000 0.00726 0.01317 0.22155 0.00023 0.00054 0.05309 0.00000 0.07810 0.24410 0.19979 0.81783 244.0 Viperidae Rokyta et al. (2012) Illumina https://doi.org/10.1186/1471-2164-13-312
Crotalus_durissus 0.00000 0.00000 0.00000 0.00000 0.04000 0.00000 0.00000 0.00000 0.88000 0.03000 0.02500 0.97500 150.0 Viperidae Boldrini-Franca et al. (2009) MegaBACE (GE) https://doi.org/10.1016/j.biochi.2009.02.001
Crotalus_horridus 0.00000 0.12519 0.00775 0.00213 0.02285 0.00119 0.01106 0.00000 0.22772 0.00111 0.58229 0.98129 121.0 Viperidae Rokyta et al. (2013) Illumina https://doi.org/10.1186/1471-2164-14-394
Crotalus_scutulatus 0.00000 0.20121 0.00481 0.16943 0.02000 0.00021 0.00925 0.00000 0.18171 0.14363 0.11692 0.84717 137.3 Viperidae Strickland et al. (2018) Illumina https://dx.doi.org/10.3390/toxins10040135
Crotalus_simus 0.12700 0.13400 0.00000 0.02400 0.00000 0.00150 0.03800 0.00000 0.12100 0.06500 0.28100 0.79150 150.0 Viperidae Durban et al. (2011) Roche GS FLX https://doi.org/10.1186/1471-2164-12-259
Dispholidus_typus 0.06200 0.00000 0.05400 0.00400 0.00100 0.00150 0.00000 0.07600 0.00000 0.74640 0.05500 0.99990 135.0 Colubridae Pla et al. (2017) Illumina https://doi.org/10.1016/j.bbagen.2017.01.020
Drysdalia_coronoides 0.86000 0.00000 0.02800 0.00000 0.01400 0.00000 0.00000 0.00240 0.00000 0.00000 0.09200 0.99640 40.0 Elapidae Chatrath et al. (2011) ABI PRISM https://dx.doi.org/10.1021/pr1008916
Echis_carinatus 0.00000 0.00000 0.01934 0.17602 0.00580 0.00193 0.01161 0.00000 0.08317 0.64603 0.04836 0.99226 60.0 Viperidae Casewell et al. (2009) ABI 3730 https://doi.org/10.1186/1471-2164-10-564
Echis_coloratus 0.00000 0.00000 0.05426 0.10543 0.00620 0.00000 0.04031 0.00000 0.05426 0.68527 0.04496 0.99069 62.7 Viperidae Casewell et al. (2009) ABI 3730 https://doi.org/10.1186/1471-2164-10-564
Echis_ocellatus 0.00000 0.00000 0.00297 0.10979 0.02671 0.00193 0.01187 0.00000 0.08605 0.73591 0.02374 0.99897 40.0 Viperidae Casewell et al. (2009) ABI 3730 https://doi.org/10.1186/1471-2164-10-564
Echis_pyramidum 0.00000 0.00000 0.00000 0.24011 0.00582 0.00000 0.02639 0.00000 0.20580 0.50000 0.01979 0.99791 65.0 Viperidae Casewell et al. (2009) ABI 3730 https://doi.org/10.1186/1471-2164-10-564
Gloydius_intermedius 0.00000 0.04900 0.04900 0.04900 0.14600 0.00000 0.02400 0.00000 0.12200 0.17000 0.22000 0.82900 53.3 Viperidae Yang et al. (2015) MegaBACE (GE) https://doi.org/10.1016/j.toxicon.2015.08.010
Hypsiglena_torquata 0.00070 0.07359 0.16684 0.06176 0.00012 0.00056 0.00000 0.00000 0.00000 0.68666 0.00000 0.99023 48.0 Colubridae McGivern et al. (2014) Illumina https://doi.org/10.1186/1471-2164-15-1061
Lachesis_muta 0.01600 0.73200 0.00000 0.05800 0.00300 0.00000 0.03700 0.00000 0.04700 0.05900 0.03500 0.98700 225.0 Viperidae Lnacio et al. (2007) ABI 3100 https://doi.org/10.1534/genetics.106.056515
Micrurus_altirostris 0.86400 0.00000 0.00000 0.01300 0.00000 0.03800 0.00000 0.04900 0.00000 0.00600 0.00000 0.97000 64.0 Elapidae Correa-Netto et al. (2011) ABI 3100 https://doi.org/10.1016/j.jprot.2011.04.003
Micrurus_corallinus 0.20300 0.29530 0.00000 0.06700 0.01830 0.00150 0.00840 0.37680 0.00000 0.02100 0.00100 0.99230 72.5 Elapidae Aird et al. (2017) Illumina https://dx.doi.org/10.3390/toxins9060187
Micrurus_fulvius 0.25160 0.00000 0.00000 0.00799 0.01600 0.02154 0.00880 0.64947 0.00000 0.02947 0.00000 0.98487 55.0 Elapidae Margres et al. (2013) Illumina https://doi.org/10.1186/1471-2164-14-531
Micrurus_lemniscatus 0.67220 0.00200 0.00240 0.00480 0.01420 0.07150 0.00730 0.17200 0.00000 0.01120 0.00100 0.95860 120.0 Elapidae Aird et al. (2017) Illumina https://dx.doi.org/10.3390/toxins9060187
Micrurus_spixii 0.62640 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.36390 0.00000 0.00000 0.00000 0.99030 75.0 Elapidae Aird et al. (2017) Illumina https://dx.doi.org/10.3390/toxins9060187
Micrurus_surinamensis 0.60240 0.02670 0.00000 0.00340 0.07880 0.01070 0.00000 0.26230 0.00000 0.00140 0.00000 0.98570 68.0 Elapidae Aird et al. (2017) Illumina https://dx.doi.org/10.3390/toxins9060187
Naja_atra 0.95800 0.00000 0.00700 0.00600 0.00100 0.00100 0.00100 0.01200 0.00000 0.00400 0.00000 0.99000 150.0 Elapidae Jiang et al. (2011) ABI 3100 https://doi.org/10.1186/1471-2164-12-1
Naja_kaouthia 0.82000 0.00000 0.00000 0.00500 0.00300 0.00100 0.00200 0.13600 0.00000 0.01100 0.00000 0.97800 150.0 Elapidae Xu et al. (2017) Illumina https://doi.org/10.1016/j.jprot.2017.02.018
Opheodrys_aestivus 0.13300 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.03600 0.00000 0.00000 0.00000 0.16900 91.5 Colubridae Hargreaves et al. (2014) Illumina https://doi.org/10.1016/j.toxicon.2014.10.004
Ophiophagus_hannah 0.43000 0.00000 0.08700 0.00000 0.00800 0.01000 0.05700 0.04000 0.00000 0.24400 0.00700 0.88300 425.0 Elapidae Tan et al. (2015) Illumina https://dx.doi.org/10.1186/s12864-015-1828-2
Ovophis_okinavensis 0.00000 0.00000 0.00000 0.00579 0.00008 0.00000 0.00623 0.00000 0.00647 0.04220 0.86176 0.92253 55.0 Viperidae Aird et al. (2013) Illumina https://doi.org/10.1186/1471-2164-14-790
Pantherophis_guttatus 0.01000 0.00000 0.02200 0.00000 0.00000 0.00000 0.00000 0.12900 0.00000 0.00000 0.00000 0.16100 150.0 Colubridae Hargreaves et al. (2014) Illumina https://doi.org/10.1016/j.toxicon.2014.10.004
Phalotris_mertensi 0.01800 0.00140 0.00330 0.21600 0.00040 0.43900 0.02300 0.00000 0.00000 0.22040 0.00220 0.92370 100.0 Colubridae Campos et al. (2016) ABI 3100 https://doi.org/10.1093/gbe/evw149
Philodryas_olfersii 0.00000 0.22400 0.10400 0.05800 0.00000 0.00000 0.00000 0.00000 0.00000 0.52400 0.08900 0.99900 125.0 Colubridae Ching et al. (2006) ABI 3100 https://doi.org/10.1016/j.febslet.2006.07.010
Protobothrops_elegans 0.00000 0.01400 0.00001 0.00200 0.00600 0.00000 0.00500 0.00000 0.77100 0.08100 0.10400 0.98301 120.0 Viperidae Aird et al. (2015) Illumina https://doi.org/10.1186/s12864-015-1832-6
Protobothrops_flavoviridis 0.00000 0.00000 0.03854 0.15356 0.00960 0.00000 0.09071 0.00000 0.32136 0.26975 0.11054 0.99406 137.0 Viperidae Aird et al. (2013) Illumina https://doi.org/10.1186/1471-2164-14-790
Protobothrops_mucrosquamatus 0.00000 0.00000 0.01018 0.01018 0.00647 0.00000 0.01200 0.00000 0.25746 0.39966 0.28076 0.97671 116.0 Viperidae Aird et al. (2017) Illumina https://doi.org/10.1093/gbe/evx199
Pseudonaja_textilis 0.10633 0.00000 0.03846 0.00000 0.00452 0.01810 0.00000 0.18778 0.00000 0.00000 0.42534 0.78053 150.0 Elapidae Reeks et al. (2016) Roche GS FLX https://doi.org/10.1016/j.jprot.2015.11.019
Sistrurus_catenatus 0.00830 0.00280 0.06670 0.01390 0.06670 0.00000 0.03610 0.00000 0.28060 0.12220 0.37500 0.97230 73.5 Viperidae Pahari et al. (2007) ABI 3100 https://doi.org/10.1186/1471-2199-8-115
Thamnodynastes_strigatus 0.00200 0.00000 0.01900 0.05400 0.00000 0.00000 0.00000 0.00000 0.00000 0.84400 0.00700 0.92600 70.0 Colubridae Ching et al. (2011) ABI 3100 https://dx.doi.org/10.1021/pr200876c


Scaled version of this dataset was used in all subsequent analysis

#normalize data
scaled.dat <- dat
scaled.dat[,2:12] <- t(scale(t(scaled.dat[,2:12])))

Removing low abundance toxins.

all_dat<-read.csv("./data/scaleddata/All_25_comps.csv", header = T)

Supplementary Figure 1: Average compositon of each toxin in the venom of 52 species in our study

Supplementary Figure 2: Tiles plot showing presence/absensce of each toxin in 52 snake species

Coloured tiles indicate presence of toxin, while white tiles indicate absence. As oberserved from the plot we included components that are found in at least 50% of snakes (Red dashed line). This ensured low retention of structural zeros in our complete dataset.

DTT Analysis

Higher between clade disparity in phenotypes (disparity between large diversifying lineages) is an indicator of ‘early burst’ evolution (Simpson 1953.; Tran Lucy A. P. 2014). Disparity through time measures instances when disparity was highest within versus among subclades (Harmon et al. 2003; Tran Lucy A. P. 2014). It also provides us an estimate of whether phenotypic diversity occurred early or late in the diversification process. We performed a disparity through time (DTT) analysis to get estimates of when overall disparity in toxin levels arose and the subsequent trends they experienced. Since the DTT analysis estimates overall disparity throughout a phylogeny, it could be susceptible to sampling bias, leading to interpretations that are not biologically meaningful. For this reason, we restricted the DTT analysis to the supplementary. Readers interested in the within/between clade disparity estimates can refer to the online supplementary material, however, we advise caution during its interpretation.

Relative disparity through time (DTT) for each of the 10 toxins in our dataset. DTT gives an estimate of within versus between clade diversity, and when phenotypic diversity arose. Abundant toxins families (with names highlighted in red) such as SVMP, TFTx, and ePLA2 show high disparity early in the history of snakes. CRISP showed high disparity early on followed by a decline, and a sudden spike in disparity in modern lineages. This implies that diversity in expression levels for these toxin families was higher in the past, with a reduction in modern lineages. Most other less abundant toxins display an oscillatory nature in disparity. The empirical curve represented by solid black line was compared to the median curve of a null simulation (Brownian model), represented by the fine dashed line. Points where the empirical curve deviated significantly from the null model are represented by red dots accompanied by their respective p-values estimated using the rank envelope test (Murrell 2018). 95% confidence intervals of simulated disparity is represented by thick dashed line. The DTT plots show that the prediction of early disparity in phenotypes is not a general trend in gene expression of most snake venom toxins.

Load data

tree<-read.newick("./Final_tree.tre")
trait<-read.csv("./data/scaleddata/scaled_data_fot_dtt.csv",row.names = 1, header = T)

BPP

The following script was repeated for all 11 toxins

trait<-as.matrix(trait)[,'BPP.CNP']
venomdtt <-dtt(tree, trait, index = 'avg.sq', nsim = 5000,
               CI = 0.95, plot = F)

Rank envelope function taken from (Murrell 2018)

rank_env_dtt(venomdtt)

## Rank envelope test
##  p-value of the test: 0.00389922 (ties method: midrank)
##  p-interval         : (0.00019996, 0.00759848)



CRISP

## Rank envelope test
##  p-value of the test: 0.00389922 (ties method: midrank)
##  p-interval         : (0.00019996, 0.00759848)



CTL

## Rank envelope test
##  p-value of the test: 0.01939612 (ties method: midrank)
##  p-interval         : (0.01679664, 0.0219956)

GF

## Rank envelope test
##  p-value of the test: 0.00359928 (ties method: midrank)
##  p-interval         : (0.00019996, 0.0069986)



KSPI

## Rank envelope test
##  p-value of the test: 0.00339932 (ties method: midrank)
##  p-interval         : (0.00019996, 0.00659868)



LAAO

## Rank envelope test
##  p-value of the test: 0.00359928 (ties method: midrank)
##  p-interval         : (0.00019996, 0.0069986)



SVMP

## Rank envelope test
##  p-value of the test: 0.4585083 (ties method: midrank)
##  p-interval         : (0.4569086, 0.460108)



SVSP

## Rank envelope test
##  p-value of the test: 0.00429914 (ties method: midrank)
##  p-interval         : (0.00019996, 0.00839832)



TFTx

## Rank envelope test
##  p-value of the test: 0.00479904 (ties method: midrank)
##  p-interval         : (0.00019996, 0.00939812)



vPLA2

## Rank envelope test
##  p-value of the test: 0.01159768 (ties method: midrank)
##  p-interval         : (0.0079984, 0.01519696)



ePLA2

## Rank envelope test
##  p-value of the test: 0.1921616 (ties method: midrank)
##  p-interval         : (0.1903619, 0.1939612)

Bayesian analysis of macroevolutionary mixtures (BAMM)

We estimated evolutionary rate dynamics under the Bayesian analysis of macroevolutionary mixtures (BAMM), to observe shifts in rates of gene expression evolution that took place at different points throughout the history of snake venom evolution. Parameters for phenotypic evolution suggested in the online BAMM documentation (http://bamm-project.org/) were used to run the model. The model was run for 108 generations with MCMC write frequency of 104. Priors were obtained using the setBAMMpriors function in BAMMTools (Rabosky et al. 2014).

Supplementary Figure 3: Convergence diagnosis of MCMC chiains

Convergence of MCMC chains were determined by visual inspection, and plotting effective sample size of log-likelihood and number of shifts in each sample, both of which well exceeded the recommended value of 200

BPP

mcmcout<- read.csv("./Analysis_scaled/ControlFiles_scaleddata/BPP/mcmc_out.txt", header = T)
plot(mcmcout$logLik ~ mcmcout$generation, type="l")

burnstart <- floor(0.8 * nrow(mcmcout))
postburn <- mcmcout[burnstart:nrow(mcmcout),]

effectiveSize(postburn$N_shifts)
##     var1 
## 1873.928
effectiveSize(postburn$logLik)
##     var1 
## 1743.364
plot(postburn$logLik ~ postburn$generation, type="l")

plotPrior(mcmcout, expectedNumberOfShifts = 1)



CRISP

## var1 
## 8001
##     var1 
## 6696.861



CTL

##     var1 
## 7961.472
## var1 
## 9001



GF

## var1 
## 9001
##     var1 
## 8567.747



KUN

##     var1 
## 6775.222
##     var1 
## 5244.185



LAAO

## var1 
## 3001
## var1 
## 3001



ePLA2

##    var1 
## 2353.79
##     var1 
## 1715.655



vPLA2

##     var1 
## 1750.436
##     var1 
## 1777.745



SVMP

##     var1 
## 8014.687
##     var1 
## 7752.112



SVSP

## var1 
## 8501
##     var1 
## 7195.953



TFTx

##     var1 
## 3196.441
##    var1 
## 739.902

Computing Bayes factor

Bayes factor estimates support the occurrence of at least one rate shift in all toxin families, indicating that toxins families have experienced changes in their evolutionary rates since becoming a part of the venom arsenal. A Bayes factor of 20 and above means that we have enough certainty that the number of rate shifts above the lowest number of shifts in likely.

Supplementary Figure 4: Bayes factor estimates for shift configurations

BPP

mcmcout<- read.csv("./Analysis_scaled/ControlFiles_scaleddata/BPP/mcmc_out.txt", header = T)

bfmat<-computeBayesFactors(mcmcout,expectedNumberOfShifts = 1, burnin = 0.8)

# Bayes factor matrix
bfmat
##            3         4         5         6         7         8         9
## 3   1.000000 0.3031136 0.1872172 0.1277006 0.1004248 0.1136676 0.1261433
## 4   3.299094 1.0000000 0.6176471 0.4212963 0.3313107 0.3750000 0.4161585
## 5   5.341390 1.6190476 1.0000000 0.6820988 0.5364078 0.6071429 0.6737805
## 6   7.830816 2.3736264 1.4660633 1.0000000 0.7864078 0.8901099 0.9878049
## 7   9.957704 3.0183150 1.8642534 1.2716049 1.0000000 1.1318681 1.2560976
## 8   8.797583 2.6666667 1.6470588 1.1234568 0.8834951 1.0000000 1.1097561
## 9   7.927492 2.4029304 1.4841629 1.0123457 0.7961165 0.9010989 1.0000000
## 10  4.640483 1.4065934 0.8687783 0.5925926 0.4660194 0.5274725 0.5853659
## 11  3.093656 0.9377289 0.5791855 0.3950617 0.3106796 0.3516484 0.3902439
## 12  3.093656 0.9377289 0.5791855 0.3950617 0.3106796 0.3516484 0.3902439
## 14  6.187311 1.8754579 1.1583710 0.7901235 0.6213592 0.7032967 0.7804878
## 15 12.374622 3.7509158 2.3167421 1.5802469 1.2427184 1.4065934 1.5609756
##           10        11        12        14         15
## 3  0.2154948 0.3232422 0.3232422 0.1616211 0.08081055
## 4  0.7109375 1.0664062 1.0664063 0.5332031 0.26660156
## 5  1.1510417 1.7265625 1.7265625 0.8632813 0.43164063
## 6  1.6875000 2.5312500 2.5312500 1.2656250 0.63281250
## 7  2.1458333 3.2187500 3.2187500 1.6093750 0.80468750
## 8  1.8958333 2.8437500 2.8437500 1.4218750 0.71093750
## 9  1.7083333 2.5625000 2.5625000 1.2812500 0.64062500
## 10 1.0000000 1.5000000 1.5000000 0.7500000 0.37500000
## 11 0.6666667 1.0000000 1.0000000 0.5000000 0.25000000
## 12 0.6666667 1.0000000 1.0000000 0.5000000 0.25000000
## 14 1.3333333 2.0000000 2.0000000 1.0000000 0.50000000
## 15 2.6666667 4.0000000 4.0000000 2.0000000 1.00000000

CRISP

##            0           1           2            3            4           5
## 0     1.0000 0.002663116 0.001402525 0.0009784736 0.0008116883 0.000697026
## 1   375.5000 1.000000000 0.526647966 0.3674168297 0.3047889610 0.261733271
## 2   713.0000 1.898801598 1.000000000 0.6976516634 0.5787337662 0.496979554
## 3  1022.0000 2.721704394 1.433380084 1.0000000000 0.8295454545 0.712360595
## 4  1232.0000 3.280958722 1.727910238 1.2054794521 1.0000000000 0.858736059
## 5  1434.6667 3.820683533 2.012155213 1.4037834312 1.1645021645 1.000000000
## 6  1562.6667 4.161562361 2.191678354 1.5290280496 1.2683982684 1.089219331
## 7  1664.0000 4.431424767 2.333800842 1.6281800391 1.3506493506 1.159851301
## 8  1664.0000 4.431424767 2.333800842 1.6281800391 1.3506493506 1.159851301
## 9  1877.3333 4.999556147 2.633006078 1.8369210698 1.5238095238 1.308550186
## 10 1450.6667 3.863293387 2.034595605 1.4194390085 1.1774891775 1.011152416
## 11 1536.0000 4.090545939 2.154277700 1.5029354207 1.2467532468 1.070631970
## 12 1024.0000 2.727030626 1.436185133 1.0019569472 0.8311688312 0.713754647
## 13  682.6667 1.818020417 0.957456755 0.6679712981 0.5541125541 0.475836431
## 14 1365.3333 3.636040834 1.914913511 1.3359425962 1.1082251082 0.951672862
##               6            7            8            9           10
## 0  0.0006399317 0.0006009615 0.0006009615 0.0005326705 0.0006893382
## 1  0.2402943686 0.2256610577 0.2256610577 0.2000177557 0.2588465074
## 2  0.4562713311 0.4284855769 0.4284855769 0.3797940341 0.4914981618
## 3  0.6540102389 0.6141826923 0.6141826923 0.5443892045 0.7045036765
## 4  0.7883959044 0.7403846154 0.7403846154 0.6562500000 0.8492647059
## 5  0.9180887372 0.8621794872 0.8621794872 0.7642045455 0.9889705882
## 6  1.0000000000 0.9391025641 0.9391025641 0.8323863636 1.0772058824
## 7  1.0648464164 1.0000000000 1.0000000000 0.8863636364 1.1470588235
## 8  1.0648464164 1.0000000000 1.0000000000 0.8863636364 1.1470588235
## 9  1.2013651877 1.1282051282 1.1282051282 1.0000000000 1.2941176471
## 10 0.9283276451 0.8717948718 0.8717948718 0.7727272727 1.0000000000
## 11 0.9829351536 0.9230769231 0.9230769231 0.8181818182 1.0588235294
## 12 0.6552901024 0.6153846154 0.6153846154 0.5454545455 0.7058823529
## 13 0.4368600683 0.4102564103 0.4102564103 0.3636363636 0.4705882353
## 14 0.8737201365 0.8205128205 0.8205128205 0.7272727273 0.9411764706
##              11           12          13           14
## 0  0.0006510417 0.0009765625 0.001464844 0.0007324219
## 1  0.2444661458 0.3666992188 0.550048828 0.2750244141
## 2  0.4641927083 0.6962890625 1.044433594 0.5222167969
## 3  0.6653645833 0.9980468750 1.497070313 0.7485351563
## 4  0.8020833333 1.2031250000 1.804687500 0.9023437500
## 5  0.9340277778 1.4010416667 2.101562500 1.0507812500
## 6  1.0173611111 1.5260416667 2.289062500 1.1445312500
## 7  1.0833333333 1.6250000000 2.437500000 1.2187500000
## 8  1.0833333333 1.6250000000 2.437500000 1.2187500000
## 9  1.2222222222 1.8333333333 2.750000000 1.3750000000
## 10 0.9444444444 1.4166666667 2.125000000 1.0625000000
## 11 1.0000000000 1.5000000000 2.250000000 1.1250000000
## 12 0.6666666667 1.0000000000 1.500000000 0.7500000000
## 13 0.4444444444 0.6666666667 1.000000000 0.5000000000
## 14 0.8888888889 1.3333333333 2.000000000 1.0000000000

CTL

##             1           2           3           4           5           6
## 1     1.00000  0.04133545 0.007011866 0.003128008 0.001852908 0.001314725
## 2    24.19231  1.00000000 0.169633225 0.075673725 0.044826112 0.031806230
## 3   142.61538  5.89507154 1.000000000 0.446102021 0.264253136 0.187500000
## 4   319.69231 13.21462639 2.241639698 1.000000000 0.592360319 0.420307443
## 5   539.69231 22.30842607 3.784250270 1.688161694 1.000000000 0.709546926
## 6   760.61538 31.44038156 5.333333333 2.379210780 1.409350057 1.000000000
## 7   836.92308 34.59459459 5.868392665 2.617901829 1.550741163 1.100323625
## 8   866.46154 35.81558029 6.075512406 2.710298364 1.605473204 1.139158576
## 9  1048.61538 43.34499205 7.352750809 3.280076997 1.942987457 1.378640777
## 10  836.92308 34.59459459 5.868392665 2.617901829 1.550741163 1.100323625
## 11  886.15385 36.62957075 6.213592233 2.771896054 1.641961231 1.165048544
## 12  512.00000 21.16375199 3.590075512 1.601539942 0.948688712 0.673139159
## 13  236.30769  9.76788553 1.656957929 0.739172281 0.437856328 0.310679612
## 14  472.61538 19.53577107 3.313915858 1.478344562 0.875712657 0.621359223
## 15  945.23077 39.07154213 6.627831715 2.956689124 1.751425314 1.242718447
## 16  630.15385 26.04769475 4.418554477 1.971126083 1.167616876 0.828478964
##              7           8            9          10          11          12
## 1  0.001194853 0.001154119 0.0009536385 0.001194853 0.001128472 0.001953125
## 2  0.028906250 0.027920810 0.0230707160 0.028906250 0.027300347 0.047250601
## 3  0.170404412 0.164595170 0.1360035211 0.170404412 0.160937500 0.278545673
## 4  0.381985294 0.368963068 0.3048708920 0.381985294 0.360763889 0.624399038
## 5  0.644852941 0.622869318 0.5146713615 0.644852941 0.609027778 1.054086538
## 6  0.908823529 0.877840909 0.7253521127 0.908823529 0.858333333 1.485576923
## 7  1.000000000 0.965909091 0.7981220657 1.000000000 0.944444444 1.634615385
## 8  1.035294118 1.000000000 0.8262910798 1.035294118 0.977777778 1.692307692
## 9  1.252941176 1.210227273 1.0000000000 1.252941176 1.183333333 2.048076923
## 10 1.000000000 0.965909091 0.7981220657 1.000000000 0.944444444 1.634615385
## 11 1.058823529 1.022727273 0.8450704225 1.058823529 1.000000000 1.730769231
## 12 0.611764706 0.590909091 0.4882629108 0.611764706 0.577777778 1.000000000
## 13 0.282352941 0.272727273 0.2253521127 0.282352941 0.266666667 0.461538462
## 14 0.564705882 0.545454545 0.4507042254 0.564705882 0.533333333 0.923076923
## 15 1.129411765 1.090909091 0.9014084507 1.129411765 1.066666667 1.846153846
## 16 0.752941176 0.727272727 0.6009389671 0.752941176 0.711111111 1.230769231
##             13          14          15          16
## 1  0.004231771 0.002115885 0.001057943 0.001586914
## 2  0.102376302 0.051188151 0.025594076 0.038391113
## 3  0.603515625 0.301757813 0.150878906 0.226318359
## 4  1.352864583 0.676432292 0.338216146 0.507324219
## 5  2.283854167 1.141927083 0.570963542 0.856445313
## 6  3.218750000 1.609375000 0.804687500 1.207031250
## 7  3.541666667 1.770833333 0.885416667 1.328125000
## 8  3.666666667 1.833333333 0.916666667 1.375000000
## 9  4.437500000 2.218750000 1.109375000 1.664062500
## 10 3.541666667 1.770833333 0.885416667 1.328125000
## 11 3.750000000 1.875000000 0.937500000 1.406250000
## 12 2.166666667 1.083333333 0.541666667 0.812500000
## 13 1.000000000 0.500000000 0.250000000 0.375000000
## 14 2.000000000 1.000000000 0.500000000 0.750000000
## 15 4.000000000 2.000000000 1.000000000 1.500000000
## 16 2.666666667 1.333333333 0.666666667 1.000000000

GF

##           1            2            3            4           5            6
## 1     1.000 0.0006912442 0.0003205128 0.0002042484 0.000161221 0.0001233553
## 2  1446.667 1.0000000000 0.4636752137 0.2954793028 0.233233018 0.1784539474
## 3  3120.000 2.1566820276 1.0000000000 0.6372549020 0.503009458 0.3848684211
## 4  4896.000 3.3843317972 1.5692307692 1.0000000000 0.789337919 0.6039473684
## 5  6202.667 4.2875576037 1.9880341880 1.2668845316 1.000000000 0.7651315789
## 6  8106.667 5.6036866359 2.5982905983 1.6557734205 1.306964746 1.0000000000
## 7  7850.667 5.4267281106 2.5162393162 1.6034858388 1.265692175 0.9684210526
## 8  9344.000 6.4589861751 2.9948717949 1.9084967320 1.506448839 1.1526315789
## 9  7594.667 5.2497695853 2.4341880342 1.5511982571 1.224419604 0.9368421053
## 10 6144.000 4.2470046083 1.9692307692 1.2549019608 0.990541702 0.7578947368
## 11 3072.000 2.1235023041 0.9846153846 0.6274509804 0.495270851 0.3789473684
## 12 3413.333 2.3594470046 1.0940170940 0.6971677560 0.550300946 0.4210526316
## 13 1365.333 0.9437788018 0.4376068376 0.2788671024 0.220120378 0.1684210526
## 14 2730.667 1.8875576037 0.8752136752 0.5577342048 0.440240757 0.3368421053
## 15 5461.333 3.7751152074 1.7504273504 1.1154684096 0.880481513 0.6736842105
##               7            8            9           10           11
## 1  0.0001273777 0.0001070205 0.0001316713 0.0001627604 0.0003255208
## 2  0.1842730978 0.1548230594 0.1904845506 0.2354600694 0.4709201389
## 3  0.3974184783 0.3339041096 0.4108146067 0.5078125000 1.0156250000
## 4  0.6236413043 0.5239726027 0.6446629213 0.7968750000 1.5937500000
## 5  0.7900815217 0.6638127854 0.8167134831 1.0095486111 2.0190972222
## 6  1.0326086957 0.8675799087 1.0674157303 1.3194444444 2.6388888889
## 7  1.0000000000 0.8401826484 1.0337078652 1.2777777778 2.5555555556
## 8  1.1902173913 1.0000000000 1.2303370787 1.5208333333 3.0416666667
## 9  0.9673913043 0.8127853881 1.0000000000 1.2361111111 2.4722222222
## 10 0.7826086957 0.6575342466 0.8089887640 1.0000000000 2.0000000000
## 11 0.3913043478 0.3287671233 0.4044943820 0.5000000000 1.0000000000
## 12 0.4347826087 0.3652968037 0.4494382022 0.5555555556 1.1111111111
## 13 0.1739130435 0.1461187215 0.1797752809 0.2222222222 0.4444444444
## 14 0.3478260870 0.2922374429 0.3595505618 0.4444444444 0.8888888889
## 15 0.6956521739 0.5844748858 0.7191011236 0.8888888889 1.7777777778
##              12           13           14           15
## 1  0.0002929688 0.0007324219 0.0003662109 0.0001831055
## 2  0.4238281250 1.0595703125 0.5297851563 0.2648925781
## 3  0.9140625000 2.2851562500 1.1425781250 0.5712890625
## 4  1.4343750000 3.5859375000 1.7929687500 0.8964843750
## 5  1.8171875000 4.5429687500 2.2714843750 1.1357421875
## 6  2.3750000000 5.9375000000 2.9687500000 1.4843750000
## 7  2.3000000000 5.7500000000 2.8750000000 1.4375000000
## 8  2.7375000000 6.8437500000 3.4218750000 1.7109375000
## 9  2.2250000000 5.5625000000 2.7812500000 1.3906250000
## 10 1.8000000000 4.5000000000 2.2500000000 1.1250000000
## 11 0.9000000000 2.2500000000 1.1250000000 0.5625000000
## 12 1.0000000000 2.5000000000 1.2500000000 0.6250000000
## 13 0.4000000000 1.0000000000 0.5000000000 0.2500000000
## 14 0.8000000000 2.0000000000 1.0000000000 0.5000000000
## 15 1.6000000000 4.0000000000 2.0000000000 1.0000000000

KUN

##              1           2            3            4          5            6
## 1     1.000000   0.1021127  0.005484115  0.001942658 0.00108468 0.0007602768
## 2     9.793103   1.0000000  0.053706505  0.019024652 0.01062238 0.0074454698
## 3   182.344828  18.6197183  1.000000000  0.354233655 0.19778576 0.1386325503
## 4   514.758621  52.5633803  2.822995461  1.000000000 0.55834829 0.3913590604
## 5   921.931034  94.1408451  5.055975794  1.790996785 1.00000000 0.7009228188
## 6  1315.310345 134.3098592  7.213313162  2.555198285 1.42669060 1.0000000000
## 7  1600.000000 163.3802817  8.774583964  3.108252947 1.73548773 1.2164429530
## 8  2193.655172 224.0000000 12.030257186  4.261521972 2.37941352 1.6677852349
## 9  2162.758621 220.8450704 11.860816944  4.201500536 2.34590066 1.6442953020
## 10 2560.000000 261.4084507 14.039334342  4.973204716 2.77678037 1.9463087248
## 11 2754.206897 281.2394366 15.104387292  5.350482315 2.98743268 2.0939597315
## 12 3177.931034 324.5070423 17.428139183  6.173633441 3.44703770 2.4161073826
## 13 3248.551724 331.7183099 17.815431165  6.310825295 3.52363854 2.4697986577
## 14 3107.310345 317.2957746 17.040847201  6.036441586 3.37043686 2.3624161074
## 15 2259.862069 230.7605634 12.393343419  4.390139335 2.45122681 1.7181208054
## 16 2259.862069 230.7605634 12.393343419  4.390139335 2.45122681 1.7181208054
## 17 6779.586207 692.2816901 37.180030257 13.170418006 7.35368043 5.1543624161
## 18 4519.724138 461.5211268 24.786686838  8.780278671 4.90245362 3.4362416107
##             7            8            9          10           11           12
## 1  0.00062500 0.0004558602 0.0004623724 0.000390625 0.0003630809 0.0003146701
## 2  0.00612069 0.0044642857 0.0045280612 0.003825431 0.0035556891 0.0030815972
## 3  0.11396552 0.0831237425 0.0843112245 0.071228448 0.0662059295 0.0573784722
## 4  0.32172414 0.2346579477 0.2380102041 0.201077586 0.1868990385 0.1619791667
## 5  0.57620690 0.4202716298 0.4262755102 0.360129310 0.3347355769 0.2901041667
## 6  0.82206897 0.5995975855 0.6081632653 0.513793103 0.4775641026 0.4138888889
## 7  1.00000000 0.7293762575 0.7397959184 0.625000000 0.5809294872 0.5034722222
## 8  1.37103448 1.0000000000 1.0142857143 0.856896552 0.7964743590 0.6902777778
## 9  1.35172414 0.9859154930 1.0000000000 0.844827586 0.7852564103 0.6805555556
## 10 1.60000000 1.1670020121 1.1836734694 1.000000000 0.9294871795 0.8055555556
## 11 1.72137931 1.2555331992 1.2734693878 1.075862069 1.0000000000 0.8666666667
## 12 1.98620690 1.4486921529 1.4693877551 1.241379310 1.1538461538 1.0000000000
## 13 2.03034483 1.4808853119 1.5020408163 1.268965517 1.1794871795 1.0222222222
## 14 1.94206897 1.4164989940 1.4367346939 1.213793103 1.1282051282 0.9777777778
## 15 1.41241379 1.0301810865 1.0448979592 0.882758621 0.8205128205 0.7111111111
## 16 1.41241379 1.0301810865 1.0448979592 0.882758621 0.8205128205 0.7111111111
## 17 4.23724138 3.0905432596 3.1346938776 2.648275862 2.4615384615 2.1333333333
## 18 2.82482759 2.0603621730 2.0897959184 1.765517241 1.6410256410 1.4222222222
##              13           14           15           16           17
## 1  0.0003078295 0.0003218217 0.0004425049 0.0004425049 0.0001475016
## 2  0.0030146060 0.0031516335 0.0043334961 0.0043334961 0.0014444987
## 3  0.0561311141 0.0586825284 0.0806884766 0.0806884766 0.0268961589
## 4  0.1584578804 0.1656605114 0.2277832031 0.2277832031 0.0759277344
## 5  0.2837975543 0.2966974432 0.4079589844 0.4079589844 0.1359863281
## 6  0.4048913043 0.4232954545 0.5820312500 0.5820312500 0.1940104167
## 7  0.4925271739 0.5149147727 0.7080078125 0.7080078125 0.2360026042
## 8  0.6752717391 0.7059659091 0.9707031250 0.9707031250 0.3235677083
## 9  0.6657608696 0.6960227273 0.9570312500 0.9570312500 0.3190104167
## 10 0.7880434783 0.8238636364 1.1328125000 1.1328125000 0.3776041667
## 11 0.8478260870 0.8863636364 1.2187500000 1.2187500000 0.4062500000
## 12 0.9782608696 1.0227272727 1.4062500000 1.4062500000 0.4687500000
## 13 1.0000000000 1.0454545455 1.4375000000 1.4375000000 0.4791666667
## 14 0.9565217391 1.0000000000 1.3750000000 1.3750000000 0.4583333333
## 15 0.6956521739 0.7272727273 1.0000000000 1.0000000000 0.3333333333
## 16 0.6956521739 0.7272727273 1.0000000000 1.0000000000 0.3333333333
## 17 2.0869565217 2.1818181818 3.0000000000 3.0000000000 1.0000000000
## 18 1.3913043478 1.4545454545 2.0000000000 2.0000000000 0.6666666667
##              18
## 1  0.0002212524
## 2  0.0021667480
## 3  0.0403442383
## 4  0.1138916016
## 5  0.2039794922
## 6  0.2910156250
## 7  0.3540039063
## 8  0.4853515625
## 9  0.4785156250
## 10 0.5664062500
## 11 0.6093750000
## 12 0.7031250000
## 13 0.7187500000
## 14 0.6875000000
## 15 0.5000000000
## 16 0.5000000000
## 17 1.5000000000
## 18 1.0000000000

LAAO

##              0           1           2            3           4           5
## 0     1.000000   0.1065089  0.01294964  0.006276151 0.004136029 0.002826633
## 1     9.388889   1.0000000  0.12158273  0.058926081 0.038832721 0.026538945
## 2    77.222222   8.2248521  1.00000000  0.484658298 0.319393382 0.218278894
## 3   159.333333  16.9704142  2.06330935  1.000000000 0.659007353 0.450376884
## 4   241.777778  25.7514793  3.13093525  1.517433752 1.000000000 0.683417085
## 5   353.777778  37.6804734  4.58129496  2.220362622 1.463235294 1.000000000
## 6   392.888889  41.8461538  5.08776978  2.465829847 1.625000000 1.110552764
## 7   394.666667  42.0355030  5.11079137  2.476987448 1.632352941 1.115577889
## 8   440.888889  46.9585799  5.70935252  2.767085077 1.823529412 1.246231156
## 9   369.777778  39.3846154  4.78848921  2.320781032 1.529411765 1.045226131
## 10  284.444444  30.2958580  3.68345324  1.785216179 1.176470588 0.804020101
## 11  341.333333  36.3550296  4.42014388  2.142259414 1.411764706 0.964824121
## 12  341.333333  36.3550296  4.42014388  2.142259414 1.411764706 0.964824121
## 13  227.555556  24.2366864  2.94676259  1.428172943 0.941176471 0.643216080
## 15 1820.444444 193.8934911 23.57410072 11.425383543 7.529411765 5.145728643
##              6           7           8           9          10          11
## 0  0.002545249 0.002533784 0.002268145 0.002704327 0.003515625 0.002929687
## 1  0.023897059 0.023789414 0.021295363 0.025390625 0.033007812 0.027506510
## 2  0.196549774 0.195664414 0.175151210 0.208834135 0.271484375 0.226236979
## 3  0.405542986 0.403716216 0.361391129 0.430889423 0.560156250 0.466796875
## 4  0.615384615 0.612612613 0.548387097 0.653846154 0.850000000 0.708333333
## 5  0.900452489 0.896396396 0.802419355 0.956730769 1.243750000 1.036458333
## 6  1.000000000 0.995495495 0.891129032 1.062500000 1.381250000 1.151041667
## 7  1.004524887 1.000000000 0.895161290 1.067307692 1.387500000 1.156250000
## 8  1.122171946 1.117117117 1.000000000 1.192307692 1.550000000 1.291666667
## 9  0.941176471 0.936936937 0.838709677 1.000000000 1.300000000 1.083333333
## 10 0.723981900 0.720720721 0.645161290 0.769230769 1.000000000 0.833333333
## 11 0.868778281 0.864864865 0.774193548 0.923076923 1.200000000 1.000000000
## 12 0.868778281 0.864864865 0.774193548 0.923076923 1.200000000 1.000000000
## 13 0.579185520 0.576576577 0.516129032 0.615384615 0.800000000 0.666666667
## 15 4.633484163 4.612612613 4.129032258 4.923076923 6.400000000 5.333333333
##             12          13           15
## 0  0.002929688 0.004394531 0.0005493164
## 1  0.027506510 0.041259766 0.0051574707
## 2  0.226236979 0.339355469 0.0424194336
## 3  0.466796875 0.700195313 0.0875244141
## 4  0.708333333 1.062500000 0.1328125000
## 5  1.036458333 1.554687500 0.1943359375
## 6  1.151041667 1.726562500 0.2158203125
## 7  1.156250000 1.734375000 0.2167968750
## 8  1.291666667 1.937500000 0.2421875000
## 9  1.083333333 1.625000000 0.2031250000
## 10 0.833333333 1.250000000 0.1562500000
## 11 1.000000000 1.500000000 0.1875000000
## 12 1.000000000 1.500000000 0.1875000000
## 13 0.666666667 1.000000000 0.1250000000
## 15 5.333333333 8.000000000 1.0000000000

vPLA2

##            1          2         3          4          5          6          7
## 1   1.000000 0.03659081 0.0196737 0.01790393 0.01856884 0.02628205 0.03416667
## 2  27.329268 1.00000000 0.5376679 0.48930131 0.50747283 0.71826923 0.93375000
## 3  50.829268 1.85988398 1.0000000 0.91004367 0.94384058 1.33589744 1.73666667
## 4  55.853659 2.04373048 1.0988484 1.00000000 1.03713768 1.46794872 1.90833333
## 5  53.853659 1.97054886 1.0595010 0.96419214 1.00000000 1.41538462 1.84000000
## 6  38.048780 1.39223561 0.7485605 0.68122271 0.70652174 1.00000000 1.30000000
## 7  29.268293 1.07095047 0.5758157 0.52401747 0.54347826 0.76923077 1.00000000
## 8  24.195122 0.88531905 0.4760077 0.43318777 0.44927536 0.63589744 0.82666667
## 9  14.048780 0.51405622 0.2763916 0.25152838 0.26086957 0.36923077 0.48000000
## 10  9.365854 0.34270415 0.1842610 0.16768559 0.17391304 0.24615385 0.32000000
## 11  6.243902 0.22846943 0.1228407 0.11179039 0.11594203 0.16410256 0.21333333
## 12 12.487805 0.45693887 0.2456814 0.22358079 0.23188406 0.32820513 0.42666667
##             8          9        10        11         12
## 1  0.04133065 0.07118056 0.1067708 0.1601562 0.08007813
## 2  1.12953629 1.94531250 2.9179688 4.3769531 2.18847656
## 3  2.10080645 3.61805556 5.4270833 8.1406250 4.07031250
## 4  2.30846774 3.97569444 5.9635417 8.9453125 4.47265625
## 5  2.22580645 3.83333333 5.7500000 8.6250000 4.31250000
## 6  1.57258065 2.70833333 4.0625000 6.0937500 3.04687500
## 7  1.20967742 2.08333333 3.1250000 4.6875000 2.34375000
## 8  1.00000000 1.72222222 2.5833333 3.8750000 1.93750000
## 9  0.58064516 1.00000000 1.5000000 2.2500000 1.12500000
## 10 0.38709677 0.66666667 1.0000000 1.5000000 0.75000000
## 11 0.25806452 0.44444444 0.6666667 1.0000000 0.50000000
## 12 0.51612903 0.88888889 1.3333333 2.0000000 1.00000000

ePLA2

##             1            2            3            4            5            6
## 1      1.0000   0.02857143  0.001447876 0.0005215577 0.0002741228 0.0001897773
## 2     35.0000   1.00000000  0.050675676 0.0182545202 0.0095942982 0.0066422065
## 3    690.6667  19.73333333  1.000000000 0.3602225313 0.1893274854 0.1310728745
## 4   1917.3333  54.78095238  2.776061776 1.0000000000 0.5255847953 0.3638663968
## 5   3648.0000 104.22857143  5.281853282 1.9026425591 1.0000000000 0.6923076923
## 6   5269.3333 150.55238095  7.629343629 2.7482614743 1.4444444444 1.0000000000
## 7   7050.6667 201.44761905 10.208494208 3.6773296245 1.9327485380 1.3380566802
## 8   8746.6667 249.90476190 12.664092664 4.5618915160 2.3976608187 1.6599190283
## 9  10453.3333 298.66666667 15.135135135 5.4520166898 2.8654970760 1.9838056680
## 10 10069.3333 287.69523810 14.579150579 5.2517385257 2.7602339181 1.9109311741
## 11 10581.3333 302.32380952 15.320463320 5.5187760779 2.9005847953 2.0080971660
## 12 10922.6667 312.07619048 15.814671815 5.6968011127 2.9941520468 2.0728744939
## 13 13653.3333 390.09523810 19.768339768 7.1210013908 3.7426900585 2.5910931174
## 14  9557.3333 273.06666667 13.837837838 4.9847009736 2.6198830409 1.8137651822
## 15  8192.0000 234.05714286 11.861003861 4.2726008345 2.2456140351 1.5546558704
## 16  5461.3333 156.03809524  7.907335907 2.8484005563 1.4970760234 1.0364372470
## 17 10922.6667 312.07619048 15.814671815 5.6968011127 2.9941520468 2.0728744939
##               7            8            9           10           11
## 1  0.0001418306 0.0001143293 9.566327e-05 9.931144e-05 9.450605e-05
## 2  0.0049640696 0.0040015244 3.348214e-03 3.475900e-03 3.307712e-03
## 3  0.0979576399 0.0789634146 6.607143e-02 6.859110e-02 6.527218e-02
## 4  0.2719364599 0.2192073171 1.834184e-01 1.904131e-01 1.811996e-01
## 5  0.5173978820 0.4170731707 3.489796e-01 3.622881e-01 3.447581e-01
## 6  0.7473524962 0.6024390244 5.040816e-01 5.233051e-01 4.979839e-01
## 7  1.0000000000 0.8060975610 6.744898e-01 7.002119e-01 6.663306e-01
## 8  1.2405446293 1.0000000000 8.367347e-01 8.686441e-01 8.266129e-01
## 9  1.4826021180 1.1951219512 1.000000e+00 1.038136e+00 9.879032e-01
## 10 1.4281391831 1.1512195122 9.632653e-01 1.000000e+00 9.516129e-01
## 11 1.5007564297 1.2097560976 1.012245e+00 1.050847e+00 1.000000e+00
## 12 1.5491679274 1.2487804878 1.044898e+00 1.084746e+00 1.032258e+00
## 13 1.9364599092 1.5609756098 1.306122e+00 1.355932e+00 1.290323e+00
## 14 1.3555219365 1.0926829268 9.142857e-01 9.491525e-01 9.032258e-01
## 15 1.1618759455 0.9365853659 7.836735e-01 8.135593e-01 7.741935e-01
## 16 0.7745839637 0.6243902439 5.224490e-01 5.423729e-01 5.161290e-01
## 17 1.5491679274 1.2487804878 1.044898e+00 1.084746e+00 1.032258e+00
##              12           13           14           15           16
## 1  9.155273e-05 7.324219e-05 0.0001046317 0.0001220703 0.0001831055
## 2  3.204346e-03 2.563477e-03 0.0036621094 0.0042724609 0.0064086914
## 3  6.323242e-02 5.058594e-02 0.0722656250 0.0843098958 0.1264648438
## 4  1.755371e-01 1.404297e-01 0.2006138393 0.2340494792 0.3510742188
## 5  3.339844e-01 2.671875e-01 0.3816964286 0.4453125000 0.6679687500
## 6  4.824219e-01 3.859375e-01 0.5513392857 0.6432291667 0.9648437500
## 7  6.455078e-01 5.164063e-01 0.7377232143 0.8606770833 1.2910156250
## 8  8.007813e-01 6.406250e-01 0.9151785714 1.0677083333 1.6015625000
## 9  9.570313e-01 7.656250e-01 1.0937500000 1.2760416667 1.9140625000
## 10 9.218750e-01 7.375000e-01 1.0535714286 1.2291666667 1.8437500000
## 11 9.687500e-01 7.750000e-01 1.1071428571 1.2916666667 1.9375000000
## 12 1.000000e+00 8.000000e-01 1.1428571429 1.3333333333 2.0000000000
## 13 1.250000e+00 1.000000e+00 1.4285714286 1.6666666667 2.5000000000
## 14 8.750000e-01 7.000000e-01 1.0000000000 1.1666666667 1.7500000000
## 15 7.500000e-01 6.000000e-01 0.8571428571 1.0000000000 1.5000000000
## 16 5.000000e-01 4.000000e-01 0.5714285714 0.6666666667 1.0000000000
## 17 1.000000e+00 8.000000e-01 1.1428571429 1.3333333333 2.0000000000
##              17
## 1  9.155273e-05
## 2  3.204346e-03
## 3  6.323242e-02
## 4  1.755371e-01
## 5  3.339844e-01
## 6  4.824219e-01
## 7  6.455078e-01
## 8  8.007812e-01
## 9  9.570313e-01
## 10 9.218750e-01
## 11 9.687500e-01
## 12 1.000000e+00
## 13 1.250000e+00
## 14 8.750000e-01
## 15 7.500000e-01
## 16 5.000000e-01
## 17 1.000000e+00

SVMP

##              0           1           2           3           4           5
## 0     1.000000   0.4204545 0.004693049 0.002073061 0.001344477 0.001008944
## 1     2.378378   1.0000000 0.011161847 0.004930524 0.003197674 0.002399651
## 2   213.081081  89.5909091 1.000000000 0.441730166 0.286482558 0.214986911
## 3   482.378378 202.8181818 2.263825469 1.000000000 0.648546512 0.486692845
## 4   743.783784 312.7272727 3.490613902 1.541909458 1.000000000 0.750436300
## 5   991.135135 416.7272727 4.651445967 2.054683998 1.332558140 1.000000000
## 6  1101.837838 463.2727273 5.170979198 2.284177499 1.481395349 1.111692845
## 7  1217.729730 512.0000000 5.714865550 2.524428507 1.637209302 1.228621291
## 8  1300.756757 546.9090909 6.104515474 2.696548633 1.748837209 1.312390925
## 9  1383.783784 581.8181818 6.494165398 2.868668758 1.860465116 1.396160558
## 10 1217.729730 512.0000000 5.714865550 2.524428507 1.637209302 1.228621291
## 11  830.270270 349.0909091 3.896499239 1.721201255 1.116279070 0.837696335
## 12  774.918919 325.8181818 3.636732623 1.606454505 1.041860465 0.781849913
## 13 1328.432432 558.5454545 6.234398782 2.753922008 1.786046512 1.340314136
## 14  885.621622 372.3636364 4.156265855 1.835948005 1.190697674 0.893542757
## 15  885.621622 372.3636364 4.156265855 1.835948005 1.190697674 0.893542757
##               6            7            8            9           10          11
## 0  0.0009075746 0.0008212003 0.0007687832 0.0007226563 0.0008212003 0.001204427
## 1  0.0021585557 0.0019531250 0.0018284574 0.0017187500 0.0019531250 0.002864583
## 2  0.1933869702 0.1749822443 0.1638131649 0.1539843750 0.1749822443 0.256640625
## 3  0.4377943485 0.3961292614 0.3708444149 0.3485937500 0.3961292614 0.580989583
## 4  0.6750392465 0.6107954545 0.5718085106 0.5375000000 0.6107954545 0.895833333
## 5  0.8995290424 0.8139204545 0.7619680851 0.7162500000 0.8139204545 1.193750000
## 6  1.0000000000 0.9048295455 0.8470744681 0.7962500000 0.9048295455 1.327083333
## 7  1.1051805338 1.0000000000 0.9361702128 0.8800000000 1.0000000000 1.466666667
## 8  1.1805337520 1.0681818182 1.0000000000 0.9400000000 1.0681818182 1.566666667
## 9  1.2558869702 1.1363636364 1.0638297872 1.0000000000 1.1363636364 1.666666667
## 10 1.1051805338 1.0000000000 0.9361702128 0.8800000000 1.0000000000 1.466666667
## 11 0.7535321821 0.6818181818 0.6382978723 0.6000000000 0.6818181818 1.000000000
## 12 0.7032967033 0.6363636364 0.5957446809 0.5600000000 0.6363636364 0.933333333
## 13 1.2056514914 1.0909090909 1.0212765957 0.9600000000 1.0909090909 1.600000000
## 14 0.8037676609 0.7272727273 0.6808510638 0.6400000000 0.7272727273 1.066666667
## 15 0.8037676609 0.7272727273 0.6808510638 0.6400000000 0.7272727273 1.066666667
##             12           13          14          15
## 0  0.001290458 0.0007527669 0.001129150 0.001129150
## 1  0.003069196 0.0017903646 0.002685547 0.002685547
## 2  0.274972098 0.1604003906 0.240600586 0.240600586
## 3  0.622488839 0.3631184896 0.544677734 0.544677734
## 4  0.959821429 0.5598958333 0.839843750 0.839843750
## 5  1.279017857 0.7460937500 1.119140625 1.119140625
## 6  1.421875000 0.8294270833 1.244140625 1.244140625
## 7  1.571428571 0.9166666667 1.375000000 1.375000000
## 8  1.678571429 0.9791666667 1.468750000 1.468750000
## 9  1.785714286 1.0416666667 1.562500000 1.562500000
## 10 1.571428571 0.9166666667 1.375000000 1.375000000
## 11 1.071428571 0.6250000000 0.937500000 0.937500000
## 12 1.000000000 0.5833333333 0.875000000 0.875000000
## 13 1.714285714 1.0000000000 1.500000000 1.500000000
## 14 1.142857143 0.6666666667 1.000000000 1.000000000
## 15 1.142857143 0.6666666667 1.000000000 1.000000000

SVSP

##            2           3          4          5          6           7
## 2    1.00000  0.09135802 0.02617925 0.01359628 0.01003254 0.008575402
## 3   10.94595  1.00000000 0.28655660 0.14882411 0.10981562 0.093865884
## 4   38.19820  3.48971193 1.00000000 0.51935326 0.38322487 0.327564895
## 5   73.54955  6.71934156 1.92547170 1.00000000 0.73788865 0.630716934
## 6   99.67568  9.10617284 2.60943396 1.35521803 1.00000000 0.854758962
## 7  116.61261 10.65349794 3.05283019 1.58549731 1.16992046 1.000000000
## 8  113.87387 10.40329218 2.98113208 1.54826066 1.14244396 0.976514215
## 9  107.81982  9.85020576 2.82264151 1.46594806 1.08170644 0.924598269
## 10  85.33333  7.79588477 2.23396226 1.16021558 0.85610991 0.731767614
## 11  85.33333  7.79588477 2.23396226 1.16021558 0.85610991 0.731767614
## 12  36.90090  3.37119342 0.96603774 0.50171485 0.37020969 0.316440049
## 13  64.57658  5.89958848 1.69056604 0.87800098 0.64786696 0.553770087
## 14  18.45045  1.68559671 0.48301887 0.25085742 0.18510484 0.158220025
## 15  73.80180  6.74238683 1.93207547 1.00342969 0.74041938 0.632880099
##              8           9         10         11         12         13
## 2  0.008781646 0.009274733 0.01171875 0.01171875 0.02709961 0.01548549
## 3  0.096123418 0.101520722 0.12827280 0.12827280 0.29663086 0.16950335
## 4  0.335443038 0.354278075 0.44763514 0.44763514 1.03515625 0.59151786
## 5  0.645886076 0.682152406 0.86190878 0.86190878 1.99316406 1.13895089
## 6  0.875316456 0.924465241 1.16807432 1.16807432 2.70117188 1.54352679
## 7  1.024050633 1.081550802 1.36655405 1.36655405 3.16015625 1.80580357
## 8  1.000000000 1.056149733 1.33445946 1.33445946 3.08593750 1.76339286
## 9  0.946835443 1.000000000 1.26351351 1.26351351 2.92187500 1.66964286
## 10 0.749367089 0.791443850 1.00000000 1.00000000 2.31250000 1.32142857
## 11 0.749367089 0.791443850 1.00000000 1.00000000 2.31250000 1.32142857
## 12 0.324050633 0.342245989 0.43243243 0.43243243 1.00000000 0.57142857
## 13 0.567088608 0.598930481 0.75675676 0.75675676 1.75000000 1.00000000
## 14 0.162025316 0.171122995 0.21621622 0.21621622 0.50000000 0.28571429
## 15 0.648101266 0.684491979 0.86486486 0.86486486 2.00000000 1.14285714
##            14        15
## 2  0.05419922 0.0135498
## 3  0.59326172 0.1483154
## 4  2.07031250 0.5175781
## 5  3.98632813 0.9965820
## 6  5.40234375 1.3505859
## 7  6.32031250 1.5800781
## 8  6.17187500 1.5429688
## 9  5.84375000 1.4609375
## 10 4.62500000 1.1562500
## 11 4.62500000 1.1562500
## 12 2.00000000 0.5000000
## 13 3.50000000 0.8750000
## 14 1.00000000 0.2500000
## 15 4.00000000 1.0000000

TFTx

##            2          3         4          5          6          7          8
## 2   1.000000  0.2439024 0.1100811 0.06384409 0.04576108 0.03971572 0.03029337
## 3   4.100000  1.0000000 0.4513326 0.26176075 0.18762042 0.16283445 0.12420281
## 4   9.084211  2.2156611 1.0000000 0.57997312 0.41570328 0.36078595 0.27519133
## 5  15.663158  3.8202824 1.7242178 1.00000000 0.71676301 0.62207358 0.47448980
## 6  21.852632  5.3299101 2.4055620 1.39516129 1.00000000 0.86789298 0.66198980
## 7  25.178947  6.1412067 2.7717265 1.60752688 1.15221580 1.00000000 0.76275510
## 8  33.010526  8.0513479 3.6338355 2.10752688 1.51059730 1.31103679 1.00000000
## 9  36.378947  8.8729140 4.0046350 2.32258065 1.66473988 1.44481605 1.10204082
## 10 43.115789 10.5160462 4.7462341 2.75268817 1.97302505 1.71237458 1.30612245
## 11 41.768421 10.1874198 4.5979143 2.66666667 1.91136802 1.65886288 1.26530612
## 12 24.252632  5.9152760 2.6697567 1.54838710 1.10982659 0.96321070 0.73469388
## 13 21.557895  5.2580231 2.3731170 1.37634409 0.98651252 0.85618729 0.65306122
## 14 21.557895  5.2580231 2.3731170 1.37634409 0.98651252 0.85618729 0.65306122
## 15 43.115789 10.5160462 4.7462341 2.75268817 1.97302505 1.71237458 1.30612245
## 16 43.115789 10.5160462 4.7462341 2.75268817 1.97302505 1.71237458 1.30612245
##             9         10         11         12         13         14         15
## 2  0.02748843 0.02319336 0.02394153 0.04123264 0.04638672 0.04638672 0.02319336
## 3  0.11270255 0.09509277 0.09816028 0.16905382 0.19018555 0.19018555 0.09509277
## 4  0.24971065 0.21069336 0.21748992 0.37456597 0.42138672 0.42138672 0.21069336
## 5  0.43055556 0.36328125 0.37500000 0.64583333 0.72656250 0.72656250 0.36328125
## 6  0.60069444 0.50683594 0.52318548 0.90104167 1.01367188 1.01367188 0.50683594
## 7  0.69212963 0.58398437 0.60282258 1.03819444 1.16796875 1.16796875 0.58398438
## 8  0.90740741 0.76562500 0.79032258 1.36111111 1.53125000 1.53125000 0.76562500
## 9  1.00000000 0.84375000 0.87096774 1.50000000 1.68750000 1.68750000 0.84375000
## 10 1.18518519 1.00000000 1.03225806 1.77777778 2.00000000 2.00000000 1.00000000
## 11 1.14814815 0.96875000 1.00000000 1.72222222 1.93750000 1.93750000 0.96875000
## 12 0.66666667 0.56250000 0.58064516 1.00000000 1.12500000 1.12500000 0.56250000
## 13 0.59259259 0.50000000 0.51612903 0.88888889 1.00000000 1.00000000 0.50000000
## 14 0.59259259 0.50000000 0.51612903 0.88888889 1.00000000 1.00000000 0.50000000
## 15 1.18518519 1.00000000 1.03225806 1.77777778 2.00000000 2.00000000 1.00000000
## 16 1.18518519 1.00000000 1.03225806 1.77777778 2.00000000 2.00000000 1.00000000
##            16
## 2  0.02319336
## 3  0.09509277
## 4  0.21069336
## 5  0.36328125
## 6  0.50683594
## 7  0.58398438
## 8  0.76562500
## 9  0.84375000
## 10 1.00000000
## 11 0.96875000
## 12 0.56250000
## 13 0.50000000
## 14 0.50000000
## 15 1.00000000
## 16 1.00000000

Plotting BAMM data

BPP

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/BPP/event_dataBPP.txt", burnin = 0.1, type= 'trait')
plot.bammdata(edata, lwd = 1.5, legend = T,labels = T)

distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 87 plots

cred_Shift$number.distinct
## [1] 91
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  91
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.22738631  0.2273863          3
##          2 0.09245377  0.3198401          3
##          3 0.09245377  0.4122939          4
##          4 0.05897051  0.4712644          4
##          5 0.03648176  0.5077461          4
##          6 0.03098451  0.5387306          4
##          7 0.03048476  0.5692154          4
##          8 0.02948526  0.5987006          5
##          9 0.02898551  0.6276862          3
## 
## ...omitted 82 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5, method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

CRISP

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/CRISP/edata_CRISP.txt", burnin = 0.1, type= 'trait')
plot.bammdata(edata, lwd = 1.5, legend = T,labels = T)

distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 8 plots

cred_Shift$number.distinct
## [1] 12
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  12
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.37607799  0.3760780          1
##          2 0.16535433  0.5414323          1
##          3 0.12773403  0.6691664          1
##          4 0.09423822  0.7634046          2
##          5 0.06411699  0.8275216          1
##          6 0.03212098  0.8596425          2
##          7 0.02887139  0.8885139          2
##          8 0.02737158  0.9158855          2
##          9 0.01449819  0.9303837          2
## 
## ...omitted 3 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks", method = 'polar') + addBAMMshifts(best, cex = 2.5, method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

CTL

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/CTL/event_dataCTL.txt", burnin = 0.1, type= 'trait')
plot.bammdata(edata, lwd = 1.5, legend = T,labels = T)

distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 139 plots

cred_Shift$number.distinct
## [1] 143
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  143
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.18986779  0.1898678          3
##          2 0.12143095  0.3112987          2
##          3 0.07076991  0.3820687          3
##          4 0.06354849  0.4456172          3
##          5 0.04655038  0.4921675          3
##          6 0.04243973  0.5346073          2
##          7 0.03188535  0.5664926          4
##          8 0.02210865  0.5886013          3
##          9 0.01733141  0.6059327          4
## 
## ...omitted 134 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5,method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

GF

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/GF/event_dataGF.txt", burnin = 0.1, type= 'trait')
plot.bammdata(edata, lwd = 1.5, legend = T,labels = T)

distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 110 plots

cred_Shift$number.distinct
## [1] 114
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  114
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.25452727  0.2545273          2
##          2 0.17298078  0.4275081          2
##          3 0.05799356  0.4855016          2
##          4 0.05754916  0.5430508          2
##          5 0.04221753  0.5852683          2
##          6 0.03910677  0.6243751          4
##          7 0.02377514  0.6481502          3
##          8 0.02377514  0.6719253          3
##          9 0.02321964  0.6951450          3
## 
## ...omitted 105 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5,method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

KSPI

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/KUN/event_dataKUN.txt", burnin = 0.1, type= 'trait')
plot.bammdata(edata, lwd = 1.5, legend = T,labels = T)

distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 317 plots

cred_Shift$number.distinct
## [1] 321
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  321
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.23447069  0.2344707          3
##          2 0.09123860  0.3257093          3
##          3 0.05436820  0.3800775          4
##          4 0.03237095  0.4124484          3
##          5 0.03187102  0.4443195          4
##          6 0.02824647  0.4725659          4
##          7 0.02287214  0.4954381          3
##          8 0.01987252  0.5153106          3
##          9 0.01424822  0.5295588          3
## 
## ...omitted 312 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5,method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

LAAO

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/LAAO/event_dataLAAO.txt", burnin = 0.1, type= 'trait')
plot.bammdata(edata, lwd = 1.5, legend = T,labels = T)

distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 42 plots

cred_Shift$number.distinct
## [1] 46
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  46
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.23058980  0.2305898          2
##          2 0.11262912  0.3432189          2
##          3 0.11262912  0.4558481          2
##          4 0.06864379  0.5244918          0
##          5 0.05598134  0.5804732          3
##          6 0.04798401  0.6284572          2
##          7 0.04265245  0.6711096          1
##          8 0.02765745  0.6987671          2
##          9 0.02465845  0.7234255          3
## 
## ...omitted 37 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5,method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

vPLA2

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/vPLA2ePLA2/event_data_vPLA2.txt", burnin = 0.1, type= 'trait')
distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 9 plots

cred_Shift$number.distinct
## [1] 13
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  13
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.56805399  0.5680540          2
##          2 0.16222972  0.7302837          3
##          3 0.04911886  0.7794026          3
##          4 0.04424447  0.8236470          1
##          5 0.02612173  0.8497688          2
##          6 0.02124734  0.8710161          4
##          7 0.02037245  0.8913886          1
##          8 0.01687289  0.9082615          2
##          9 0.01399825  0.9222597          3
## 
## ...omitted 4 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5,method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

ePLA2

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/vPLA2ePLA2/event_data_ePLA2.txt", burnin = 0.1, type= 'trait')
distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 352 plots

cred_Shift$number.distinct
## [1] 356
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  356
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.14885639  0.1488564          3
##          2 0.10323710  0.2520935          3
##          3 0.06249219  0.3145857          4
##          4 0.05311836  0.3677040          3
##          5 0.02987127  0.3975753          3
##          6 0.02849644  0.4260717          3
##          7 0.02412198  0.4501937          4
##          8 0.02374703  0.4739408          4
##          9 0.01887264  0.4928134          3
## 
## ...omitted 347 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5,method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

SVMP

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/SVMP/event_dataSVMP.txt", burnin = 0.1, type= 'trait')
distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 0 plots

cred_Shift$number.distinct
## [1] 4
summary(cred_Shift)
## 
##  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.84672391  0.8467239          2
##          2 0.05599341  0.9027173          3
##          3 0.03623103  0.9389484          3
##          4 0.02929067  0.9682390          3
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5,method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

SVSP

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/SVSP/event_dataLAAO.txt", burnin = 0.1, type= 'trait')
distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 36 plots

cred_Shift$number.distinct
## [1] 40
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  40
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.18821315  0.1882132          4
##          2 0.16303964  0.3512528          3
##          3 0.10692860  0.4581814          4
##          4 0.09398894  0.5521703          3
##          5 0.04246559  0.5946359          2
##          6 0.04070109  0.6353370          3
##          7 0.03093754  0.6662746          5
##          8 0.03034937  0.6966239          5
##          9 0.02623221  0.7228561          4
## 
## ...omitted 31 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T,breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5,method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

TFTx

edata<-getEventData(tree, eventdata="./Analysis_scaled/ControlFiles_scaleddata/TFTx/event_dataTFTx.txt", burnin = 0.1, type= 'trait')
distShiftconfigs<-distinctShiftConfigurations(edata, threshold = 0.01, expectedNumberOfShifts = 1)
plot(distShiftconfigs$frequency)

cred_Shift = credibleShiftSet(edata, expectedNumberOfShifts = 1, threshold = 5, set.limit = 0.95)
plot.credibleshiftset(cred_Shift, lwd=1.7, plotmax = 4)
## Omitted 139 plots

cred_Shift$number.distinct
## [1] 143
summary(cred_Shift)
## 
##  95 % credible set of rate shift configurations sampled with BAMM
## 
## Distinct shift configurations in credible set:  143
## 
## Frequency of 9 shift configurations with highest posterior probability:
## 
## 
##    rank     probability cumulative  Core_shifts
##          1 0.17520620  0.1752062          2
##          2 0.17270682  0.3479130          3
##          3 0.07273182  0.4206448          2
##          4 0.05723569  0.4778805          4
##          5 0.04098975  0.5188703          4
##          6 0.03924019  0.5581105          3
##          7 0.02624344  0.5843539          3
##          8 0.02549363  0.6098475          3
##          9 0.01599600  0.6258435          4
## 
## ...omitted 134 additional distinct shift configurations
## from the credible set. You can access the full set from your 
## credibleshiftset object
dev.off()
## null device 
##           1
best<-getBestShiftConfiguration(edata, expectedNumberOfShifts = 1)
## 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.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks",method = 'polar') + addBAMMshifts(best, cex = 2.5,method = 'polar')
## Error in plot.bammdata(best, lwd = 2, legend = T, breaksmethod = "jenks", : non-numeric argument to binary operator

Rate through time

Since the early burst model states that branching events account for most of the phenotypic diversity amongst species (Simpson 1953.; Ingram et al. 2012), we expect the rate of trait evolution to be high immediately after major lineages diverge. To identify these increases we estimated the evolutionary rate of each toxin family after the split of the three venomous families. Our estimates of ‘rate through time’ revealed that toxin families show unique evolutionary rates and rate dynamics in each venomous snake family (Fig. 4)

BPP

st <- max(branching.times(tree))
plotRateThroughTime(edataBPP, intervalCol="blue", avgCol="blue", start.time=st, node=97, ylim=c(0,.08),cex.axis=1)+
text(x=30, y= 0.008, label="Colubrids only", font=4, cex=2.0, pos=4)
plotRateThroughTime(edataBPP, intervalCol="red", avgCol="red", start.time=st, node=54, ylim=c(0,.08),cex.axis=1)+
text(x=30, y= 0.008, label="Vipers only", font=4, cex=2.0, pos=4)
plotRateThroughTime(edataBPP, intervalCol="green", avgCol="green", start.time=st, node=85, ylim=c(0,.08),cex.axis=1)+
text(x=30, y= 0.008, label="Elapids only", font=4, cex=2.0, pos=4)

## integer(0)

## integer(0)

## integer(0)

CRISP

## integer(0)

## integer(0)

## integer(0)

CTL

## integer(0)

## integer(0)

## integer(0)

GF

## integer(0)

## integer(0)

## integer(0)

KSPI

## integer(0)

## integer(0)

## integer(0)

LAAO

## integer(0)

## integer(0)

## integer(0)

vPLA2

## integer(0)

## integer(0)

## integer(0)

ePLA2

## integer(0)

## integer(0)

## integer(0)

SVMP

## integer(0)

## integer(0)

SVSP

## integer(0)

## integer(0)

TFTx

## integer(0)

## integer(0)

Analytical considerations

A major issue with these types of analysis is species sampling. While our dataset spans the three families of venomous snakes, it does undersample Colubrids and most Australian elapids. Despite this, the nature of our analysis helps alleviate some sample bias. For instance, if we had a larger phylogeny and trait data, we would be able to identify more potential rate shifts dating further back in time and along more branches in the phylogeny, however, we would also obtain the same shifts configurations we observe in our data (although the relative abundance in the posterior could be lower). Rather than being interested in the precise locations of shifts in evolutionary dynamics, we are interested to identify evidence of non-zero rate shifts. Even using the most conservative estimate (expectedshifts = 1, see methods) Bayes factor tests shows that we always obtain shifts that are higher than expected by chance, therefore, it is unlikely that we would get a lower number of shifts that we get from our estimates if we use a larger phylogeny. Furthermore, the modelling approach developed by Landis and Schraiber (Landis and Schraiber 2017) has low false positive rates for selecting jump processes, even for subsets of small trees that are different from a true tree. Therefore, while a much larger dataset would undoubtedly improve the resolution of our macroevolutionary dynamics, we believe our analysis provides evidence of similar trends using a concise dataset.

In the estimation of evolutionary rate dynamics, two rate configurations are potentially misrepresented in the phylorate plots. One of them is that of TFTx in vipers (at the ancestral branch of Crotalus simus and Crotalus durissus) and the other is the occurrence of a shift for ePLA2 at the same branch. The shift for TFTx can be explained by the primary source of data that reports levels of TFTx in the C. simus venom transcriptome (Durban et al. 2011). Despite being sister species, the gene expression profiles of venom toxins in C. simus and C. durissus vary extensively from each other. This is not unusual, and has been previously reported in many other closely related venomous snakes (Chippaux, Williams, and White 1991). When observing the untransformed data we find that C. durissus has more structural zeros as compared to C. simus (Supplementary Table 1). After normalization, the difference between toxins levels becomes more prominent. While this is appropriate for toxins that are indeed differentially expressed, for toxins that are absent, their normalized value gets inflated to represent a quantitative difference between the two based on the values of other toxins. We believe that the shift in ePLA2 at this branch is an artefact of normalization. While this would cause BAMM to place a rate shift along the branch to explain this stark difference, since the normalization is done within species, this would not affect estimates in other species and locations of other shifts. If the normalization was causing spurious differences throughout the dataset we would see it in other closely related species as well; Protobothops, which like Crotalus do not have ePLA2, or Naja, and Micrurus, which do not have vPLA2. Since we do not see spurious shifts in these lineages, we can be fairly certain that the normalization only had an effect in the two Crotalus species. Thus, this would not affect our interpretations as shifts in evolutionary rates are significant and common to all toxin families.

Chippaux, J. P., V. Williams, and J. White. 1991. “Snake Venom Variability: Methods of Study, Results and Interpretation.” Toxicon: Official Journal of the International Society on Toxinology 29 (11): 1279–1303.

Durban, Jordi, Paula Juárez, Yamileth Angulo, Bruno Lomonte, Marietta Flores-Diaz, Alberto Alape-Girón, Mahmood Sasa, et al. 2011. “Profiling the Venom Gland Transcriptomes of Costa Rican Snakes by 454 Pyrosequencing.” BMC Genomics 12 (May): 259.

Landis, Michael J., and Joshua G. Schraiber. 2017. “Pulsed Evolution Shaped Modern Vertebrate Body Sizes.” Proceedings of the National Academy of Sciences of the United States of America 114 (50): 13224–29.