# Likelihood of a dataset under models with biogeography fit to a subgroup. Computes the likelihood of a dataset under either the linear or exponential diversity dependent model with specified `sigma2` and slope values and with a `geography.object` formed using `CreateGeoObject`. ```r likelihood_subgroup_model(data,phylo,geography.object,model=c("MC","DDexp","DDlin"), par,return.z0=FALSE,maxN=NULL,error=NULL) ``` ## Arguments - `phylo`: an object of type 'phylo' (see ape documentation) produced as "map" from CreateGeobyClassObject. NB: the length of this object need not match number of items in data, since map may include tips outside of group with some part of their branch in the group - `data`: a named vector of continuous data for a subgroup of interest with names corresponding to `phylo$tip.label` - `geography.object`: a list of sympatry/group membership through time created using `CreateGeobyClassObject` - `model`: model chosen to fit trait data, `"DDlin"` is the diversity-dependent linear model, and `"DDexp"` is the diversity-dependent exponential model of Weir & Mursleen 2013. - `par`: a vector listing a value for `log(sig2)` (see Note) and either `b` (for the linear diversity dependent model) or `r` (for the exponential diversity dependent model), in that order. - `return.z0`: logical indicating whether to return an estimate of the trait value at the root given the parameter values (if `TRUE`, function returns root value rather than negative log-likelihood) - `maxN`: when fitting `DDlin` model, it is necessary to specify the maximum number of sympatric lineages to ensure that the rate returned does not correspond to negative sig2 values at any point in time (see Details). - `error`: A named vector with standard errors (SE) of trait values for each species (with names matching `"phylo$tip.label"`). The default is NULL, in this case potential error is ignored in the fit. If set to NA, the SE is estimated from the data (to be used when there are no error measurements, a nuisance parameter is estimated). Note: When standard errors are provided, a nuisance parameter is also estimated. ## Details When specifying `par`, `log(sig2)` (see Note) must be listed before the slope parameter (`b` or `r`). maxN can be calculated using `maxN=max(vapply(geo.object$geography.object,function(x)max(rowSums(x)),1))`, where geo.object is the output of `CreateGeoObject` ## Note To stabilize optimization, this function exponentiates the input `sig2` value, thus the user must input the log(sig2) value to compute the correct log likelihood (see example). ## Returns The negative log-likelihood value of the dataset (accordingly, the negative of the output should be recorded as the likelihood), given the phylogeny, `sig2` and slope values, and `geography.object`. If `return.z0=TRUE`, the estimated root value for the par values is returned instead of the negative log-likelihood. ## References Drury, J., Clavel, J., Manceau, M., and Morlon, H. 2016. Estimating the effect of competition on trait evolution using maximum likelihood inference. **Systematic Biology** doi 10.1093/sysbio/syw020 Weir, J. & Mursleen, S. 2012. Diversity-dependent cladogenesis and trait evolution in the adaptive radiation of the auks (Aves: Alcidae). **Evolution** 67:403-416. ## Author(s) Jonathan Drury [jonathan.p.drury@gmail.com](mailto:jonathan.p.drury@gmail.com) Julien Clavel ## See Also `fit_t_comp` `CreateGeoObject` `likelihood_t_DD` ## Examples ```r data(BGB.examples) Canidae.phylo<-BGB.examples$Canidae.phylo dummy.group<-c(rep("B",3),rep("A",12),rep("B",2),rep("A",6),rep("B",5),rep("A",6)) names(dummy.group)<-Canidae.phylo$tip.label Canidae.simmap<-phytools::make.simmap(Canidae.phylo, dummy.group) set.seed(123) Canidae.data<-rnorm(length(Canidae.phylo$tip.label)) names(Canidae.data)<-Canidae.phylo$tip.label Canidae.A<-Canidae.data[which(dummy.group=="A")] Canidae.geobyclass.object<-CreateGeobyClassObject(phylo=Canidae.phylo, simmap=Canidae.simmap, trim.class="A", ana.events=BGB.examples$Canidae.ana.events, clado.events=BGB.examples$Canidae.clado.events,stratified=FALSE, rnd=5) par <- c(log(0.01),-0.000005) maxN<-max(vapply(Canidae.geobyclass.object$geo.object$geography.object, function(x)max(rowSums(x)),1)) lh <- -likelihood_subgroup_model(data=Canidae.A, phylo=Canidae.geobyclass.object$map, geography.object=Canidae.geobyclass.object$geo.object, model="DDlin", par=par, return.z0=FALSE, maxN=maxN) ```

likelihood_subgroup_model function