More enhancements to marray* functions

9 files changed, 127 insertions, 82 deletions

diff --git a/.local_checks/@making_IPMbook_package.R b/.local_checks/@making_IPMbook_package.R
index a624970..ca531a6 100644
--- a/.local_checks/@making_IPMbook_package.R
+++ b/.local_checks/@making_IPMbook_package.R
@@ -20,7 +20,7 @@ devtools::load_all("D:/Github/IPMbook_package/IPMbook")
 # ==========================
 unlink(list.files(pattern="Rplots.pdf", recursive=TRUE))
 system("R CMD build IPMbook")  # Produces the .tar.gz
-pkg <- "IPMbook_0.1.2.9003.tar.gz"  # <-- fix version number here
+pkg <- "IPMbook_0.1.2.9004.tar.gz"  # <-- fix version number here
 
 # Pick one to check:
 ## on desktop
diff --git a/DESCRIPTION b/DESCRIPTION
index 80fd109..3382ff2 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,8 +1,8 @@
 Package: IPMbook
 Type: Package
 Title: Functions and Data for the Book 'Integrated Population Models'
-Version: 0.1.2.9003
-Date: 2021-09-29
+Version: 0.1.2.9004
+Date: 2021-10-01
 Depends: R (>= 2.10)
 Imports: stats, abind
 Suggests:
diff --git a/NEWS b/NEWS
index 6185fbf..9882d28 100644
--- a/NEWS
+++ b/NEWS
@@ -1,3 +1,9 @@
+CHANGES in version 0.1.2.9004 (2021-10-01)
+
+  * `marray*` functions now deal correctly with trap losses if only 1 capture.
+
+  * `marray*` functions deal with with groups specified as columns in the `freq` argument.
+
 CHANGES in version 0.1.2.9003 (2021-09-29)
 
   * `marray*` functions now accept data frames for capture histories.
@@ -16,6 +22,6 @@ CHANGES in version 0.1.2.9000 (2021-09-25)
 
 CHANGES in version 0.1.2 (2021-09-14)
 
-  * Corrected publication date of Schaube and Kéry (2021) to (2022) in all help pages, though it will ship in Q3 2021.
+  * Corrected publication date of Schaub and Kéry (2021) to (2022) in all help pages, though it will ship in Q3 2021.
 
   * Changed the book URL to the permanent site at Swiss Ornithological Institute.
diff --git a/R/marray.R b/R/marray.R
index 4998238..1a666c7 100644
--- a/R/marray.R
+++ b/R/marray.R
@@ -23,47 +23,51 @@
 marray <- function(ch, unobs = 0, freq = 1){
   if(is.data.frame(ch))
     ch <- as.matrix(ch)
+  if(!is.matrix(ch))
+    ch <- matrix(ch, nrow=1)
   ch <- round(ch)
   stopifNegative(ch, allowNA=TRUE, allowZero=TRUE)
   unobs <- round(unobs)
   stopifNegative(unobs, allowNA=FALSE, allowZero=TRUE)
-  freq <- round(freq)
+  if(is.data.frame(freq))
+    freq <- as.matrix(freq)
   if(length(freq) == 1)
     freq <- rep(freq, nrow(ch))
+  if(!is.matrix(freq))
+    freq <- matrix(freq, ncol=1)
+  freq <- round(freq)
+  absfreq <- abs(freq)
 
-  if(!is.matrix(ch))
-    ch <- matrix(ch, nrow=1)
-  ns <- length(table(ch)) - 1 + unobs # number of states, excluding 0
+  ns <- length(table(ch)) - 1 + unobs # number of states, excluding 0, can be 1
   no <- ncol(ch)                      # number of observations (replicates, years, surveys, ...)
-
-  # Remove capture histories of individuals that are marked at last occasion
-  first <- getFirst(ch)
-  last <- which(first==no)
-  if (length(last) > 0) {
-    ch <- ch[-last,]
-    freq <- freq[-last]
-  }
-  absfreq <- abs(freq)
+  ng <- ncol(freq)                    # number of groups, can be 1
 
   # Check for trap losses
   traploss <- freq < 0
   if(all(traploss == FALSE))
-    traploss <- is.na(rowSums(ch))
+    traploss <- matrix(is.na(rowSums(ch)), nrow=nrow(ch), ncol=ng)
 
   # Create empty m-array, add dimnames
-  out <- matrix(0, ncol = ns*(no-1)+1, nrow = ns*(no-1))
+  out <- array(0, dim = c(ns*(no-1), ns*(no-1)+1, ng))
+
+  gNames <- colnames(freq)
+  if(is.null(gNames))
+    gNames <- paste0("G", 1:ng)
   if(ns == 1) {
     dimnames(out) <- list(released = paste0("Y", 1:(no-1)),
-        recaptured = c(paste0("Y", 2:no), "never"))
+        recaptured = c(paste0("Y", 2:no), "never"), group=gNames)
   } else {
     YStmp <- expand.grid(paste0("S", 1:ns), paste0("Y", 1:no))
     Y.S <- paste(YStmp[,2], YStmp[,1], sep=".")
     dimnames(out) <- list(released = Y.S[1:(ns*(no-1))],
-        recaptured = c(Y.S[-(1:ns)], "never"))
+        recaptured = c(Y.S[-(1:ns)], "never"), group=gNames)
   }
 
   # Insert values in m-array
+  first <- getFirst(ch)
   for (i in 1:nrow(ch)){
+    if(first[i] == no)   # Skip cases where first capture is last occasion
+      next
     cap.occ <- which(ch[i,]!=0)
     state <- ch[i,cap.occ]
     if(ns == 1) {
@@ -71,19 +75,22 @@ marray <- function(ch, unobs = 0, freq = 1){
     } else {
       capID <- paste0("Y", cap.occ, ".S", state)
     }
-    if (length(state) == 1) {  # never recaptured
-       out[capID[1], 'never'] <- out[capID[1], 'never'] + absfreq[i]
-    }
-    if (length(state) > 1) {  # recaptured at least once
-      for (t in 2:length(cap.occ)){
-        out[capID[t-1], capID[t]] <- out[capID[t-1], capID[t]] + absfreq[i]
+
+    for(g in 1:ng) {
+      if (length(state) == 1 && !traploss[i, g]) {  # no recaptures
+         out[capID[1], 'never', g] <- out[capID[1], 'never', g] + absfreq[i, g]
       }
-      if (max(cap.occ) < no && !traploss[i]){  # never recaptured after last release
-        out[capID[t], 'never'] <- out[capID[t], 'never'] + absfreq[i]
+      if (length(state) > 1) {  # recaptured at least once
+        for (t in 2:length(cap.occ)){
+          out[capID[t-1], capID[t], g] <- out[capID[t-1], capID[t], g] + absfreq[i, g]
+        }
+        if (max(cap.occ) < no && !traploss[i, g]){  # never recaptured after last release
+          out[capID[t], 'never', g] <- out[capID[t], 'never', g] + absfreq[i, g]
+        }
       }
     }
   }
 
-  return(out)
+  return(drop(out))
 }
 
diff --git a/R/marrayAge.R b/R/marrayAge.R
index 0c46189..d3e5ecc 100644
--- a/R/marrayAge.R
+++ b/R/marrayAge.R
@@ -23,40 +23,46 @@
 marrayAge <- function(ch, age, mAge = 1, freq = 1){
   if(is.data.frame(ch))
     ch <- as.matrix(ch)
+  if (!is.matrix(ch))
+    ch <- matrix(ch, nrow = 1)
   ch <- round(ch)
   stopifNegative(ch, allowNA=TRUE, allowZero=TRUE)
   age <- round(age)
   stopifNegative(age, allowNA=FALSE, allowZero=FALSE)
   mAge <- round(mAge[1])
   stopifNegative(mAge, allowNA=FALSE, allowZero=FALSE)
-  freq <- round(freq)
+  if(is.data.frame(freq))
+    freq <- as.matrix(freq)
   if(length(freq) == 1)
     freq <- rep(freq, nrow(ch))
+  if(!is.matrix(freq))
+    freq <- matrix(freq, ncol=1)
+  freq <- round(freq)
+  absfreq <- abs(freq)
 
-  if (!is.matrix(ch))
-    ch <- matrix(ch, nrow = 1)
   maxAge <- max(age, mAge)
   nind <- nrow(ch)
   stopifnotLength(age, nind, allow1=FALSE)
-  stopifnotLength(freq, nind, allow1=FALSE)
+  # stopifnotLength(freq, nind, allow1=FALSE)  ## FIXME
   n.occasions <- ncol(ch)
-
+  ng <- ncol(freq)   # number of groups, can be 1
+  
   # Remove capture histories of individuals that are marked at last occasion
   first <- getFirst(ch)
-  last <- which(first == n.occasions)
-  if (length(last) > 0) {
-    ch <- ch[-last,]
-    age <- age[-last]
-    freq <- freq[-last]
-    first <- getFirst(ch)
-    nind <- nrow(ch)
-  }
+  # last <- which(first == n.occasions)
+  # if (length(last) > 0) {
+    # ch <- ch[-last,]
+    # age <- age[-last]
+    # freq <- freq[-last]
+    # first <- getFirst(ch)
+    # nind <- nrow(ch)
+  # }
   absfreq <- abs(freq)
 
-    # Check for trap losses
+  # Check for trap losses
   traploss <- freq < 0
   if(all(traploss == FALSE))
-    traploss <- is.na(rowSums(ch))
+    traploss <- matrix(is.na(rowSums(ch)), nrow=nrow(ch), ncol=ng)
 
   age.matrix <- matrix(0, ncol = n.occasions, nrow = nind)
   for (i in 1:nind){
@@ -64,29 +70,41 @@ marrayAge <- function(ch, age, mAge = 1, freq = 1){
   }
   age.matrix[age.matrix > maxAge] <- maxAge
 
-  marr <- array(0, dim = c(n.occasions-1, n.occasions, maxAge))
+  # Create empty m-array, add dimnames
+  marr <- array(0, dim = c(n.occasions-1, n.occasions, maxAge, ng))
+  
+  gNames <- colnames(freq)
+  if(is.null(gNames))
+    gNames <- paste0("G", 1:ng)
   dimnames(marr) <- list(released = paste0("Y", 1:(n.occasions-1)),
         recaptured = c(paste0("Y", 2:n.occasions), "never"),
-        age = 1:maxAge)
+        age = 1:maxAge,
+        group = gNames)
 
+  # Insert values in m-array
+  first <- getFirst(ch)
   for (i in 1:nind){
+    if(first[i] == n.occasions)   # Skip cases where first capture is last occasion
+      next
     cap.occ <- which(ch[i,]!=0)
     capID <- paste0("Y", cap.occ)
     cap.age <- age.matrix[i, cap.occ]
-
-    if (length(capID) == 1) {  # never recaptured
-       marr[capID[1], 'never', cap.age[1]] <- marr[capID[1], 'never', cap.age[1]] + absfreq[i]
-    }
-    if (length(capID) > 1) {  # recaptured at least once
-      for (t in 2:length(cap.occ)){
-        marr[capID[t-1], capID[t], cap.age[t-1]] <- marr[capID[t-1], capID[t], cap.age[t-1]] + absfreq[i]
+    
+    for(g in 1:ng) {
+      if (length(capID) == 1 && !traploss[i,g]) {  # no recaptures
+         marr[capID[1], 'never', cap.age[1], g] <- marr[capID[1], 'never', cap.age[1],g] + absfreq[i,g]
       }
-      if (max(cap.occ) < n.occasions && !traploss[i]){  # never recaptured after last release
-        marr[capID[t], 'never', cap.age[t]] <- marr[capID[t], 'never', cap.age[t]]  + absfreq[i]
+      if (length(capID) > 1) {  # recaptured at least once
+        for (t in 2:length(cap.occ)){
+          marr[capID[t-1], capID[t], cap.age[t-1], g] <- marr[capID[t-1], capID[t], cap.age[t-1], g] + absfreq[i, g]
+        }
+        if (max(cap.occ) < n.occasions && !traploss[i,g]){  # never recaptured after last release
+          marr[capID[t], 'never', cap.age[t], g] <- marr[capID[t], 'never', cap.age[t], g]  + absfreq[i, g]
+        }
       }
     }
   }
 
-  return(marr)
+  return(drop(marr))
 }
 
diff --git a/R/marrayDead.R b/R/marrayDead.R
index 6527374..188f085 100644
--- a/R/marrayDead.R
+++ b/R/marrayDead.R
@@ -12,32 +12,46 @@ marrayDead <- function(MR, freq = 1){
     MR <- matrix(MR, nrow=1)
   if(any(rowSums(MR) > 2))
     stop("The rows of MR may not have more that two 1's", call.=FALSE)
-  freq <- round(freq)
+  if(is.data.frame(freq))
+    freq <- as.matrix(freq)
   if(length(freq) == 1)
     freq <- rep(freq, nrow(MR))
+  if(!is.matrix(freq))
+    freq <- matrix(freq, ncol=1)
+  freq <- round(freq)
+  absfreq <- abs(freq)
 
   nind <- nrow(MR)
   n.occasions <- ncol(MR)
-  out <- matrix(0, ncol=n.occasions, nrow=n.occasions-1)
+  ng <- ncol(freq)                    # number of groups, can be 1
+
+  # Create empty m-array, add dimnames
+  out <- array(0, dim=c(n.occasions-1, n.occasions, ng))
+  
+  gNames <- colnames(freq)
+  if(is.null(gNames))
+    gNames <- paste0("G", 1:ng)
   dimnames(out) <- list(released = paste0("Y", 1:(n.occasions-1)),
-      recovered = c(paste0("Y", 2:n.occasions), "never"))
+      recovered = c(paste0("Y", 2:n.occasions), "never"),
+      gNames)
 
   # Create vector with occasion of marking
   f <- getFirst(MR)  # year of release
   f.fact <- factor(f, levels=1:n.occasions)
-  # Calculate the number of released individuals at each time period
-  released <- tapply(freq, f.fact, sum)
-  released[is.na(released)] <- 0  # tapply returns NA if a value is missing from f.
+  for(g in 1:ng) {
+    # Calculate the number of released individuals at each time period
+    released <- tapply(freq[,g], f.fact, sum)
+    released[is.na(released)] <- 0  # tapply returns NA if a value is missing from f.
 
-  # Fill m-array with recovered individuals
-  rec.ind <- which(apply(MR, 1, sum)==2)  # which were recovered dead
-  rec <- getFirst(rmFirst(MR[rec.ind, ])) # year of recovery
-  for (i in seq_along(rec.ind)){
-    out[f[rec.ind[i]],rec[i]-1] <- out[f[rec.ind[i]],rec[i]-1] + abs(freq[rec.ind[i]])
+    # Fill m-array with recovered individuals
+    rec.ind <- which(apply(MR, 1, sum)==2)  # which were recovered dead
+    rec <- getFirst(rmFirst(MR[rec.ind, ])) # year of recovery
+    for (i in seq_along(rec.ind)){
+      out[f[rec.ind[i]],rec[i]-1, g] <- out[f[rec.ind[i]],rec[i]-1, g] + abs(freq[rec.ind[i], g])
+    }
+    # Calculate the number of individuals that are never recovered
+    out[ ,n.occasions, g] <- released[-n.occasions] - rowSums(out[,,g])
   }
-  # Calculate the number of individuals that are never recovered
-  out[ ,n.occasions] <- released[-n.occasions] - rowSums(out)
-
-  return(out)
+  return(drop(out))
 }
 
diff --git a/man/marray.Rd b/man/marray.Rd
index 356f5e2..36cbef5 100644
--- a/man/marray.Rd
+++ b/man/marray.Rd
@@ -6,27 +6,27 @@
 Converts capture-histories to an m-array for one age class
 }
 \description{
-Creates an m-array for a single- or multistate capture-recapture data with one age class.
+Creates an m-array for a single- or multistate capture-recapture data with one age class and optionally more than one group.
 }
 \usage{
 marray(ch, unobs = 0, freq = 1)
 }
 \arguments{
   \item{ch}{
-an individuals x time matrix with single- or multistate capture histories (0: not captured; 1...X: captured in the 1...X states). This can be a matrix of unique capture histories accompanied by a vector, \code{freq}, specifying the number of animals with each capture history. Trap losses can be indicated either by negative values for \code{freq}, or by filling the row with NA after the last capture.
+an individuals x time matrix with single- or multistate capture histories (0: not captured; 1...X: captured in the 1...X states). This can be a matrix of unique capture histories accompanied by a vector or matrix, \code{freq}, specifying the number of animals with each capture history. Trap losses can be indicated either by negative values for \code{freq}, or by filling the row with NA after the last capture.
 }
   \item{unobs}{
 number of unobserved states (default is 0, needs to be given only in specific cases).
 }
   \item{freq}{
-a vector with the number of animals with each capture history. If a single value is supplied, it will be used for all rows in the capture history. The default is to assume each row corresponds to a single animal.
+a vector with the number of animals with each capture history, or a matrix with a column for each group. If a single value is supplied, it will be used for all rows in the capture history; the default is to assume each row corresponds to a single animal.
 }
 
 }
 \value{
-For single-state capture recapture data, an m-array which is a (years-1) x years matrix, where element [i, j] contains the number of individuals released in year i and recaptured in year j+1 (by definition no recaptures can occur in year 1). The last column contains the number of individuals released in year i and never recaptured.
+For single-state capture recapture data, an m-array which is a (years-1) x years x groups array, where element [i, j, g] contains the number of individuals in group g released in year i and recaptured in year j+1 (by definition no recaptures can occur in year 1). If no groups are specified, this will be a (years-1) x years matrix. The last column contains the number of individuals released in year i and never recaptured.
 
-For multi-state capture-recapture data with s states (including potential unobservable states), an m-array which is a (years-1)*s x (years-1)*2 + 1 matrix. An element [i, j] contains the number of individuals released in year t and state n (i = (t-1)*s + n) and recaptured in year k+1 in state m (j = (k-1)*s + m). The last column contains the number of individuals released in year t and state s and never recaptured. The labeling of unobserved states starts with the number of observed states + 1.
+For multi-state capture-recapture data with s states (including potential unobservable states), an m-array which is a (years-1)*s x (years-1)*2 + 1 x groups array. An element [i, j, g] contains the number of individuals in group g released in year t and state n (i = (t-1)*s + n) and recaptured in year k+1 in state m (j = (k-1)*s + m). If no groups are specified, this will be a matrix. The last column contains the number of individuals released in year t and state s and never recaptured. The labeling of unobserved states starts with the number of observed states + 1.
 }
 \author{
 Michael Schaub
diff --git a/man/marrayAge.Rd b/man/marrayAge.Rd
index 753c317..86a449b 100644
--- a/man/marrayAge.Rd
+++ b/man/marrayAge.Rd
@@ -13,7 +13,7 @@ marrayAge(ch, age, mAge = 1, freq = 1)
 }
 \arguments{
   \item{ch}{
-an individuals x time matrix with capture histories (0: not captured; 1: captured). Alternativel, this can be a matrix of unique capture histories accompanied by a vector, \code{freq}, specifying the number of animals with each combination of capture history and age. Trap losses can be indicated either by negative values for \code{freq}, or by filling the row with NA after the last capture.
+an individuals x time matrix with capture histories (0: not captured; 1: captured). Alternatively, this can be a matrix of unique capture histories accompanied by a vector or matrix, \code{freq}, specifying the number of animals with each combination of capture history and age. Trap losses can be indicated either by negative values for \code{freq}, or by filling the row with NA after the last capture.
 }
   \item{age}{
 vector with the age class at first capture for each individual.
@@ -22,12 +22,12 @@ vector with the age class at first capture for each individual.
 maximum number of age classes for which m-arrays are constructed; ignored if \code{max(age) > mAge}. Only required if the \code{age} vector has fewer age classes than we want to separate (e.g. capture histories \code{ch} contains only individuals marked as juveniles, and we want 2 age classes).
 }
   \item{freq}{
-a vector with the number of animals with each capture history. If a single value is supplied, it will be used for all rows in the capture history. The default is to assume each row corresponds to a single animal.
+a vector with the number of animals with each capture history, or a matrix with a column for each group. If a single value is supplied, it will be used for all rows in the capture history; the default is to assume each row corresponds to a single animal.
 }
 
 }
 \value{
-A 3-d array, (years-1) x years x age classes, where element [i, j, k] contains the number of individuals of age class k released in year i and recaptured in year j+1 (by definition no recaptures can occur in year 1). The last column contains the number of individuals released in year i and never recaptured.
+A 4-d array, (years-1) x years x age classes x groups, where element [i, j, k, g] contains the number of individuals in group g of age class k released in year i and recaptured in year j+1 (by definition no recaptures can occur in year 1). if no groups are specified, this will be a 3-d array, (years-1) x years x age classes. The last column contains the number of individuals released in year i and never recaptured.
 }
 \author{
 Michael Schaub
diff --git a/man/marrayDead.Rd b/man/marrayDead.Rd
index 0d6d494..d8c5324 100644
--- a/man/marrayDead.Rd
+++ b/man/marrayDead.Rd
@@ -16,11 +16,11 @@ marrayDead(MR, freq = 1)
 an individuals x time matrix with 1 denoting either the time of marking or the time of recovery; otherwise 0.
 }
   \item{freq}{
-a vector with the number of animals with each capture history. If a single value is supplied, it will be used for all rows in the capture history. The default is to assume each row corresponds to a single animal.
+a vector with the number of animals with each capture history, or a matrix with a column for each group. If a single value is supplied, it will be used for all rows in the capture history; the default is to assume each row corresponds to a single animal.
 }
 }
 \value{
-An m-array, a (years-1) x years matrix, where element [i, j] contains the number of individuals marked in year i and recovered in year j+1. The last column contains the number of individuals marked in year i and never recovered.
+An m-array, a (years-1) x years x groups array, where element [i, j, g] contains the number of individuals in group g marked in year i and recovered in year j+1. The last column contains the number of individuals marked in year i and never recovered. If no groups are specified, this will be a (years-1) x years matrix.
 }
 \author{
 Michael Schaub