Fit two models to two separate segments to an evolutionary sequence (time-series).
Source:R/fit.mode.shift.R
fit.mode.shift.RdWrapper function to find maximum likelihood solutions to two models to an evolutionary sequence.
Arguments
- y
an univariate evoTS object.
- model1
the model fitted to the first segment. Options are Stasis, URW, GRW, OU.
- model2
the model fitted to the second segment. Options are Stasis, URW, GRW, OU.
- fit.all
logical indicating whether to fit all pairwise combinations of the four models to the evolutionary sequence (time-series).
- minb
the minimum number of samples within a segment to consider
- shift.point
The sample that split the time-series into two segments. The samples are passed to the argument as a vector. Default is NULL, which means all possible shift points will be assessed constrained by how minb is defined.
- pool
logical indicating whether to pool variances across samples
- silent
if TRUE, less information is printed to the screen as the model is fit
- hess
logical, indicating whether to calculate standard errors from the Hessian matrix.
#'
Value
the function returns a list of all investigated models and their highest log-likelihood (and their corresponding AICc and AICc weight).
- logL
the log-likelihood of the optimal solution
- AICc
AIC with a correction for small sample sizes
- parameters
parameter estimates
- modelName
abbreviated model name
- method
Joint consideration of all samples
- K
number of parameters in the model
- n
the number of observations/samples
- all.logl
log-likelihoods for all tested partitions of the series into segments. Will return a single value if shift points have been given
- GG
matrix of indices of initial samples of each tested segment configuration; each column of GG corresponds to the elements of all.logl
In addition, if fit.all=TRUE the function also returns a list of all investigated models and their highest log-likelihood (and their corresponding AICc and AICc weight).
Note
The models have been implemented to be compatible with the joint parameterization routine in the package paleoTS. The optimization is therefore fit using the actual sample values, with the autocorrelation among samples accounted for in the log-likelihood function. The joint distribution of sample means is multivariate normal, with means and variance-covariances determined by evolutionary parameters and sampling errors.
References
Hunt, G. 2006. Fitting and comparing models of phyletic evolution: random walks and beyond. Paleobiology 32:578–601
Hunt, G., Bell, M. A. & Travis, M. P. Evolution towards a new adaptive optimum: Phenotypic evolution in a fossil stickleback lineage. Evolution 62:700–710 (2008)
Examples
##Generate a paleoTS object.
x <- paleoTS::sim.GRW(30)
## Fit a mode-shift model without defining a shift point (the example may take > 5 seconds to run)
fit.mode.shift(x, model1="URW", model2="Stasis")
#> [1] "Searching all possible shift points in the evolutionary sequence"
#> Total # hypotheses: 17
#> 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17
#>
#> paleoTSfit object [n = 30 , K = 5 ]
#>
#> Model: URW-Stasis
#> Method: Joint
#> log-likelihood = -22.07077
#> AICc = 56.64154
#>
#> Parameter estimates:
#> anc vstep theta omega shift1
#> 0.24968675 0.07526602 1.89963879 2.70278953 24.00000000
#>
#> Log-likelihoods of all tested shift-points
#>
#> shift all.logl
#> 8 -32.70246
#> 9 -31.68724
#> 10 -31.07975
#> 11 -29.54590
#> 12 -28.64491
#> 13 -30.72209
#> 14 -32.03039
#> 15 -30.62360
#> 16 -29.62931
#> 17 -28.58026
#> 18 -27.68512
#> 19 -27.36638 *
#> 20 -26.65537 *
#> 21 -25.48556 *
#> 22 -24.37638 *
#> 23 -23.11386 *
#> 24 -22.07077 **
#>
#> Shift occurs immediately AFTER listed sample number
#> ** = maximum-likelihood shift point
#> * = additional shift points in CI [ within 5.535 logL units; Chi-sq P = 0.95 , df = 5 ]
#>
#>
#> Additional elements not printed: convergence logLFunction