Likelihood ratio test (LRT)

Likelihood ratio test (LRT)
• Always compares an unconstrained to a constrained
model
• Constrained model must be nested within the
unconstrained model
• Parameter(s) take on their maximum likelihood
estimates (MLEs) in the unconstrained model
• At least some parameters that are estimated in the
unconstrained model are set to some particular value
of interest in the constrained model
• Unconstrained model must therefore equal or exceed
the constrained model in its fit to the data (as measured
by the maximized likelihood)
Paul O. Lewis ~ Phylogenetics, Spring 2014
Likelihood Ratio Test Statistic (LRT)
is the MLE
(possibly multidimensional)
is some other value
Coin-flipping example:
Data:
Constrained model:
Unconstrained model:
6 heads out of 10 flips
fair coin (θ = 0.5)
biased coin (θ = )
Example of likelihood calculation for case of θ = 0.6
Paul O. Lewis ~ Phylogenetics, Spring 2014
Likelihood Ratio Test
Coin-flipping example:
Data:
Constrained model:
Unconstrained model:
6 heads out of 10 flips
fair coin (θ = 0.5)
biased coin (θ = )
L(0.6)
0.251
LRT = 2 log
= 2 log
= 0.404886
L(0.5)
0.205
Not significant: P = 0.527
This means that the
simpler, constrained
model cannot be
rejected
LRT approximates a chi-square random variable with d.f. equal to the
difference in the number of free parameters between the two models
Paul O. Lewis ~ Phylogenetics, Spring 2014
Likelihood Ratio Test
Chi-squared
Test relies on fact that LR
test statistic approximates
the chi-squared test statistic
55 heads/100 flips
Paul O. Lewis ~ Phylogenetics, Spring 2014
LRT
Examples of unconstrained vs.
constrained model comparisons
Model
Constrained
Unconstrained
GTR+G
shape = ∞
shape = MLE
K80
κ = 1.0
κ = MLE
HKY+I+G
p
p
HKY+I+G
p
shape = ∞
p
shape = MLE
GY94 (codon)
ω = 1.0
ω = MLE
Paul O. Lewis ~ Phylogenetics, Spring 2014
Problems at the border
0.4
Cases in which the constrained model involves setting a parameter to the
edge of its valid range require special consideration (see Ota et al. 2000)
0.0
0.1
0.2
0.3
The problem is that the theory supporting LRT ≈ Χ2
depends on the assumption that the likelihood is
approximately normal, which isn’t true when a
parameter value is on the border
-4
-2
0
2
4
Ignoring this causes the test to be conservative
(i.e., the simpler model is rejected less often
than it should be)
Note: Modeltest implements this correction
Ota, R., P. J. Waddell, M. Hasegawa, H. Shimodaira, and H. Kishino. 2000.
Appropriate likelihood ratio tests and marginal distributions for evolutionary
tree models with constraints on parameters. Molecular Biology and Evolution
17:798-803.
Paul O. Lewis ~ Phylogenetics, Spring 2014
In fact it is half-normal
and half zeros
Testing the molecular clock
Unconstrained model: need to estimate 2n-3 = 11 branch lengths
Constrained model: need to estimate n-1 = 6 divergence times
1
2
3 4
5 6
7
t1
t2
t3
t4
t5
t6
Likelihood
ratio test thus has
(2n-3) - (n-1) = n-2 d.f.
n = 7 taxa
Felsenstein, J. 1983. Statistical inference of phylogenies. Journal of the Royal Statistical Society A 146:246-272.
Paul O. Lewis ~ Phylogenetics, Spring 2014
Akaike Information Criterion
•
•
•
•
•
AIC = -2 max(lnL) + 2K
K is number of free model parameters
Measures expected relative distance to true model
Model with smallest AIC wins
Advantage over LRT: non-nested models
Example: 6 heads/10 flips revisited
Unconstrained model:θ = 0.6, AIC = -2(-1.383) + 2(1) = 4.766
Constrained model: θ = 0.5, AIC = -2(-1.584) + 2(0) = 3.168 (best)
Paul O. Lewis ~ Phylogenetics, Spring 2014
Akaike, H. 1973. Information theory as an extension of the maximum likelihood
principle. Pages 267-281 in B. N. Petrov and F. Csaki (eds.), Second International
Symposium on Information Theory. Akademiai Kiado, Budapest.
Akaike Information Criterion
(AIC)
true
Calculate AIC for each model:
AIC = 2k
2 log(Lmax )
AICfree = 2(3)
2( 43.1) = 92.2
AICequal = 2(0)
2( 44.4) = 88.8
The constrained model ("equal") is
a better choice than the
unconstrained model ("free")
according to AIC
92.2 = twice
expected (relative)
K-L divergence
from free model to
true model
88.8 = twice
expected (relative)
K-L divergence
from equal model
to true model
(K-L stands for
Kullback-Leibler)
equal
Paul O. Lewis ~ Phylogenetics, Spring 2014
Bayesian Information Criterion
•
•
•
•
•
•
BIC = -2 max(lnL) + K log(n)
K is number of free model parameters
n is the sample size
Model with smallest BIC wins
Advantage over LRT: non-nested models
Considered superior to both AIC and LRT
Example: 6 heads/10 flips one more time. Note: log(10) ≈ 2.3
Unconstrained model:θ = 0.6, BIC = -2(-1.383) + (2.3)(1) = 5.066
Constrained model: θ = 0.5, BIC = -2(-1.584) + 0 = 3.168 (best)
(We will discuss BIC more after an introduction to Bayesian statistics.)
Paul O. Lewis ~ Phylogenetics, Spring 2014
Decision-Theoretic Approach
• Rationale: “If a simple model is returning estimates of branch
lengths that are nearly identical to those from a more complex
model, there will be little difference in phylogenetic
estimation under the two models.” — Minin et al. (2003:676)
• Therefore, you might as well choose the simpler model.
• Suppose there are 3 models in contention: Model 1: GTR Model 2: HKY85 Model 3: JC
• Suppose further that we have a measure, P(M|D), of the
probability that a model M is the true model given the data D
Paul O. Lewis ~ Phylogenetics, Spring 2014
Minin, V., Z. Abdo, P. Joyce, and J. Sullivan. 2003. Performance-based selection of
likelihood models for phylogeny estimation. Systematic Biology 52:674–683.
Decision-Theoretic Approach
Risk (R) associated with each model:
GTR:
R1 = d12 P (M2 |D) + d13 P (M3 |D)
HKY:
R2 = d12 P (M1 |D) + d23 P (M3 |D)
JC:
R3 = d13 P (M1 |D) + d23 P (M2 |D)
R1 = 0.30
1.8
0.5
GTR
0.3
0.4
HKY
R2 = 0.26
best
Paul O. Lewis ~ Phylogenetics, Spring 2014
0.1
JC
R3 = 1.34
worst
1.1
HKY wins according to DT
even though it does not fit
the data as well as GTR
Automating model
selection using
JModelTest
• Download JModelTest from http://darwin.uvigo.es/
!
• Unzip it, and double-click the JModelTest.jar file
!
• Use File > Load DNA Alignment to read in a data
file (note: it may choke if your Nexus-formatted
data file is complicated)
!
• Use Analysis > Compute likelihood scores to
compute maximized log-likelihoods for each
model (JModelTest uses a built-in version of
PhyML for this step)
!
• Use Analysis > Do AIC calculations... to do
model comparison using AIC (can also choose
BIC, DT, or hLRT to use these criteria)
JModelTest: Posada D. 2008. jModelTest: Phylogenetic Model Averaging. Molecular Biology and Evolution 25: 1253-1256.
PhyML: Guindon S, Gascuel O. 2003. A simple, fast, and accurate algorithm to estimate large phylogenies by maximum likelihood. Systematic Biology 52:
696-704.
Paul O. Lewis ~ Phylogenetics, Spring 2014

Download Report