Comparative methods in the genomic era

Every modern biologist should know about phylogenetic comparative methods, even if they are not familiar with the term. The idea is that we cannot compare biological traits without taking into account the evolutionary relations between the species/populations used as "sample points". This is because individuals are not independent samples of whatever variables we are trying to correlate (they are "pseudoreplicates").

However, not everybody is an up-to-date biologist, as we see in the figure below: birds with similar body masses would also have similar flight speeds not because these two traits are related, but because these birds are closely related. That is, not long ago they were a single species (with a single value for each trait). The solution to avoid these spurious correlations is, in abstract terms, to think about what trait values do we expect for their ancestors. Hopefully we will have this in mind in upcoming large-scale genomic analyses, where sophisticated statistical algorithms for ever-increasing amounts of data might obfuscate their biological appropriateness.

from doi:10.1063/1.4886855

Here I include a few links to good comments and papers discussing comparative methods, starting with Felsenstein's reaction to the paper above (which is about airplanes, by the way). I also include at the bottom the slides of an impromptu journal club I presented a few years ago about things we can model over a tree, with strong emphasis on phylogenetic contrasts.

Physicists and engineers decide how to analyze evolution (by Joe Felsenstein):

They make allometric plots of features of new airplane models, log-log plots over many orders of magnitude. The airplanes show allometry: did you know that a 20-foot-long airplane won’t have 100-foot-long wings? That you need more fuel to carry a bigger load?

But permit me a curmudgeonly point: This paper would have been rejected in any evolutionary biology journal. Most of its central citations to biological allometry are to 1980s papers on allometry that failed to take the the phylogeny of the organisms into account. The points plotted in those old papers are thus not independently sampled, a requirement of the statistics used. (More precisely, their error residuals are correlated).

Steps towards understanding comparative methods – Wainwright Lab:

Felsenstein elaborates on his method of calculating standardized contrasts (phylogenetically independent contrasts) to help overcome the non-independence of character traits.  These contrasts are basically the differences between trait values of species pairs weighted by the evolutionary change separating them; they are estimates of the rate of change over time. A common use of standardized contrasts is to look for correlation in this rate between two traits; if standardized contrasts of traits X and Y are compared in a regression analysis, a linear trend suggests correlated rates of evolution between the traits.

Strong phylogenetic inertia on genome size and transposable element content among 26 species of flies | Biology Letters doi:10.1098/rsbl.2016.0407:

To date, quantifying the importance of phylogenetic inertia in TE content distribution remains a key question as the dynamic of TE accumulation is still poorly understood. Here, we analysed the evolution of genome size and genomic TE content in 26 Drosophila, using a phylogenetic framework. We estimated genomic TE content using a de novo TE assembly approach, tested the correlation between TE content and genome size among closely related species and finally estimated the phylogenetic inertia.
(...)
Comparative analyses were performed using ape [21], nlme [22] and phytools [23] packages in R. Ancestral trait reconstruction of genome size was calculated using phylogenetic independent contrasts. We tested the phylogenetic signal using Pagel's λ [24]. Best-fitting model to the trait evolution and its covariance structure was tested among (i) absence of phylogenetic signal, (ii) neutral Brownian motion and (iii) constrained evolution Ornstein–Uhlenbeck (OU) models using generalized least squares (GLS) and selected according to minimum Akaike information criterion (AIC).

Controlling for non-independence in comparative analysis of patterns across populations within species | Philosophical Transactions of the Royal Society B: Biological Sciences doi:10.1098/rstb.2010.0311:

How do we quantify patterns (such as responses to local selection) sampled across multiple populations within a single species? Key to this question is the extent to which populations within species represent statistically independent data points in our analysis. Comparative analyses across species and higher taxa have long recognized the need to control for the non-independence of species data that arises through patterns of shared common ancestry among them (phylogenetic non-independence), as have quantitative genetic studies of individuals linked by a pedigree.

The Unsolved Challenge to Phylogenetic Correlation Tests for Categorical Characters | Systematic Biology doi:10.1093/sysbio/syu070:

When comparative biologists observe that animal species living in caves also tend to have reduced eyes, they may see such correlation as evidence that the traits are adaptively or functionally linked: for instance, selection to maintain eye function is relaxed when light is unavailable.
(...)
However, the last few decades have taught us that among-species correlative tests should take into account evolutionary relationships (Felsenstein 1985; Ridley 1989; Harvey and Pagel 1991). If phylogeny is not taken into account, an interpreted correlation may have a trivial explanation different from the biological relationship we claim. There is a correlation among species in the distribution of fur and bones in the middle ear—species with fur also have three bones in the middle ear, and vice versa. These two traits are characteristics of mammals, and absent outside the mammals. Using their shared distribution as evidence of an interesting biological relationship between fur and middle ear bones would be considered a mistake, however, for reasons understood long ago by Darwin (1872)

Walking on a tree from Leonardo de Oliveira Martins

Bar charts must start at zero (or something)

The other day I was mentioning to a colleague that a bar chart should start at zero, and I may have given the impression that it was just my personal taste. It is not. It is a universal standard in statistical visualisation. However, since it is very easy to get distracted or to misinterpret it, I'll summarise it here very quickly, and include some links below if you want to read more.

1. It is OK to set the Y axis to whatever range you want (min an max are good choices). However, in this case, do not use bar charts. Use mean-and-error plots, boxplots, lines, scatterplots, violin plots, etc.
2. Bar charts start at the bottom of the plot for a reason: you want to easily compare the heights between bars. Our eyes try to find the proportion of increase/decrease between them. Then, for example, if your Y axis starts at  a higher value, a 0.1% increase may look twice the height.
3. The Y axis doesn't always need to start at zero, but it must be a value that makes sense as a baseline. Therefore, if you are comparing ratios maybe one is a better baseline, or for a time series the baseline can be the value at day one (as in climate change graphs). The idea is that there are "natural" baseline values, such that the proportions between bars are meaningful.
4. In all cases, ask yourself if a bar chart is really necessary, and if you can justify it against other options. This may help you pay attention to the interpretation of the baseline. 

I believe that the same rules apply to an impulse plot or a stem plot -- like the ones used in ACF plots -- but I am not so sure. Anyway, here are some links:

It’s OK not to start your y-axis at zero — Quartz

Of course column and bar charts should always have zeroed axes, since that is the only way for the visualization to accurately represent the data. Bar and column charts rely on bars that stretch to zero to accurately mirror the ratios between data points. Truncating the axis breaks the relationship between the size of the rectangle and the value of the data. There is no debating this one (except for a few exceptions).

(The exception linked above is about a ratio, where the natural baseline is one instead of zero.)

Bar Chart Baselines Start at Zero | FlowingData

The main argument for bar charts without a zero baseline is this: There’s no point in extending the range of the value axis if the range of the data never includes zero. Ok.
Now instead of weight, let’s look at height. I think we can agree that it’s difficult to find people who are zero inches tall. I’m 70 inches tall, and my son is half my height at 35 inches. The bar chart on the left shows the comparison with a zero baseline, and as expected, the bar for me is twice the length of the bar for my son. On the right, I take it to the extreme and set the baseline to 35 inches, and the bar for my son disappears.

Maybe the latter communicates that I’m much taller, but the magnitude is infinitely exaggerated.

(The link above also provides an elegant solution if you insist on using barplots: explicitly model your baseline such that you can compare all values to it.)

Kick the bar chart habit : Nature Methods : Nature Research

Instead of comparing to an abstract zero level, scientists often compare multiple experimental samples to one another. Because the samples are usually generated from populations with a potentially large and irregular underlying variation, graphing their means using bar charts misleadingly assigns importance to the distance of the means from zero and poorly represents the distribution of the data used to calculate the means. Instead of bar charts, mean-and-error plots and box plots should be used for statistical sample data.

Zero is zero - Statistical Modeling, Causal Inference, and Social Science

The idea is that the area of the bar represents “how many” or “how much.” The bar has to go down to 0 for that to work. You don’t have to have your y-axis go to zero, but if you want the axis to go anywhere else, don’t use a bar graph, use a line graph. Usually line graphs are better anyway.  I’m sure this is all in a book somewhere.

10 Ways Charts Can Lie, Cheat & Lead Astray | QVDesign

Here we’re looking at exactly the same basic dataset covering Visitor numbers per year, the chart on the left seems very placid and visitor numbers appear to be relatively steady year on year; not bad but not great, whilst the chart on the right looks much more positive; visitor numbers took a huge leap in 2010; cue extra investment and pay rises all round! The charts show identical data so what’s going on?

The ONLY difference between these 2 charts is that the one on the left has the Properties > Axis > ‘Forced 0’ checked for the Y-Axis whilst the other doesn’t; one check box and we’re completely changing the message of this chart. Now of course if you spend a little time looking at this chart and read off the numbers it becomes clear that the increase isn’t all that marked, however; (...)  To accentuate the effect even further you can use the ‘Static Min/Max’ settings to zero in even more on the difference making the smaller value seem even smaller and the increase all the larger.

How to summarise a collection of trees that came from a Bayesian analysis

After running a Bayesian phylogenetic analysis we are usually left with a large collection of trees, that came from the posterior distribution of the model given our data. Then if we want to work with a single tree - that is, to have a point estimate of this posterior distribution of trees - the most usual ways are to calculate the consensus tree or to select the most frequent tree. There are other ways, but let's fix on those by now.

We might not be aware of it, but when we choose for one or another summary we are in fact deciding for the tree estimate that minimizes its distance to all other trees in the set, and in expectation this will be the closest to the true tree under this distance metric (the so called Bayes estimator). This depends on what exactly do we mean by "distance" between trees, and that's what the article "Bayes Estimators for Phylogenetic Reconstruction" (doi 10.1093/sysbio/syr021) is about.

For example, the majority-rule consensus tree is the best we can get if we assume that the Robinson-Foulds distance (RF distance) is a good way of penalizing trees far away from the true one (I won't dwell into the meaning of "truth"; for us, the True tree^® is the one that originated the data). To be more explicit, the consensus tree is the one whose RF distance to all trees in the sample is the shortest possible. This will be the closest we can get to the true tree for this sample, if by "close" we mean "with a small RF distance".

Now suppose I don't like the RF metric because I can only count to two: if the trees are the same the distance is zero, but if they are different then the distance is not zero, and I don't care how different they are (think of apples and oranges). In this case the best representative of my sample is the one that appears more often, known as modal value or Maximum A Posteriori (MAP) value, since our sample comes from a posterior distribution. Is it the closest I can get to the true tree for this distribution? Yes, for this particular definition of distance: the MAP tree is the tree that maximizes the expected coincidence with the true tree.

In the article they also mention that if you want to find the tree that minimizes the expected quartet distance to the true value, then the quartet puzzling method will find this tree for you. But the quartet puzzling tree is not as easy to calculate as the consensus or MAP tree, and there is no straightforward way to find the tree that minimizes other distances in general (e.g. the dSPR, the geodesic distance or the Gene Tree Parsimony). Therefore the authors offer the well-known hill-climbing heuristics for finding the best tree, and use the squared path difference as an example of distance.

Below you can find the presentation I gave to my group last week about this paper, it contains basically some background information and a summary of their method. One thing that is absent from the slides are the results, which I briefly summarize below:

• their method (called "Bayes" in the figures or "BE") always used the path difference as distance measure; this is the overall distance they were trying to minimize.
• they simulated many data sets with several levels of sequence divergence, and reconstructed the phylogeny using Maximum Likelihood, Neighbor-Joining, and Bayesian analysis. From the Bayesian posterior distribution they elected as point estimates the consensus tree, MAP tree, and used their method to find the BE under the path difference.
• Figures 3 and 5 show the distance between the inferred and the true trees, where on figure 3 this distance is the path difference and in figure 5 it is the RF distance. As expected,  the Bayes estimator is better than any other measure at minimizing the path difference distance to the true tree, while the consensus tree wins if we want the closest in terms of RF distance.
• this result is rephrased in figures 8 and 9, which now look specifically at the distances between BE or MAP trees and the true tree. What they plot is distance(BE, true) - distance(MAP, true) for a different definition of distance(,) in each case. The MAP tree is correlated to the consensus tree (if the MAP frequency is larger than 50% they are equal, for instance). Therefore it should come as no surprise that if we define closeness to the true tree in terms of RF distance, the MAP tree will be closer than the BE as shown in figure 9. Because BE assumes that closeness to true is calculated in terms of the path difference, which is reinforced in figure 8.
• The authors wisely avoid offering the "best" Bayes estimator, since it depends on your judgment of how to penalize trees different from the true one.

Journal Club @ UVigo 2011.07.22

View more presentations from Leonardo de Oliveira Martins

OBS: This was my first time using beamer for Latex (after all these years, I know), so the slides are not prime time material. This is also my first submission to slideshare, and I like the idea of an embedded presentation within the blog post. I use latex a lot, and I think it would be easier for me to prepare a post with figures, equations and text within a presentation, and then simply embed it here with a minimum of extra text.

Maybe I'll try this next time, a presentation but with much more text than the recommended - in real life presentations the slides should support and complement but not replace the lecturer. Then you tell me if you would like to read on such a format or if you prefer a more traditional article-ish post.

Reference:
Huggins, P., Li, W., Haws, D., Friedrich, T., Liu, J., & Yoshida, R. (2011). Bayes Estimators for Phylogenetic Reconstruction Systematic Biology, 60 (4), 528-540 DOI: 10.1093/sysbio/syr021