Tree incongruence: genes whose own alignment rejects the species tree

Test every gene alignment against the species tree and find the genes whose own sites reject it, with the false discovery rate controlled across genes.

What it doesCopy Anchor Link

Choose the gene alignments, one aligned FASTA file per gene with one record per species, and give the species tree in Newick with its tips named as the records are. For each gene, the species tree and every tree one nearest-neighbour interchange away from it are fitted to the gene's alignment by maximum likelihood, and the Shimodaira–Hasegawa test asks whether the species tree is significantly worse than the best of them. One false discovery rate then covers all the genes. The result lists the genes that reject the species tree, smallest p first, and the whole table of genes, each with its number of sites, how far the species tree falls below its best tree, its p and its q. Beside it are the alternative trees that genes prefer, and for each internal branch of the species tree how many genes prefer each of its two rearrangements.

When you’d use thisCopy Anchor Link

• Finding the genes that tell a different story. Among closely spaced species a share of genes follows another tree because ancestral variation sorted differently along each lineage. Genes that significantly reject the species tree are the ones to read first.
• Asking whether discordance is only lineage sorting. Incomplete lineage sorting alone makes the two rearrangements of a branch equally common. A significant excess of one is what gene flow between the lineages it joins predicts, and also a bias shared by those lineages' data.
• Checking a set of orthologue alignments. A paralog taken for an orthologue, or a stretch that two species' gene models share and the others lack, makes a gene reject the species tree strongly.

What this computesCopy Anchor Link

For every gene:

• sites and site patterns of the alignment
• the fitted GTR + Γ model: six exchangeabilities, four base frequencies and the Gamma shape of rate variation among sites, fitted with the branch lengths on the species tree
• the log-likelihood of the species tree and of each nearest-neighbour rearrangement (two per internal branch), each with its branch lengths refitted under the gene's model
• Δ ln L, how far each tree falls below the gene's best tree
• the Shimodaira–Hasegawa p of every tree, by RELL resampling of the sites. Each replicate is measured from its best tree, as Shimodaira and Hasegawa define the test, or from its worst tree, as phangorn computes it. Trees within 0.001 log-likelihood units of the best count as tied with it (p 1)
• the q-value of the species tree's p, adjusted across the genes by Benjamini–Hochberg (independent or positively dependent tests) or Benjamini–Yekutieli (any dependence, less power), and the rejection call at the false discovery rate you set
• whether the gene's best tree is the species tree, and which tree it is if not

Across the genes:

• the number of genes that reject the species tree and the number that are discordant (prefer another tree, rejected or not)
• each alternative tree some gene prefers, with how many genes prefer it and how many of those reject the species tree
• per internal branch of the species tree, the counts of genes preferring each of its two rearrangements and the exact binomial p of their equality

Every result states, in words generated from the run, the counts around each branch and whether they depart from the equality lineage sorting predicts, the causes a rejection can have, which genes have fewer than 300 sites (too few for the test to have power), and, below 1,000 replicates, the smallest p the test can give.

Worked exampleCopy Anchor Link

122 human genes and their one-to-one orthologues in chimpanzee, gorilla, orangutan and rhesus macaque (Ensembl Compara coding alignments, a seeded sample of protein-coding genes), against the Open Tree of Life species tree ((human, chimpanzee), gorilla), orangutan, macaque. 1,000 replicates, from the best tree, Benjamini–Hochberg at q ≤ 0.05:

• 48 genes prefer another tree, and 5 reject the species tree: MAP4K2, DUSP18, UCHL1, GLRX2 and DNAJB11, each with the species tree 57 to 275 log-likelihood units below its best tree.
• Around the human–chimpanzee branch, 21 genes prefer gorilla with chimpanzee and 20 prefer gorilla with human (binomial p 1): the two lineage-sorting alternatives, equally often.
• Read the five rejections before naming a cause. In DUSP18 gorilla and chimpanzee share a 213-column gap the other three do not, and in DNAJB11 gorilla and orangutan carry a 45-column block the other three lack: differences between gene models, not between genealogies.

ControlsCopy Anchor Link

List decides how much of the result the table shows: Rejecting or Every gene.

LimitationsCopy Anchor Link

• A rejection is not gene flow or transfer. Lineage sorting, hidden paralogy, alignment error and selection all make a gene fit another tree better.
• Only trees one nearest-neighbour interchange from the species tree are tested. A gene whose best tree is further away still rejects the species tree when one neighbour fits better, but its own best tree is not reported.
• The model is fitted on the species tree and kept for the rearrangements, whose branch lengths alone are refitted, so a rearrangement is never favoured by a model fitted to it.
• The SH test is conservative, and more so the more trees are compared; a gene with few sites rarely rejects anything.
• Measured from the worst tree, as phangorn computes it, trees whose likelihoods tie (a gene with no informative site) get p 0 and are reported as rejecting. The default, from the best tree, counts a tie as p 1.
• Coding genes sort more completely than the genome around them: in the gorilla genome, 30% of the sequence is closer to human or chimpanzee than those two are to each other, and fewer coding genes are.

Where the numbers come fromCopy Anchor Link

Every gene's log-likelihoods are checked against phangorn's pml and optim.pml on the 122 great-ape genes, fitted the same way, and the SH arithmetic against phangorn's SH.test on phangorn's own per-site likelihoods. The form measured from the best tree is the one IQ-TREE computes.

AutomationCopy Anchor Link

The same test is a workflow node and a tool the assistant can call, with the same settings and defaults as the panel; each alignment is the name of a FASTA file in the project.