library(lifelihood)
library(tidyverse)
#> ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
#> ✔ dplyr 1.1.4 ✔ readr 2.1.5
#> ✔ forcats 1.0.0 ✔ stringr 1.5.1
#> ✔ ggplot2 3.5.2 ✔ tibble 3.3.0
#> ✔ lubridate 1.9.4 ✔ tidyr 1.3.1
#> ✔ purrr 1.0.4
#> ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
#> ✖ dplyr::filter() masks stats::filter()
#> ✖ dplyr::lag() masks stats::lag()
#> ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
If you haven’t already check it, have a look at:
Load libraries
Data preparation
Load the dataset from .csv file:
Prepare input parameters for the lifelihoodData() function:
# name of the columns of the clutchs into a single vector
clutchs <- c(
"clutch_start1",
"clutch_end1",
"clutch_size1",
"clutch_start2",
"clutch_end2",
"clutch_size2"
)Note: If you have a large number of clutches, it is easier to generate this vector programmatically, particularly if your dataset contains a large number of clutches. See the Generate clutch names vignette.
Create the lifelihoodData object
lifelihoodData() creates a lifelihoodData object, which is a list containing all the information needed to run the lifelihood program of a given dataset of individual life history.
This function mostly takes as input your dataset, your column names.
It also has the model_specs argument, which is a vector of characters with the name of the statistical law to use. Must be of length 3 and each element must be one of "wei" (Weibull law), "exp" (Exponential law), "gam" (Gamma law) or "lgn" (Log-normal law). The first one is used for maturity, the second one is used for clutchs and the third one for death.
dataLFH <- lifelihoodData(
df = df,
sex = "sex",
sex_start = "sex_start",
sex_end = "sex_end",
maturity_start = "mat_start",
maturity_end = "mat_end",
clutchs = clutchs,
death_start = "death_start",
death_end = "death_end",
covariates = c("geno", "type"),
model_specs = c("gam", "lgn", "wei")
)Get the results
All default parameters
Once you have created your lifelihoodData object with lifelihoodData(), you can call the lifelihood() function to run the lifelihood program.
It returns a lifelihoodResults object, which is a list containing all the results of the analysis.
Here it’s a minimalist usage of the function, where we only specify the lifelihoodData object, the path to the configuration file and the seeds to use. The raise_estimation_warning argument will be the focus of the next vignette.
results <- lifelihood(
lifelihoodData = dataLFH,
path_config = get_config_path("config"),
seeds = c(1, 2, 3, 4),
raise_estimation_warning = FALSE
)
summary(results)
#> LIFELIHOOD RESULTS
#>
#> Likelihood:
#> [1] -829.0578
#>
#> Effects:
#> name estimation stderror
#> 1 int_expt_death -14.95337701 0
#> 2 eff_expt_death_geno_1 2.27064527 0
#> 3 eff_expt_death_type_1 -9.55371757 0
#> 4 eff_expt_death_type_2 0.53315647 0
#> 5 int_survival_shape 10.28792377 0
#> 6 eff_survival_shape_geno_1 -1.93963186 0
#> 7 eff_survival_shape_type_1 7.57026237 0
#> 8 eff_survival_shape_type_2 4.74998927 0
#> 9 eff_survival_shape_type_1:geno_1 -0.24718760 0
#> 10 eff_survival_shape_type_2:geno_1 3.52427897 0
#> 11 int_ratio_expt_death -1.45299283 0
#> 12 eff_ratio_expt_death_geno_1 3.02680049 0
#> 13 int_prob_death -17.15815657 0
#> 14 eff_prob_death_geno_1 0.02066385 0
#> 15 int_expt_maturity -0.07781130 0
#> 16 eff_expt_maturity_geno_1 -0.34729932 0
#> 17 eff_expt_maturity_type_1 -0.76380036 0
#> 18 eff_expt_maturity_type_2 -0.75476253 0
#> 19 int_maturity_shape -10.02261060 0
#> 20 eff_maturity_shape_geno_1 -3.30865726 0
#> 21 int_ratio_expt_maturity 1.74264178 0
#> 22 eff_ratio_expt_maturity_geno_1 3.36233997 0
#> 23 eff_ratio_expt_maturity_type_1 -2.79189482 0
#> 24 eff_ratio_expt_maturity_type_2 -4.87384990 0
#> 25 int_expt_reproduction -1.79104715 0
#> 26 eff_expt_reproduction_geno_1 -0.89759305 0
#> 27 int_pontn -2.96400432 0
#> 28 eff_pontn_geno_1 -1.03064224 0
#> 29 int_increase_death_hazard 3.32326622 0
#> 30 eff_increase_death_hazard_geno_1 2.18470307 0
#> 31 eff_increase_death_hazard_type_1 -0.38970742 0
#> 32 eff_increase_death_hazard_type_2 -3.84713473 0
#> 33 eff_increase_death_hazard_type_1:geno_1 -1.13985166 0
#> 34 eff_increase_death_hazard_type_2:geno_1 -3.69784224 0
#> 35 int_tof_reduction_date -1.27913674 0
#> 36 eff_tof_reduction_date_geno_1 -3.51198283 0
#> 37 int_increase_tof_n_offspring 0.24026695 0
#> 38 eff_increase_tof_n_offspring_geno_1 -0.32428130 0
#> 39 int_lin_decrease_hazard 0.81074061 0
#> 40 eff_lin_decrease_hazard_geno_1 0.19505972 0
#> 41 eff_lin_decrease_hazard_type_1 5.50927699 0
#> 42 eff_lin_decrease_hazard_type_2 -0.85842716 0
#> 43 eff_lin_decrease_hazard_type_1:geno_1 3.98397092 0
#> 44 eff_lin_decrease_hazard_type_2:geno_1 0.69812447 0
#> 45 int_quad_senescence -0.49711337 0
#> 46 eff_quad_senescence_geno_1 0.47173154 0
#> 47 int_quad_decrease_hazard -0.07370047 0
#> 48 eff_quad_decrease_hazard_geno_1 -0.14027753 0
#> 49 int_quad_change_n_offspring 0.07862816 0
#> 50 eff_quad_change_n_offspring_geno_1 0.00640363 0
#> 51 int_tof_n_offspring -0.21320110 0
#> parameter kind event
#> 1 expt_death intercept mortality
#> 2 expt_death coefficient mortality
#> 3 expt_death coefficient mortality
#> 4 expt_death coefficient mortality
#> 5 survival_shape intercept mortality
#> 6 survival_shape coefficient mortality
#> 7 survival_shape coefficient mortality
#> 8 survival_shape coefficient mortality
#> 9 survival_shape coefficient mortality
#> 10 survival_shape coefficient mortality
#> 11 ratio_expt_death intercept mortality
#> 12 ratio_expt_death coefficient mortality
#> 13 prob_death intercept mortality
#> 14 prob_death coefficient mortality
#> 15 expt_maturity intercept maturity
#> 16 expt_maturity coefficient maturity
#> 17 expt_maturity coefficient maturity
#> 18 expt_maturity coefficient maturity
#> 19 maturity_shape intercept maturity
#> 20 maturity_shape coefficient maturity
#> 21 ratio_expt_maturity intercept maturity
#> 22 ratio_expt_maturity coefficient maturity
#> 23 ratio_expt_maturity coefficient maturity
#> 24 ratio_expt_maturity coefficient maturity
#> 25 expt_reproduction intercept reproduction
#> 26 expt_reproduction coefficient reproduction
#> 27 pontn intercept reproduction
#> 28 pontn coefficient reproduction
#> 29 increase_death_hazard intercept reproduction
#> 30 increase_death_hazard coefficient reproduction
#> 31 increase_death_hazard coefficient reproduction
#> 32 increase_death_hazard coefficient reproduction
#> 33 increase_death_hazard coefficient reproduction
#> 34 increase_death_hazard coefficient reproduction
#> 35 tof_reduction_date intercept reproduction
#> 36 tof_reduction_date coefficient reproduction
#> 37 increase_tof_n_offspring intercept reproduction
#> 38 increase_tof_n_offspring coefficient reproduction
#> 39 lin_decrease_hazard intercept reproduction
#> 40 lin_decrease_hazard coefficient reproduction
#> 41 lin_decrease_hazard coefficient reproduction
#> 42 lin_decrease_hazard coefficient reproduction
#> 43 lin_decrease_hazard coefficient reproduction
#> 44 lin_decrease_hazard coefficient reproduction
#> 45 quad_senescence intercept reproduction
#> 46 quad_senescence coefficient reproduction
#> 47 quad_decrease_hazard intercept reproduction
#> 48 quad_decrease_hazard coefficient reproduction
#> 49 quad_change_n_offspring intercept reproduction
#> 50 quad_change_n_offspring coefficient reproduction
#> 51 tof_n_offspring intercept reproductionGet specific results
The lifelihoodResults object is a list containing all the results of the analysis. We can get specific results by calling the list element.
Get estimations:
results$effectsGet maximum likelihood found:
results$likelihood
#> [1] -829.0578Next step
Now that you have seen how to use the package, you can go further and customise your parameter boundaries and deal with estimation warnings.