Computational protocol: A Phylogenomic Approach to Vertebrate Phylogeny Supports a Turtle-Archosaur Affinity and a Possible Paraphyletic Lissamphibia

Similar protocols

Protocol publication

[…] Markers used in this study are single-copy, orthologous, protein-coding genes . The single-copy nature of the markers was checked both during marker development and after data collection. During marker development, markers were pre-screened using BLAST to compare with the high-coverage, well-annotated chicken genome and gene families were identified in Ensembl and Metazome. After data collection, to identify and remove paralogous genes, preliminary phylogenetic analyses using RAxML were performed for each gene. Each of the trees was analyzed individually by eye for erroneous phylogenetic relationships (e.g. mammal species more closely related to fish) and signatures of gene duplications (i.e. replicated tree topologies within the larger tree). Sequences with erroneous positions were removed, and when gene duplications were detected, the sub-tree that included newly sequenced data was retained for subsequent analyses. We obtained new sequences according to the methods of Fong & Fujita (2011) . Briefly, we used cDNA preparations as template to amplify the protein-coding genes using conserved primers. Amplicons were sequenced using ABI3730 chemistry, and sequences were edited using Geneious 5 (Biomatters Ltd.) and aligned using MUSCLE . [...] Sequences were combined into two main categories of datasets: individual genes and concatenations. Individual datasets for the 75 genes consisted of orthologous sequences from online genomes and the new samples. Combining individual genes using a Perl script (available upon request) produced the concatenated datasets. We compiled seven different concatenated datasets: 1) All taxa-75 genes, 2) 16 taxa-75 genes, 3) All taxa-31 genes (turtle), 4) All taxa-26 genes (Lissamphibia), 5) 16 taxa-31 genes (turtle), 6) 16 taxa-26 genes (Lissamphibia), 7) slow genes (removal of fastest 25% of genes). For the 16-taxon datasets, the vertebrate group Crocodilia is represented by two individuals of the species Alligator mississippiensis (from the EST database and a new sample). We combined the data from both individuals to minimize missing data; this approach is justified, as when there were data from both individuals for a marker, data were identical. The dataset of reduced loci for all taxa was used when evaluating the specific phylogenetic questions (turtle and Lissamphibia). Loci without representatives of all the focal groups were removed, leaving 31 genes for the turtle analysis and 26 genes for the lissamphibian analysis.The standard nucleotide (NUCL) dataset was transformed to three data-types using the following methods. AA was translated in Geneious 5 (Biomatters Ltd.), 3rd codon positions were removed for N12 using MacClade v4.08 , and DEGEN1 was converted using a Perl script .The rate of evolution for each of the 75 genes was calculated by computing tree length and averaging branch lengths using an R script . Based on the shape of the frequency histogram (), we drew a cut-off of average branch length of 0.04, which denoted the top 25% fastest genes (19 genes) for removal. The names of these genes are: DSCR3, EXOC2, GAPDH, GDE1, GNAS, HPD, H2AFY2, IFT57, MAT2B, OAT, OAZ1, PPM1A, PPTC7, PSAT1, SEC13, SGK1, TAT, UBE2J2, and XPOT () . [...] All datasets were subject to maximum likelihood analyses using RAxML , and a subset of datasets were also subject to Bayesian analyses using MrBayes .Since this study deals with a complex, multi-gene dataset, we explored heterogeneous processes of molecular evolution through partitioning the data. Tests of alternative partitioning strategies were performed on the NUCL dataset only, as the N12 dataset is a subset of the NUCL dataset, and the DEGEN1 and AA datasets have information on codon position integrated into gene partitions. For RAxML analyses of the NUCL dataset, three different partitioning strategies were tested: by gene (75 partitions), by gene and 1st+2nd and 3rd codon position (150 partitions), and by gene and codon position (225 partitions). Likelihood ratio tests selected the 150 partitions as the best partitioning strategy.RAxML nucleotide analyses used the GTRGAMMA model of evolution for tree inference and bootstrapping (1,000 replicates). RAxML amino acid analyses employed the protein gamma model of evolution and the appropriate model of protein evolution selected using ProtTest v2.4 , with empirical amino-acid frequencies for both tree inference and bootstrapping (1,000 replicates). All concatenated datasets were partitioned according to the optimal partitioning strategy. RAxML v7.2.5 and v7.2.6 , analyses were run locally and on the CIPRES portal v2.2 .For individual gene analyses, the support values of clades are generally very low, since these relatively short genes (average length is ∼450bp) were used to infer the entire vertebrate phylogeny. However, to understand and summarize the phylogenetic signal for each gene, we classified them based on preferred topology (see ) irrespective of nodal support.Bayesian analyses were only run on individual genes and 16-taxon datasets, as the computational burden for the larger datasets would require extremely long analysis times to achieve stationarity (i.e. >2000 hours). When both RAxML and MrBayes analyses were run, preferred topologies were almost identical, so results should not be compromised by reporting only inferences from RAxML. MrBayes v3.1.2 and v3.2 analyses were run locally and using the [email protected] resource at Cornell University ( All analyses were run with four chains for 10 million generations. Appropriate models of DNA substitution for each partition were selected using MrModeltest v2.3 , and amino acid substitution models were the same as those used in RAxML analyses. Burn-in of MCMC chains was evaluated using the online program AWTY, examining cumulative plots of posterior probabilities of the 20 most variable splits . [...] Approximately unbiased topology tests (AU tests) were used to test whether sub-optimal trees were significantly worse than the maximum likelihood tree. AU tests were performed to compare the five turtle and four lissamphibian alternative hypotheses for each of the different datasets. Constrained RAxML analyses were run for each of the different topologies using the GTRGAMMA model of sequence evolution, and per-site log likelihoods calculated. These per-site log likelihoods were then input into the program CONSEL . [...] Scripts for DFA analyses of the NUCL dataset were written using the R language and rely on the SeqinR and MASS libraries. These scripts are available from the authors upon request. Per-site log-likelihood scores (LLS) were calculated for each of the constrained phylogenies pertaining to a relevant hypothesis (five turtle positions, four Lissamphibia hypotheses) using RAxML. Site-wise GC content (%GC) and proportion of missing data (%missing) were computed for major clades with potential sister-relationships for turtle placement (turtles, archosaurs, crocodilians, reptiles excluding turtles, and amniotes excluding turtles) and lissamphibian relationships (amniotes, caecilians, caecilians+salamanders, frogs+salamanders, and lissamphibians). Site-specific rates of evolution (site-rates) were calculated for each of the nine constrained phylogenies using HyPhy under a GTR model of sequence evolution with model parameters estimated independently for each phylogenetic hypothesis. These rates were then averaged across the five turtle and four lissamphibian topologies. %GC and %missing were both calculated for each site in the NUCL data matrix, averaged across all taxa in each clade of interest.DFA (from the MASS library) was run with preferred topology as the predicted category, and %GC (for relevant clades), %missing (for relevant clades), and site-rates as predictor variables in one case, and without %GC in another case. Posterior probabilities from the DFA were calculated using leave-one-out cross-validation and normalized with prior probabilities (posterior/prior ratio). The prior probability of assignment to any particular topology was simply the proportion of sites in the alignment preferring that topology. Two different types of DFA were tested to maximize the predictive power of our analysis: linear discriminant analysis (LDA), and quadratic discriminant analysis (QDA). Comparisons of average posterior/prior ratios show that QDA performed best (Turtle: LDA = 1.483, QDA = 1.693; Lissamphibian: LDA: 1.219, QDA = 1.572). […]

Pipeline specifications