---
title: "Interpreting large quartet distances"
author: "Martin R. Smith"
date: "`r Sys.Date()`"
output: 
  bookdown::pdf_document2:
    includes:
      in_header: ../inst/preamble.tex
  rmarkdown::html_document:
bibliography: ../inst/REFERENCES.bib
csl: https://raw.githubusercontent.com/citation-style-language/styles/master/apa-old-doi-prefix.csl

vignette: >
  %\VignetteIndexEntry{Interpreting large quartet distances}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---
```{r Load package, include=FALSE}
require('Quartet')
require('ape') # for plot.phylo
```

Properties of the quartet distance are explored fully in Steel [-@Steel1993].  
As quartet distances of 1 can only be accomplished for small trees (five or fewer
leaves), it is perhaps more appropriate to consider whether or not trees are more
dissimilar than a pair of random trees, whose distance will be, on average, 
$\frac{2}{3}$.

## Minimum quartet distance

When there are six or more tips in a bifurcating tree, some quartets are necessarily shared between trees.

Consider the tree:
```{r example-tree}
tree_a <- ape::read.tree(text="((1, 2), (3, (4, 5)));")
```
```{r plot-trees, fig.height=1.6, fig.width=2, echo=FALSE}
par(mar=rep(0.3, 4))
plot(tree_a)
```

The only trees with no quartets in common with Tree A are symmetric with  

```{r none-in-common}
tree_b <- ape::read.tree(text="((1, 5), (3, (2, 4)));")
```
```{r plot-tree-b, fig.height=1.6, fig.width=2, echo=FALSE}
par(mar=rep(0.3, 4))
plot(tree_b)
```

Now create Tree C by adding a 6<sup>th</sup> tip as a sister to tip `3` on Tree A.

```{r tree-c}
tree_c <- ape::read.tree(text="((1, 2), ((3, 6), (4, 5)));")
```
```{r Add-tip-6-to-Tree-C, fig.height=1.6, fig.width=2, echo=FALSE}
par(mar=rep(0.3, 4))
plot(tree_c, tip.color=c(1,1,1,2,1,1))
```

There's nowhere to add tip `6` to Tree B without creating a quartet that 
exists in Tree C.

```{r Adding-tip-6-to-Tree-B-duplicates-a-quartet, fig.height=5, fig.width=2.5, echo=FALSE}
PlotApeTree <- function (text, quartet) {
  orig <- TreeTools::RenumberTips(tree_c, as.character(1:6))
  tree <- ape::read.tree(text=text)
  PlotQuartet(list(orig, TreeTools::RenumberTips(tree, as.character(1:6))), quartet, overwritePar=FALSE, cex=0.9)
}

par(mfrow=c(7, 2), mar=rep(0.4, 4), cex=0.9)
PlotApeTree("(((1, 6), 5), (3, (2, 4)));", c(1, 6, 4, 5))
PlotApeTree("((1, 5), (3, ((2, 6), 4)));", c(2, 6, 4, 5))
PlotApeTree("((1, 5), ((3, 6), (2, 4)));", c(3, 6, 4, 5))
PlotApeTree("((1, 5), (3, (2, (4, 6))));", c(4, 6, 1, 2))
PlotApeTree("((1, (5, 6)), (3, (2, 4)));", c(5, 6, 1, 2))
PlotApeTree("(((1, 5), 6), (3, (2, 4)));", c(1, 5, 3, 6))
PlotApeTree("((1, 5), (3, ((2, 4), 6)));", c(4, 2, 3, 6))
```


## Quartet distance in a pair of random trees

On average, $\frac{1}{3}$ of the quartets resolved in a pair of random trees will match. This is because there are
three quartets involving any set of four tips, each of which is equally likely to 
occur on a truly random tree.

The below code calculates the mean proportion of matching quartets for random 
trees with 4 to 20 tips, and the corresponding standard deviation.

```{r mean-proportion}
round(vapply(4:20, function (n_tip) {
 trees <- lapply(logical(56), function (X) 
   ape::rtree(n_tip, tip.label=seq_len(n_tip), br=NULL))
 results <- QuartetStatus(trees)[1, ] / choose(n_tip, 4)
 c(mean(results[-1]), sd(results[-1]))
}, double(2)), 3)
```

## References
