Merge pull request #7 from MindTheGap-ERC/dev

Dev

R/StratPal-package.R

@@ -0,0 +1,6 @@
+#' @keywords internal
+"_PACKAGE"
+
+## usethis namespace: start
+## usethis namespace: end
+NULL

R/apply_niche_pref.R

@@ -3,13 +3,13 @@ apply_niche_pref = function(x, niche_def, gc){
   #'
   #' @title apply niche preference
   #'
-  #' @param x events, e.g. fossil occurrences
-  #' @param niche_def function, specifying the niche along a gradient. should return 0 when taxon is outside of niche, and 1 when fully inside niche. Values between 0 and 1 are interpreted as probabilities.
+  #' @param x events, e.g. times/ages of fossil occurrences or their stratigraphic position.
+  #' @param niche_def function, specifying the niche along a gradient. Should return 0 when taxon is outside of niche, and 1 when fully inside niche. Values between 0 and 1 are interpreted as probabilities.
   #' @param gc function, stands for "gradient change". Specifies how the gradient changes, e.g. with time
   #'
   #' @description
-    #' models niche preferences by removing events (fossil occurrences) when they are outside of their preferred niche using the function `thin`
-    #' Combines the functions niche_def and gc ("gradient change") to determine how the taxons preferred environment changes with time. This is done by composing `niche_def` and `gc`. The result is then used as a thinning.
+    #' Models niche preferences by removing events (fossil occurrences) when they are outside of a preferred niche using the function `thin`.
+    #' Combines the functions `niche_def` and `gc` ("gradient change") to determine how the taxons' preferred environment changes with time. This is done by composing `niche_def` and `gc`. The result is then used as a thinning on the events `x`.
   #'
   #' @examples
     #' \dontrun{
@@ -36,9 +36,13 @@ apply_niche_pref = function(x, niche_def, gc){
     #'  # fossil occurrences with niche preference
     #'  hist(foss_occ_niche, xlab = "time")
     #'
+    #'  # see also
+    #'  vignette("event_data")
+    #'  # for a longer example
+    #'
     #' }
     #'
-  #' @seealso [apply_taphonomy()] to model taphonomic effects based on the sampe principle, [thin()] for the underlying mathematical procedure
+  #' @seealso [apply_taphonomy()] to model taphonomic effects based on the same principle, [thin()] for the underlying mathematical procedure.
 
   # function that returns niche preference as a function of y (typically time)
   change_in_niche = function(y) niche_def(gc(y))

R/apply_taphonomy.R

@@ -1,17 +1,18 @@
 apply_taphonomy = function(x, pres_potential, ctc){
   #' @export
   #'
-  #' @title model taphonomy
+  #' @title Model taphonomic effects
   #'
-  #' @param x numeric vector. event type data, e.g. fossil occurrences
-  #' @param pres_potential function. takes taphonomic conditions as input and returns the preservation potential (a number between 0 and 1)
-  #' @param ctc function, change in taphonomic conditions (ctc) with time.
+  #' @param x events, e.g. times/ages of fossil occurrences or their stratigraphic position.
+  #' @param pres_potential function. Takes taphonomic conditions as input and returns the preservation potential (a number between 0 and 1)
+  #' @param ctc function, change in taphonomic conditions (ctc) with time or stratigraphic position.
   #'
   #' @description
-    #' models taphonomy by combining the change in taphonomic conditions with the preservation potential as a function of taphonomic conditions to determine how preservation potential changes. This is then used to systematically remove (thin) the event data using `thin`
+    #' Models taphonomy by combining the change in taphonomic conditions with the preservation potential as a function of taphonomic conditions to determine how preservation potential changes. This is then used to systematically remove (thin) the event data using `thin`.
     #'
     #'
-  #' @seealso [apply_niche_pref()] for modeling niche preferences based on the same principle, [thin()] for the underlying mathematical procedure
+  #' @seealso [apply_niche_pref()] for modeling niche preferences based on the same principle, [thin()] for the underlying mathematical procedure.
+  #'
   #'
 
   # function that returns preservation potential as a function of input (e.g. time or position)

R/cnd.R

@@ -13,19 +13,25 @@ cnd = function(mean, sd, inc = 1, cut_neg = TRUE){
   #'
   #' @examples
     #' \dontrun{
-    #' # using water depht as niche
+    #' # using water depth as niche
     #' wd = seq(-3, 40, by = 0.5)
     #' f = cnd(mean = 10, sd = 5, inc = 15, cut_neg = FALSE)
     #' # 1 indicates high preference, 0 indicates low preference
     #' plot(wd, f(wd), xlab = "Water depth", ylab = "Env. preference")
     #' # set value at neg wd to 0 - non-terrestrial species.
     #' f = cnd(mean = 10, sd = 5, inc = 15, cut_neg = TRUE)
     #' plot(wd, f(wd), xlab = "Water depth", ylab = "Env. preference")
+    #'
+    #' # see also
+    #' vignette("event_data")
+    #' #for examples how to use it for niche modeling
     #' }
     #'
   #'
   #' @description
-    #' returns a function that defines a niche based on a capped normal distribution, i.e. a pdf of a normal distribution where all values above 1 are capped. Mathematically, this is f(x) = min( inc * pdf(x), 1)
+    #' Returns a function that defines a niche based on a capped normal distribution, i.e. a pdf of a normal distribution where all values above 1 are capped. Mathematically, this is f(x) = min( inc * pdf(x), 1).
+    #' Corresponds to the probability of collection model by Holland and Patzkowsky (2002).
+  #' @references * Holland, Steven M. and Patzkowsky, Mark E. 2002. "Stratigraphic variation in the timing of first and last occurrences." PALAIOS. https://doi.org/10.1669/0883-1351(2002)017<0134:SVITTO>2.0.CO;2
 
   fa = function(x){
     y = pmin(inc * stats::dnorm(x, mean, sd), rep(1, length(x)))

R/ornstein_uhlenbeck.R

@@ -3,16 +3,16 @@ ornstein_uhlenbeck = function(t, mu = 0, theta = 1, sigma = 1, y0 = 0){
   #'
   #' @title simulate ornstein-uhlenbeck (OU) process
   #'
-  #' @param t times at which the process is simulated
-  #' @param mu numeric, long term mean
-  #' @param theta numeric, mean reversion speed
+  #' @param t times at which the process is simulated. Can be heterodistant
+  #' @param mu scalar, long term mean
+  #' @param theta scalar, mean reversion speed
   #' @param sigma positive scalar, strength of randomness
-  #' @param y0 scalar. initial value (value of process at the first entry of t)
+  #' @param y0 scalar, initial value (value of process at the first entry of t)
   #'
   #' @description
-    #' simulates the Ornstein-Uhlenbeck process using the Euler-Maruyame method. The process is simulated on a scale of 0.25 * min(diff(t)) and then interpolated to the values of `t`.
+    #' Simulates an Ornstein-Uhlenbeck process using the Euler-Maruyama method. The process is simulated on a scale of `0.25 * min(diff(t))` and then interpolated to the values of `t`.
     #'
-  #' @returns a list with two elements: `t` and `y`. `t` is a duplicate of the input `t`, `y` are the values of the OU process at these times. Outputs are of class `timelist` and can thus be plotted directly using `plot`, see `?plot.timelist`
+  #' @returns A list with two elements: `t` and `y`. `t` is a duplicate of the input `t`, `y` are the values of the OU process at these times. Output list is of class `timelist` and can thus be plotted directly using `plot`, see `?admtools::plot.timelist`
   #'
   #' @examples
     #' \dontrun{

R/p3.R

@@ -12,19 +12,27 @@ p3 = function(rate, from, to, n = NULL){
   #' @description
     #' Simulates events in the interval `from` to `to` based on a Poisson point process with rate `rate`. If the parameter `n` is used, the number of fossils is conditioned to be `n`
     #' In the context of paleontology, these events can be interpreted as fossil occurrences or first/last occurrences of species. In this case, the rate is the average number of fossil occurrences (resp first/last occurrences) per unit
+  #'
+  #' @returns a numeric vector with timing/location of events.
+  #'
   #' @examples
     #' \dontrun{
     #' # for fossil occ.
     #' x = p3(rate = 5, from = 0, to = 1) # 5 fossil occurrences per myr on avg.
-    #' hist(x, xlab = "Time (Myr"), ylab = "Fossil Occurrences" )
+    #' hist(x, xlab = "Time (Myr)", ylab = "Fossil Occurrences" )
     #'
     #' x = p3(rate = 3, from = 0, to = 4)
     #' hist(x, main = paste0(length(x), " samples")) # no of events is random
     #'
     #' x = p3(rate = 3, from = 0, to = 4, n = 10)
     #' hist(x, main = paste0(length(x), " samples")) # no of events is fixed to n
+    #'
+    #' # see also
+    #' vignette("event_data")
+    #' # for details on usage and applications to paleontology
     #' }
     #'
+  #' @seealso [p3_var_rate()] for the variable rate implementation
 
   if (rate <= 0){ stop("Need strictly positive rate")}
   if (from >= to) {stop("parameter \"from\" needs to be smaller than \"to\".")}

R/p3_var_rate.R

@@ -18,18 +18,23 @@ p3_var_rate = function(x, y = NULL, from = 0, to = 1, f_max = 1, n = NULL){
     #' \dontrun{
     #' # assuming events are fossil occurrences
     #' # then rate is the avg rate of fossil occ. per unit
-    #' linear decrease in rate from 50 at x = 0 to 0 at x = 1
+    #' #linear decrease in rate from 50 at x = 0 to 0 at x = 1
     #' x = c(0, 1)
     #' y = c(50, 0)
     #' s = p3_var_rate(x, y, f_max = 50)
     #' hist(s, xlab = "Time (myr)", main = "Fossil Occurrences")
-    #' # conditoned to return 100 samples
+    #' # conditioned to return 100 samples
     #' s = p3_var_rate(x, y, f_max = 50, n = 100)
     #' # hand over function
     #' s = p3_var_rate(x = sin, from = 0 , to = 3 * pi, n = 50)
     #' hist(s) # note that negative values of f (sin) are ignored in sampling
+    #'
+    #' # see also
+    #' vignette("event_data")
+    #' # for details on usage and applications to paleontology
     #' }
     #'
+  #' @seealso [p3()] for the constant rate implementation, [rej_samp()] for the underlying random number generation.
 
   if (from >= to){
     stop("\"from\" must be smaller than \"to\".")

R/random_walk.R

@@ -3,15 +3,15 @@ random_walk = function(t, sigma = 1, mu = 0, y0 = 0){
   #'
   #' @title continuous time random walk
   #'
-  #' @param t numeric vector with strictly increasing elements. Times at which the random walk is evaluated
-  #' @param sigma variance parameter
-  #' @param mu directionality parameter
-  #' @param y0 starting value, i.e. value of the random walk at the first entry of `t`
+  #' @param t numeric vector with strictly increasing elements, can be heterodistant. Times at which the random walk is evaluated
+  #' @param sigma positive scalar, variance parameter
+  #' @param mu scalar, directionality parameter
+  #' @param y0 scalar, starting value (value of the random walk at the first entry of `t`)
   #'
   #' @description
-    #' simulates a random walk as a Brownian drift
+    #' Simulates a (continuous time) random walk as a Brownian drift
     #'
-  #' @returns a list with elements `t` and `y`. `t` is a duplicate of the input parameter and is the times at which the random walk is evaluated. `y` are the values of the random walk at said times
+  #' @returns A list with elements `t` and `y`. `t` is a duplicate of the input parameter and is the times at which the random walk is evaluated. `y` are the values of the random walk at said times. Output list is of class `timelist` and can thus be plotted directly using `plot`, see `?admtools::plot.timelist`
   #'
   #' @examples
     #' \dontrun{

R/rej_samp.R

@@ -4,7 +4,7 @@ rej_samp = function(f, x_min, x_max, n = 1L, f_max = 1){
   #' @title rejection sampling
   #'
   #' @description
-    #' rejection sampling from the (pseudo) pdf `f` in the interval between `x_min` and `x_max`. Returns `n` samples. Note that values of `f` below 0 are capped to zero
+    #' Rejection sampling from the (pseudo) pdf `f` in the interval between `x_min` and `x_max`. Returns `n` samples. Note that values of `f` below 0 are capped to zero
     #'
   #' @param f function. (pseudo) pdf from which the sample is drawn
   #' @param x_min scalar. lower limit of the examined interval
@@ -18,7 +18,7 @@ rej_samp = function(f, x_min, x_max, n = 1L, f_max = 1){
     #' x = rej_samp(f, 0, 3*pi, n = 100)
     #' hist(x) # note that no samples are drawn where sin is negative
     #' }
-    #'
+  #' @seealso [p3_var_rate()] for the derived variable rate Poisson point process implementation.
   x = c()
   warn = FALSE
   if (f_max <= 0) {stop("`f_max` must be positive.")}

R/scenarioA.R

@@ -1,4 +1,8 @@
-#' @title scenario A from Hohmann et al 2024
+#' @title example data, scenario A from Hohmann et al. (2024)
+#'
+#' @description
+  #' Scenario A as described in Hohmann et al. (2024), published in Hohmann et al. (2023). Contains data from a carbonate platform simulated using CarboCAT Lite (Burgess 2013, 2023)
+#'
 #'
 #' @format A list with 6 elements:
 #'
@@ -7,7 +11,14 @@
 #' * `dist_from_shore` : character vector. Distance from shore in km of locations at which the observations were made
 #' * `h_m` : matrix of size length(t_myr) x length(dist_from_shore). Accumulated sediment height in m at examined locations
 #' * `wd_m`: matrix of size length(t_myr) x length(dist_from_shore). Water depth in m at examined locations
-#' * `strat_col`: list with length(dist_from shore) elements. Represents a stratigraphic column. Each element is a list with two elements
-#'  * `bed_thickness_m`: numeric vector. Bed thickness in m
-#'  * `facies_code` : integer vector. facies code of the bed
+#' * `strat_col`: list with length(dist_from shore) elements. Represents a stratigraphic column. Each element is a list with two elements:
+#'      * `bed_thickness_m`: numeric vector. Bed thickness in m
+#'      * `facies_code` : integer vector. facies code of the bed
+#'
+#'
+#' @references
+  #' * Burgess, Peter. 2013. "CarboCAT: A cellular automata model of heterogeneous carbonate strata." Computers & Geosciences. <https://doi.org/10.1016/j.cageo.2011.08.026>.
+  #' * Burgess, Peter. 2023. "CarboCATLite v1.0.1." Zenodo. <https://doi.org/10.5281/zenodo.8402578>
+  #' * Hohmann, Niklas; Koelewijn, Joël R.; Burgess, Peter; Jarochowska, Emilia. 2024. "Identification of the mode of evolution in incomplete carbonate successions." BMC Ecology and Evolution, In Press. https://doi.org/10.1101/2023.12.18.572098.
+  #' * Hohmann, Niklas, Koelewijn, Joël R.; Burgess, Peter; Jarochowska, Emilia. 2023. “Identification of the Mode of Evolution in Incomplete Carbonate Successions - Supporting Data.” Open Science Framework. https://doi.org/10.17605/OSF.IO/ZBPWA, published under the [CC-BY 4.0](https://creativecommons.org/licenses/by/4.0/) license.
 "scenarioA"

R/stasis.R

@@ -3,13 +3,14 @@ stasis = function(t, mean = 0, sd = 1){
   #'
   #' @title simulate phenotypic stasis
   #'
-  #' @param t times at which the phenotype is determined
-  #' @param mean mean value
-  #' @param sd standard deviation
+  #' @param t times at which the traits are determined
+  #' @param mean scalar, mean trait value
+  #' @param sd strictly positive scalar, standard deviation of traits
   #'
   #' @description
-    #' simulates stasis as independent, normally distributed random variables with mean `mean` and standard deviation `sd`
+    #' Simulates stasis as independent, normally distributed random variables with mean `mean` and standard deviation `sd`
     #'
+  #' @returns A list with two elements: `t` and `y`. `t` is a duplicate of the input `t`, `y` are the corresponding trait values. Output list is of class `timelist` and can thus be plotted directly using `plot`, see `?admtools::plot.timelist`
   #' @examples
     #' \dontrun{
     #'  library("admtools") # required for plotting of results

R/thin.R

@@ -9,17 +9,18 @@ thin = function(x, thin){
   #'
   #'
   #' @description
-    #' Thins a vector of events using the function thin, meaning the probability that the ith event in x is preserved is given by _thin (x(i))_. Values of
-    #' `thin` below 0 and above 1 are ignored
-    #' Is used to model niche preferences in apply_niche_pref, where events are thinned based on how close they are to their preferred niche
-    #' Can also be used to model taphonomic effects. In this case, `thin` describes the preservation potential.
+    #' Thins a vector of events using the function thin, meaning the probability that the ith event in x is preserved is given by _thin(x(i))_. Values of
+    #' `thin` below 0 and above 1 are ignored.
+    #' Is used to model niche preferences in `apply_niche_pref` and taphonomic effects in `apply_taphonomy`.
     #'
   #' @examples
+    #' \dontrun{
     #' x = p3(rate = 100, from = 0, to = 3 * pi) # simulate Poisson point process
     #' y = thin(x, sin)
     #' hist(y) # not how negative values of sin are treated as 0
     #' yy = thin(x, function(x) 5 * sin(x))
     #' hist(yy) # note how values of 5 * sin above 1 are not affecting the thinning
+    #' }
     #'
   #' @seealso [apply_niche_pref()] and [apply_taphonomy()] for use cases with biological meaning