R/RI.R In phenology: Tools to Manage a Parametric Function that Describes Phenology and More

Documented in RI

```#' RI returns an expected remigration interval
#' @title Return an expected remigration interval.
#' @author Marc Girondot
#' @return Return a remigration interval.\cr
#' @param s Time-conditional probability of survival
#' @param t Time-conditional probability of tag retention
#' @param r Time-conditional probability of return
#' @param c Probability of return
#' @param p Annual probability of observation
#' @description Model of remigration interval\cr
#' Note that r, s and t are conditional probabilities. If c is null, then return probabilities are
#' estimated from r. r can be named vector. For example:\cr
#' r <- c(r1=0.5, r2=0.60, r3=1) is equivalent to c <- c(c1=0.5, c2=0.3, c3=0.2)\cr
#' The vector of r described the probability that a female returned after
#' 1, 2, 3 years among those who have not nested before.
#' The vector of r is the same but defining the return probability for an initial female.\cr
#' @family Model of Remigration Interval
#' @examples
#' \dontrun{
#' library(phenology)
#' # Example
#' s <- c(s1=1, s2=1, s3=1, s4=1, s5=1)
#' t <- c(t1=0.95, t2=1, t3=1, t4=1, t5=1)
#' r <- c(r1=0.1, r2=0.8, r3=0.7, r4=0.7, r5=1)
#' p <- c(p1=0.6, p2=0.6, p3=0.6, p4=0.6, p5=0.6)
#'
#' # r is equivalent to
#' c <- c(c1=0.1, c2=0.72, c3=0.126, c4=0.0378, c5=0.0162)
#' # Then the true remigration interval is:
#' ri_true <- sum(1:5*c[1:5])
#'
#' s_ri  <- NULL
#' for (sx in seq(from=0.01, to=1, by=0.01)) {
#'   s[] <- sx
#'   ri1 <- RI(s=s, t=t, r=r, p=p)
#'   s_ri  <- c(s_ri,sum(1:5*ri1)/sum(ri1))
#' }
#'
#' par(mar=c(4, 4, 1, 1)+0.4)
#'
#' plot(x=seq(from=0.01, to=1, by=0.01), y=s_ri, type="l",
#'      las=1, bty="n", ylim=c(0, 4),
#'      xlab="Annuual survival probabilities", ylab="Naive Remigration Interval",
#'     main="")
#' segments(x0=0.01, x1=1, y0=ri_true, y1=ri_true, lty=2, col="red")
#' legend("topright", legend="True remigration interval", lty=2, col="red")
#'
#' }
#' @export

RI <- function(s, t, r=NULL, c=NULL, p) {

if (names(s[1]) != "s") {
s <- s[order(as.numeric(substr(names(s), 2, nchar(s))))]
}

if (names(t[1]) != "t") {
t <- t[order(as.numeric(substr(names(t), 2, nchar(t))))]
}

if (!is.null(r)) {
if (names(r[1]) != "r") {
r <- r[order(as.numeric(substr(names(r), 2, nchar(r))))]
}
}

if (!is.null(c)) {
if (names(c[1]) != "c") {
c <- c[order(as.numeric(substr(names(c), 2, nchar(c))))]
}
}

if (names(p[1]) != "p") {
p <- p[order(as.numeric(substr(names(p), 2, nchar(p))))]
}

if (is.null(c)) {

rp <- rep(NA, length(r))
rp[1] <- r[1]
if (length(r) != 1) {
for (i in 2:length(r)) rp[i] <- r[i]*prod(1 - r[1:(i-1)])
}
r <- rp
} else {
r <- c
}

# r <- r/sum(r)
k <- matrix(data = 1, nrow=1)
N <- 1

# Note that r, s and t are conditional probabilities

n.return.vue <- NULL

for (l in 1:length(s)) {
n.survival <- N * prod(s[1:l])
n.tagretention <- n.survival * prod(t[1:l])
n.return <- 0
for (h in 1:nrow(k)) {
ts <- k[h, , drop=TRUE]
j <- 1
l <- 1
if (length(ts) != 1) {
for (n in 1:(length(ts)-1)) {
if (ts[n]==1) {
j <- j*r[l]*(1-p[n])
l <- 1
} else {
l <- l + 1
}
}
}
j <- j*r[l]*p[length(s)]
n.return <- n.return + j * n.tagretention
}
n.return.vue <- c(n.return.vue, n.return)
k2 <- k
k2[, ncol(k)] <- 0
k3 <- rbind(k, k2)
k4 <- cbind(k3, rep(1, nrow(k3)))
k <- k4
}
return(unname(n.return.vue))
}
```

Try the phenology package in your browser

Any scripts or data that you put into this service are public.

phenology documentation built on Oct. 16, 2023, 9:06 a.m.