Description Usage Arguments Details Value Author(s) References See Also Examples

This function estimates the population parameter *THETA*
from a genealogy (coded a as phylogenetic tree) under the coalescent.

1 | ```
theta.tree(phy, theta, fixed = FALSE, analytical = TRUE, log = TRUE)
``` |

`phy` |
an object of class |

`theta` |
a numeric vector. |

`fixed` |
a logical specifying whether to estimate |

`analytical` |
a logical specifying whether to use analytical
formulae to estimate |

.

`log` |
a logical specifying whether to return the likelihoods on a
log scale (the default); ignored if |

The tree `phy`

is considered as a genealogy, and therefore should
be ultrametric. By default, *THETA* is estimated by
maximum likelihood and the value given in `theta`

is used as
starting value for the minimisation function (if several values are
given as a vector the first one is used). If `fixed = TRUE`

,
then the [log-]likelihood values are returned corresponding to each
value in `theta`

.

The present implementation does a numerical optimisation of the
log-likelihood function (with `nlminb`

) with the
first partial derivative as gradient. It is possible to solve the
latter and have a direct analytical MLE of *THETA* (and
its standard-error), but this does not seem to be faster.

If `fixed = FALSE`

, a list with two elements:

`theta` |
the maximum likelihood estimate of |

`logLik` |
the log-likelihood at its maximum. |

If `fixed = TRUE`

, a numeric vector with the [log-]likelihood
values.

Emmanuel Paradis

Kingman, J. F. C. (1982) The coalescent. *Stochastic Processes
and their Applications*, **13**, 235–248.

Kingman, J. F. C. (1982) On the genealogy of large
populations. *Journal of Applied Probability*, **19A**,
27–43.

Wakeley, J. (2009) *Coalescent Theory: An Introduction.*
Greenwood Village, CO: Roberts and Company Publishers.

1 2 3 4 5 6 7 | ```
tr <- rcoal(50) # assumes theta = 1
theta.tree(tr, 10)
theta.tree(tr, 10, analytical = FALSE) # uses nlminb()
## profile log-likelihood:
THETA <- seq(0.5, 1.5, 0.01)
logLikelihood <- theta.tree(tr, THETA, fixed = TRUE)
plot(THETA, logLikelihood, type = "l")
``` |

pegas documentation built on May 29, 2017, 6:33 p.m.

Embedding an R snippet on your website

Add the following code to your website.

For more information on customizing the embed code, read Embedding Snippets.