Merge pull request #135 from rsetienne/develop

v4.0.2
rsetienne · May 26, 2021 · e607ce8 · e607ce8
2 parents 6da0e3f + 094495b
commit e607ce8
Show file tree

Hide file tree

Showing 69 changed files with 1,252 additions and 1,822 deletions.
diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml
@@ -6,6 +6,8 @@ on:
       - main
       - master
       - develop
+      - DAISIE_IW_B
+
   pull_request:
     branches:
       - main
@@ -84,13 +86,6 @@ jobs:
         run: find check -name 'testthat.Rout*' -exec cat '{}' \; || true
         shell: bash
 
-      - name: Install covr and test coverage
-        if: matrix.config.r == 'devel'
-        run: |
-          remotes::install_cran("covr")
-          covr::codecov()
-        shell: Rscript {0}
-
       - name: Upload check results
         if: failure()
         uses: actions/upload-artifact@main

diff --git a/.github/workflows/test-coverage.yaml b/.github/workflows/test-coverage.yaml
@@ -0,0 +1,51 @@
+on:
+  push:
+    branches:
+      - main
+      - master
+      - develop
+  pull_request:
+    branches:
+      - main
+      - master
+      - develop
+
+name: test-coverage
+
+jobs:
+  test-coverage:
+    runs-on: ubuntu-20.04
+
+    env:
+      GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }}
+    steps:
+      - uses: actions/checkout@v2
+
+      - uses: r-lib/actions/setup-r@v1
+
+      - uses: r-lib/actions/setup-pandoc@v1
+
+      - name: Query dependencies
+        run: |
+          install.packages('remotes')
+          saveRDS(remotes::dev_package_deps(dependencies = TRUE), ".github/depends.Rds", version = 2)
+          writeLines(sprintf("R-%i.%i", getRversion()$major, getRversion()$minor), ".github/R-version")
+        shell: Rscript {0}
+
+      - name: Restore R package cache
+        uses: actions/cache@v2
+        with:
+          path: ${{ env.R_LIBS_USER }}
+          key: ${{ runner.os }}-${{ hashFiles('.github/R-version') }}-1-${{ hashFiles('.github/depends.Rds') }}
+          restore-keys: ${{ runner.os }}-${{ hashFiles('.github/R-version') }}-1-
+
+      - name: Install dependencies
+        run: |
+          install.packages(c("remotes"))
+          remotes::install_deps(dependencies = TRUE)
+          remotes::install_cran("covr")
+        shell: Rscript {0}
+
+      - name: Test coverage
+        run: covr::codecov()
+        shell: Rscript {0}
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,8 +1,8 @@
 Package: DAISIE
 Type: Package
 Title: Dynamical Assembly of Islands by Speciation, Immigration and Extinction
-Version: 3.2.1
-Date: 2021-03-23
+Version: 4.0.2
+Date: 2021-05-26
 Depends: R (>= 3.5.0)
 biocViews:
 SystemRequirements: C++14
@@ -34,7 +34,7 @@ Suggests:
     gridExtra,
     dplyr,
     ggplot2,
-    ggtree,
+    ggtree (>= 3.0.0),
     tidytree,
     tidyr,
     purrr,

diff --git a/NEWS.md b/NEWS.md
@@ -1,14 +1,38 @@
+# DAISIE 4.0.2
+
+**N.B.: MacOS users may experience issues when installing DAISIE, especially when on MacOS Big Sur. If that is you case, please see [here](https://github.com/rsetienne/DAISIE/blob/6da0e3f65680d5f237345ef80935bda7541cf230/doc/DAISIE_macOS.md) for detailed installation instructions.**
+
+* Suggest ggtree >= 3.0.0.
+
+# DAISIE 4.0.1
+
+**N.B.: MacOS users may experience issues when installing DAISIE, especially when on MacOS Big Sur. If that is you case, please see [here](https://github.com/rsetienne/DAISIE/blob/6da0e3f65680d5f237345ef80935bda7541cf230/doc/DAISIE_macOS.md) for detailed installation instructions.**
+
+* Fix possibility of fitting CS model with IW likelihood on simulated data by setting `CS_version = 0`. Improve `CS_version` documentation.
+
+# DAISIE 4.0.0
+
+**N.B.: MacOS users may experience issues when installing DAISIE, especially when on MacOS Big Sur. If that is you case, please see [here](https://github.com/rsetienne/DAISIE/blob/6da0e3f65680d5f237345ef80935bda7541cf230/doc/DAISIE_macOS.md) for detailed installation instructions.**
+
+* Fix bug when calculating conditional probabilities, which are now correctly calculated from island age to the present.
+* Fix bug when calculating probabilities upon migration, which assumed no recolonisation was possible in the CS model. Handling recolonisation in the same manner is not possible for the IW model, so an approximation is now made. The influence of the bug and of the approximation in the IW model is expected to be minimal, particularly in cases where the colonisation rate is low. Such cases of low colonisation are the norm.
+* Add `num_cycles` argument to ML functions, allowing to specify how many cycles the optimizer should take. Defaults to 1.
+* Minor vignette corrections.
+* Improve README.md documentation.
+
 # DAISIE 3.2.1
+[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.4633441.svg)](https://doi.org/10.5281/zenodo.4633441)
 
-**N.B.: MacOS users may experience issues when installing DAISIE, especially when on MacOS Big Sur. If that is you case, please see [here](doc/DAISIE_macOS.md) for detailed installation instructions.**
+**N.B.: MacOS users may experience issues when installing DAISIE, especially when on MacOS Big Sur. If that is you case, please see [here](https://github.com/rsetienne/DAISIE/blob/6da0e3f65680d5f237345ef80935bda7541cf230/doc/DAISIE_macOS.md) for detailed installation instructions.**
 
 * Minor documentation improvements.
 
 # DAISIE 3.2.0
+[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.4630604.svg)](https://doi.org/10.5281/zenodo.4630604)
+
 
-**N.B.: MacOS users may experience issues when installing DAISIE, especially when on MacOS Big Sur. If that is you case, please see [here](doc/DAISIE_macOS.md) for detailed installation instructions.**
+**N.B.: MacOS users may experience issues when installing DAISIE, especially when on MacOS Big Sur. If that is you case, please see [here](https://github.com/rsetienne/DAISIE/blob/6da0e3f65680d5f237345ef80935bda7541cf230/doc/DAISIE_macOS.md) for detailed installation instructions.**
 
-## Changes
 * `DAISIE_loglikg_IW()` is now more efficient and numerically stable. Numerical integration is now done via C++ with package `odeint`.
 * Add relaxed rate capabilities (both inference and simulations). Relaxed rate models allow for parameters to not be static, but to be sampled by specific probability distributions.
 * Introduce `MinAge` data status in DAISIE data objects. A status containing `MinAge` sets a lower boundary for colonization in situations when the precise colonization time is unknown. This is interpreted by `DAISIE_dataprep()` so that the information is passed on to the likelihood optimization functions. See the `DAISIE_dataprep()` help page for more details. In the back-end this results in new `stac` values 8 and 9.

diff --git a/R/DAISIE_ExpEIN.R b/R/DAISIE_ExpEIN.R
@@ -45,7 +45,7 @@ DAISIE_ExpEIN <- function(t, pars, M, initEI = c(0, 0)) {
    DD <- laa + 2 * lac
    E0 <- initEI[1]
    I0 <- initEI[2]
-   if (t == Inf) {
+   if (t[1] == Inf) {
       Imm <- ga * M2 / B
       End <- DD / A * Imm
    } else {

diff --git a/R/DAISIE_ML1.R b/R/DAISIE_ML1.R
@@ -95,7 +95,8 @@ DAISIE_ML1 <- function(
   verbose = 0,
   tolint = c(1E-16, 1E-10),
   island_ontogeny = NA,
-  jitter = 0) {
+  jitter = 0,
+  num_cycles = 1) {
   # datalist = list of all data: branching times, status of clade, and numnber of missing species
   # datalist[[,]][1] = list of branching times (positive, from present to past)
   # - max(brts) = age of the island
@@ -297,7 +298,8 @@ DAISIE_ML1 <- function(
     CS_version = CS_version,
     abstolint = tolint[1],
     reltolint = tolint[2],
-    jitter = jitter
+    jitter = jitter,
+    num_cycles = num_cycles
   )
   if (out$conv != 0) {
     cat(

diff --git a/R/DAISIE_ML2.R b/R/DAISIE_ML2.R
@@ -77,7 +77,8 @@ DAISIE_ML2 <- function(
   optimmethod = "subplex",
   verbose = 0,
   tolint = c(1E-16, 1E-10),
-  jitter = 0) {
+  jitter = 0,
+  num_cycles = 1) {
 # datalist = list of all data: branching times, status of clade, and numnber of missing species
 # datalist[[,]][1] = list of branching times (positive, from present to past)
 # - max(brts) = age of the island
@@ -169,7 +170,8 @@ DAISIE_ML2 <- function(
     methode = methode,
     abstolint = tolint[1],
     reltolint = tolint[2],
-    jitter = jitter
+    jitter = jitter,
+    num_cycles = num_cycles
   )
   if (out$conv != 0) {
     cat("Optimization has not converged. Try again with different initial values.\n")

diff --git a/R/DAISIE_ML3.R b/R/DAISIE_ML3.R
@@ -67,7 +67,8 @@ DAISIE_ML3 <- function(
   CS_version = 1,
   verbose = 0,
   tolint = c(1E-16, 1E-10),
-  jitter = 0) {
+  jitter = 0,
+  num_cycles = 1) {
 # datalist = list of all data: branching times, status of clade, and numnber of missing species
 # datalist[[,]][1] = list of branching times (positive, from present to past)
 # - max(brts) = age of the island
@@ -166,7 +167,8 @@ DAISIE_ML3 <- function(
       CS_version = CS_version,
       abstolint = tolint[1],
       reltolint = tolint[2],
-      jitter = jitter
+      jitter = jitter,
+      num_cycles = num_cycles
     )
   if (out$conv != 0) {
     cat("Optimization has not converged. Try again with different initial values.\n")

diff --git a/R/DAISIE_ML4.R b/R/DAISIE_ML4.R
@@ -64,7 +64,7 @@ DAISIE_loglik_all_choosepar4 <- function(trparsopt,
 #' \item{loglik}{ gives the maximum loglikelihood}
 #' \item{df}{
 #' gives the number
-#' of estimated parameters, i.e. degrees of feedom}
+#' of estimated parameters, i.e. degrees of freedom}
 #' \item{conv}{ gives a
 #' message on convergence of optimization; conv = 0 means convergence}
 #' @keywords internal
@@ -86,7 +86,8 @@ DAISIE_ML4 <- function(
   verbose = 0,
   tolint = c(1E-16, 1E-10),
   island_ontogeny = NA,
-  jitter = 0) {
+  jitter = 0,
+  num_cycles = 1) {
 
   out2err <- data.frame(
     lambda_c = NA,
@@ -194,7 +195,8 @@ DAISIE_ML4 <- function(
     CS_version = CS_version,
     abstolint = tolint[1],
     reltolint = tolint[2],
-    jitter = jitter
+    jitter = jitter,
+    num_cycles = num_cycles
   )
   if (out$conv != 0) {
     cat(

diff --git a/R/DAISIE_ML_CS.R b/R/DAISIE_ML_CS.R
@@ -182,7 +182,8 @@ DAISIE_ML_CS <- DAISIE_ML <- function(
   CS_version = 1,
   verbose = 0,
   tolint = c(1E-16, 1E-10),
-  jitter = 0) {
+  jitter = 0,
+  num_cycles = 1) {
 
   if (datatype == "single") {
     if (is.na(island_ontogeny)) {
@@ -203,7 +204,8 @@ DAISIE_ML_CS <- DAISIE_ML <- function(
                           CS_version = CS_version,
                           verbose = verbose,
                           tolint = tolint,
-                          jitter = jitter)
+                          jitter = jitter,
+                          num_cycles = num_cycles)
       } else
       {
         out <- DAISIE_ML1(datalist = datalist,
@@ -226,7 +228,8 @@ DAISIE_ML_CS <- DAISIE_ML <- function(
                           CS_version = CS_version,
                           verbose = verbose,
                           tolint = tolint,
-                          jitter = jitter)
+                          jitter = jitter,
+                          num_cycles = num_cycles)
       }
     } else
       {
@@ -246,7 +249,8 @@ DAISIE_ML_CS <- DAISIE_ML <- function(
                           CS_version = CS_version,
                           verbose = verbose,
                           tolint = tolint,
-                          jitter = jitter)
+                          jitter = jitter,
+                          num_cycles = num_cycles)
       }
   } else
     {
@@ -265,7 +269,8 @@ DAISIE_ML_CS <- DAISIE_ML <- function(
                         optimmethod = optimmethod,
                         verbose = verbose,
                         tolint = tolint,
-                        jitter = jitter)
+                        jitter = jitter,
+                        num_cycles = num_cycles)
   }
   return(out)
 }
diff --git a/R/DAISIE_ML_IW.R b/R/DAISIE_ML_IW.R
@@ -90,12 +90,13 @@ DAISIE_ML_IW <- function(
   optimmethod = "subplex",
   verbose = 0,
   tolint = c(1E-16, 1E-14),
-  jitter = 0) {
+  jitter = 0,
+  num_cycles = 1) {
   options(warn = -1)
   out2err <- data.frame(lambda_c = NA, mu = NA, K = NA, gamma = NA, lambda_a = NA, loglik = NA, df = NA, conv = NA)
   out2err <- invisible(out2err)
   if (is.null(datalist[[1]]$brts_table)) {
-    datalist <- Add_brt_table(datalist)
+    datalist <- add_brt_table(datalist)
   }
   np <- datalist[[1]]$not_present
   if (is.null(np)) {
@@ -153,7 +154,8 @@ DAISIE_ML_IW <- function(
     methode = methode,
     abstolint = tolint[1],
     reltolint = tolint[2],
-    jitter = jitter
+    jitter = jitter,
+    num_cycles = num_cycles
   )
   if (out$conv != 0) {
     cat("Optimization has not converged. Try again with different initial values.\n")

diff --git a/R/DAISIE_MW_ML.R b/R/DAISIE_MW_ML.R
@@ -218,6 +218,7 @@ DAISIE_MW_loglik_choosepar = function(
 #' the 11th parameter is the d0 for colonization and the 12th is the d0 for
 #' anagenesis.
 #'
+#' @inheritParams default_params_doc
 #' @param datalist Data object containing information on colonisation and
 #' branching times of species for several islands or archipelagos, as well as the area,
 #' isolation and age of each of the islands/archipelagos. See data(archipelagos41) for
@@ -267,8 +268,6 @@ DAISIE_MW_loglik_choosepar = function(
 #' @param optimmethod Method used in likelihood optimization. Default is
 #' "subplex" (see subplex package). Alternative is 'simplex' which was the
 #' method in previous versions.
-#' @param CS_version For internal testing purposes only. Default is 1, the
-#' original DAISIE code.
 #' @param verbose sets whether parameters and likelihood should be printed (1)
 #' or not (0)
 #' @param tolint Vector of two elements containing the absolute and relative
@@ -304,6 +303,8 @@ DAISIE_MW_loglik_choosepar = function(
 #' machine.
 #' @param cpus Number of cpus used in parallel computing. Default is 3. Will
 #' not have an effect if parallel = 'no'.
+#' @param num_cycles The number of cycles the optimizer will go through.
+#'   Default is 1.
 #' @return The output is a dataframe containing estimated parameters and
 #' maximum loglikelihood.
 #' \item{lambda_c0}{ gives the maximum likelihood estimate of lambda^c,
@@ -385,7 +386,8 @@ DAISIE_MW_ML = function(
   distance_type = 'continent',
   distance_dep = 'power',
   parallel = 'local',
-  cpus = 3
+  cpus = 3,
+  num_cycles = 1
 )
 {
   options(warn=-1)
@@ -443,7 +445,24 @@ DAISIE_MW_ML = function(
   }
   cat("Optimizing the likelihood - this may take a while.","\n")
   utils::flush.console()
-  out = DDD::optimizer(optimmethod = optimmethod,optimpars = optimpars,fun = DAISIE_MW_loglik_choosepar,trparsopt = trparsopt,idparsopt = idparsopt,trparsfix = trparsfix,idparsfix = idparsfix,pars2 = pars2,datalist = datalist,methode = methode,CS_version = CS_version,abstolint = tolint[1],reltolint = tolint[2],distance_type = distance_type,parallel = parallel,cpus = cpus,distance_dep = distance_dep)
+  out = DDD::optimizer(optimmethod = optimmethod,
+                       optimpars = optimpars,
+                       fun = DAISIE_MW_loglik_choosepar,
+                       trparsopt = trparsopt,
+                       idparsopt = idparsopt,
+                       trparsfix = trparsfix,
+                       idparsfix = idparsfix,
+                       pars2 = pars2,
+                       datalist = datalist,
+                       methode = methode,
+                       CS_version = CS_version,
+                       abstolint = tolint[1],
+                       reltolint = tolint[2],
+                       distance_type = distance_type,
+                       parallel = parallel,
+                       cpus = cpus,
+                       distance_dep = distance_dep,
+                       num_cycles = num_cycles)
   if(out$conv != 0)
   {
     cat("Optimization has not converged. Try again with different initial values.\n")