Here, we will go over the construction of m-arrays from individual capture histories. Schaub & Kéry (2022) Integrated Population Models Chapter 5 : Introduction to integrated population models

We will use the code here and on Mike Meredith’s Github page:

https://github.com/mikemeredith/IPM_code/blob/main/IPM_05/IPM_05.2.R

First load the required package and data.

library(IPMbook)
data(woodchat5)
str(woodchat5)

## List of 4
##  $ ch   : num [1:1902, 1:20] 1 1 1 1 1 1 1 1 1 1 ...
##  $ age  : num [1:1902] 2 2 2 2 2 2 2 2 2 2 ...
##  $ repro: num [1:929, 1:3] 6 2 2 5 3 5 3 2 3 2 ...
##   ..- attr(*, "dimnames")=List of 2
##   .. ..$ : NULL
##   .. ..$ : chr [1:3] "Reproduction" "Year" "Age of mother"
##  $ count: num [1:20] 91 119 131 88 139 145 148 116 112 106 ...

Here we see that the woodchat data contains data collected over 20 years:

ch: capture histories for 1902 females x 20 years.
age: the age of each individual at first capture: 1 = juvenile, 2 = at least 1 year old (adult).
repro: reproduction information for 929 monitored nests:
- Reproduction = total number of fledglings (male + female),
- Year = year monitored,
- Age of mother = 1 = 1 year old, 2 = at least 2 years old.
count: total number of breeding pairs observed in each year.

For many models, the m-array is often used as a way to repackage monitoring data to facilitate many types of survival models.

marr <- marray(woodchat5$ch)
marr[1:5,1:5] # inspect the first few bits of data

##         recaptured
## released Y2 Y3 Y4 Y5 Y6
##       Y1 24  4  0  0  0
##       Y2  0 20  7  1  0
##       Y3  0  0 26  4  0
##       Y4  0  0  0 31  2
##       Y5  0  0  0  0 35

M-arrays look very different from the raw capture histories they are built from. In the capture histories, columns are years of the study, rows are individual birds, cell values are detections (1) or no detection (0) including the initial capture and release of the bird.

(woodchat5$ch[1:10,]) # Capture history of first 10 birds in the data set

##       [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13]
##  [1,]    1    0    0    0    0    0    0    0    0     0     0     0     0
##  [2,]    1    0    0    0    0    0    0    0    0     0     0     0     0
##  [3,]    1    1    1    0    0    0    0    0    0     0     0     0     0
##  [4,]    1    1    0    0    0    0    0    0    0     0     0     0     0
##  [5,]    1    0    1    1    0    0    0    0    0     0     0     0     0
##  [6,]    1    0    0    0    0    0    0    0    0     0     0     0     0
##  [7,]    1    0    0    0    0    0    0    0    0     0     0     0     0
##  [8,]    1    1    0    0    0    0    0    0    0     0     0     0     0
##  [9,]    1    1    0    0    0    0    0    0    0     0     0     0     0
## [10,]    1    0    0    0    0    0    0    0    0     0     0     0     0
##       [,14] [,15] [,16] [,17] [,18] [,19] [,20]
##  [1,]     0     0     0     0     0     0     0
##  [2,]     0     0     0     0     0     0     0
##  [3,]     0     0     0     0     0     0     0
##  [4,]     0     0     0     0     0     0     0
##  [5,]     0     0     0     0     0     0     0
##  [6,]     0     0     0     0     0     0     0
##  [7,]     0     0     0     0     0     0     0
##  [8,]     0     0     0     0     0     0     0
##  [9,]     0     0     0     0     0     0     0
## [10,]     0     0     0     0     0     0     0

The rows of the m-array, however, are annual summaries of recaptures of woodchats, where a recapture occurs 1 or more years after a release. The columns are the year of each recapture for woodchats released each year (each row). The releases span the years from 1 to 19, while the recaptures span years 2 to 20 because the first year cannot include any recaptures and the last year is not a release (even if birds are released in the last year, there are no follow-up observations in year 21). However there are 20 columns in the m-array because the last column is the number of releases in each year that are never seen again, the ‘never’ column. In a study that runs long after the last marking (long enough for every marked individual to die), the sum of the ‘never’ column will equal the total number of individuals in the study.

The definitions of release and recapture are specific to the m-array. A release in year t is any bird held in the hand that year - whether the capture and marking of a new bird or the capture and release of a previously marked bird. A recapture is the most recent detection of a bird a year or more after its release in year t. Every recapture of a bird in year t also counts as a release for that year, except in the last year of the study, year 20. For a given row, t, in the m-array, the columns contain counts of the next subsequently observation of (in all years >t) of each newly marked and recaptured bird released in year t.

To get these counts from the capture histories, the m-array slices all individual capture histories into unique fragments where each fragment contains a release (initial marking or a recapture) and the next subsequent recapture.

How many such ‘fragments’ are there for individual 20?

woodchat5$ch[20,]

##  [1] 1 1 0 1 0 0 0 1 1 0 0 1 0 1 1 0 0 0 0 0

Individual 20 was initially marked in year 1 and then recaptured and released again in years 2, 4, 8, 9, 12, 14, and 15 and then not seen again before the end of the study

There are 8 fragments contributed by individual 20 to the m-array. Note this the same number of 1’s in the capture history for individual 20.

Unique fragments of individual 20’s captured history that contribute to the counts in the m-array:

[1,1,-,-,-,-,-,-,-,-,-,-,-,-,-,-,-,-,-,-] initially marked and recaptured in 2nd year of study
[-,1,0,1,-,-,-,-,-,-,-,-,-,-,-,-,-,-,-,-] released again in year 2 and then recaptured again in year 4
[-,-,-,1,0,0,0,1,-,-,-,-,-,-,-,-,-,-,-,-] released again in year in year 4 and then recaptured again in year 8
[-,-,-,-,-,-,-,1,1,-,-,-,-,-,-,-,-,-,-,-] released again in year 8 and then recaptured again in year 9
[-,-,-,-,-,-,-,-,1,0,0,1,-,-,-,-,-,-,-,-] released again in year 9 and then recaptured again in year 12
[-,-,-,-,-,-,-,-,-,-,-,1,0,1,-,-,-,-,-,-] released again in year 12 and then recaptured again in year 14
[-,-,-,-,-,-,-,-,-,-,-,-,-,1,1,-,-,-,-,-] released again in year 14 and then recaptured again in year 15
[-,-,-,-,-,-,-,-,-,-,-,-,-,-,1,0,-,-,-,-] released again in year 15 and then never seen again

Notice how every recapture also counts as a release for this bird. Individuals that are only marked and never seen again contributed only 1 observation to the m-array, like individal 10 :

woodchat5$ch[10,] # marked in year 1, never seen again

##  [1] 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

The counts in the m-array are observations of unique fragments, i.e. counts of consecutive resightings of individuals organized by year of release and the year of next recapture.

Consider the first two ‘fragments’ of individual 20’s capture history. The first fragment of individual 20 (release in year 1, recaptured year 2), contributes 1 observation to the first cell of the m-array. The second fragment (release in year 2, recaptured year 4) contributes 1 observation to the third column of the second row in the m-array. Below are these two contributions to the first 3 columns of the m-array by individual 20.

	recaptured year 2	recaptured year 4
year 1 release	1
year 2 release		1
year 3 release
…

Any other birds with the exact same fragments will also contribute to these counts.

For survival analysis, only intervals in detection histories matter.

Every detection of an individual in the capture histories contributes one observation to the m-array (even if only to the ‘never’ column) except for detections in the last year (T) of the study. Detections of previously marked individuals in the last year of the study (recaptures in the last year) can only inform the recapture fates in previous years.

You can see that the total number of detections in the individual capture histories in the first 19 years (all but the last year) is equal to the total number of observations in the m-array.

cat("total detections in first 19 years of a 20 year study:", sum(woodchat5$ch[,1:19]))

## total detections in first 19 years of a 20 year study: 2543

cat("total observations in the marray:",sum(marr))

## total observations in the marray: 2543

The rows of the m-array provide us with a tally of the fate of the woodchats known to be alive at the beginning of each time step (row number) showing us exactly how many of those woodchats showed up in later years or were never seen again (column number) - it is a tally of intervals of confirmed survival over a known number of years (and also a tally of permanent disappearance from the study, in the case of the ‘never’ column, over an interval of confirmed non-detection).

How do we build the tally of fates for each cell of the m-array? We start by calculating the number of releases in each year. A release is marked as a 1 in the capture histories of each individual for that year (except the last year). It can denote either a newly marked bird (getting marked and released counts as a detection) and successful recaptures of birds that were marked in previous years. We do not distinguish a difference here.

total_released<-colSums(woodchat5$ch)
# confirm this matches the row sums of the m-array for all but the last year
cbind(total_released,c(rowSums(marr),NA))

##     total_released    
## Y1              87  87
## Y2              92  92
## Y3              99  99
## Y4             126 126
## Y5             134 134
## Y6             116 116
## Y7             104 104
## Y8             137 137
## Y9             116 116
## Y10            128 128
## Y11            135 135
## Y12            156 156
## Y13            137 137
## Y14            154 154
## Y15            167 167
## Y16            155 155
## Y17            184 184
## Y18            166 166
## Y19            150 150
##                151  NA

In some years, there were as many as 184 woodchats in-hand (newly marked or recaptured).

Now, for each of the woodchats released in each of the years from year 1 to 19, we wish to tally the first time each of those woodchats were next seen in any of the years between year 2 to 20. For a given release in year t, we are only concered with the next recapture of a bird, not all recaptures of the bird in 20 years. And obviously, a woodchat released in year t can only contribute to recaptures in years greater than t.

There are custom functions in the IPMbook library that can create such and m-array from capture history data (marray, marrayAge) but it is difficult to see what these do. The function below is hopefully a little more transparent.

mymarray<-matrix(0,ncol=ncol(woodchat5$ch),nrow=ncol(woodchat5$ch)-1) # create an empty m-array full of 0's - a container to hold our data
for(i in 1:(ncol(woodchat5$ch)-1)){ # loop through 19 times, year index
  releases<-which(woodchat5$ch[,i]==1) # identify all woodchats in the year i with a detection, each one will count a release that can be next recaptures in later years or never seen again.
  years_to_next_recapture<-table(apply(matrix(woodchat5$ch[releases,(i+1):20],ncol=20-i),1,function(row) match(1,row)),useNA=("always")) # for each woodchat releaseed in the current year i, find the next column in its capture history (years into the future)  that also has a 1 - this is the next recapture.  'table' tallies recaptures 1 year into the future, 2 years into the future,etc..., as well as those that never get recaptured which are NA's in this output table.
  mymarray[i,c((i+as.numeric(na.omit(names(years_to_next_recapture))))-1,20)]<-as.numeric(years_to_next_recapture) # insert the recapture and 'never' counts for year i into the appropriate row and columns of the m-array
}

#compare to the m-array from the IPMbook function
mymarray==marr

##         recaptured
## released   Y2   Y3   Y4   Y5   Y6   Y7   Y8   Y9  Y10  Y11  Y12  Y13  Y14  Y15
##      Y1  TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y2  TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y3  TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y4  TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y5  TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y6  TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y7  TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y8  TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y9  TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y10 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y11 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y12 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y13 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y14 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y15 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y16 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y17 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y18 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##      Y19 TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
##         recaptured
## released  Y16  Y17  Y18  Y19  Y20 never
##      Y1  TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y2  TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y3  TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y4  TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y5  TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y6  TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y7  TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y8  TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y9  TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y10 TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y11 TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y12 TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y13 TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y14 TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y15 TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y16 TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y17 TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y18 TRUE TRUE TRUE TRUE TRUE  TRUE
##      Y19 TRUE TRUE TRUE TRUE TRUE  TRUE

Working through one year at a time, this code finds each woodchat with a capture history of ‘1’ in year t. These rows in the capture history make up the total number of ‘releases’ for year t. For each of these woodchats, the code then finds which column is the next column after year t (going to the right in the capture history, forward in time) that also has a 1 - the next recapture. It then tallies the number of released woodchats in each year that were subsequently recaptured 1 year after release, 2 years after release, and so on. If a 1 is never found in any of the columns after year t, it tallies these as NA’s which are the ‘never’ recaptured birds of the released woodchats that year.

Detection probability and survival

Perhaps it is obvious that in order for an individual to leave more than two 1’s in its capture history, it must survive at least 1 year after recapture. We can see that each row of the m-array contains a set of observations of individuals that we know survived over intervales of one year or more. However, a majority of observations are in the ‘never’ column each year, meaning that most releases are never seen again. Eventually, all woodchats in the data set must end up contributing one observation to ‘never’ (if marking stops but monitoring continues). For some woodchats, a ‘never’ is the only observation they contribute to the m-array.

We might be tempted to take each row of data in the m-array and estimate survival probability as the ratio of the number of recaptures divided by the total number of released woodchats.

sum(marr[1,1:19])/sum(marr[1,1:20]) # ad hoc estimate of annual survival in year 1 based on the proportion of released animals in the first year that were successfully recaptured a year or more later

## [1] 0.3218391

This assumes the woodchats in the ‘never’ recaptured column actually died, and that the probability that any of the woodchats released in year 1 are recaptured \((\pi)\) is simply:

\(\pi = S = annual\ survival\ probability\)

The ‘never’ column only indicates disappearance and disappearance could be due to death or due to a failure to find and capture the individual, which is a probability \((q)\). The complement of this failure \((1-q)\) is our success at recapture or ‘capture probability’ denoted \(p\). I will probably use the more general term, ‘detection probability’, but it is the same. It is the probability that an individual in the population is successfully captured or recaptured by our methods. The probability that any of the released woodchats from year t actually lead to a recapture in year year t+1 is actually:

\(S * p\)

where

\(p = detection\ probability\)

What about the woodchats that were released in year 1 that were next detected in year 3? These woodchats have to roll several dice, first they have survive from year 1 to year 2, then they have to go undetected in year 2, then they have to survive again from year 2 to year 3, and then they have to be detected in year 3. The probability of this happening is:

\(S * (1-p) * S * p\)

The first two terms above apply to all the woodchats that went undetected in year 2 but were alive. All woodchats have to survive before they can be detected or be ‘missed’ by the observer. Here, i use the term ‘missed’ differently from ‘undetected’ to distinguish between an animal that is alive but unseen from a 0 in the capture history which could arise from being ‘missed’ or because an individual died. When a 0 is recorded in the capture history, we do not know whether it is due to death or capture failure. IPM’s can estimate the probability of each of these, however. Remeber that detection probability only applies to those woodchats that are alive.

At any point in time, each released woodchat has three immediate options for contributing to the m-array. They can die at rate \((1-S)\) and contribute their lives to the ‘never’ column. They can live but be missed by the observer at rate \(S*(1-p)\) and eventually contribute to the ‘never’ column or one of the other recapture columns in future years. Finally, they can live and be detected at rate \(S*p\) and contribute to the recapture column for year following release.

If an individual dies, it has run out of options. If it survives a year, it faces all three options again, the next year.

Note that:

\([1-S] + [S*(1-p)] + [S*p] = 1\)

These three outcomes (die, missed, detected) soak up all the probability for what can happen in one year and there is no other option available to a woodchat, as far as we’re concerned here (emigration is an option in some systems), One of these options produces a 1 in the capture histories, the other two both produce a 0.

What about the woodchats that are released but are not detected by the end of the study? Can you derive their probabilities? Let take a woodchat marked and released in year 1 and then never detected again. One possibility is that it died, never making it to year 2, at probability

\(1-S\)

Another possibility is that it lived one year (to year 2), went undetected, and then died the following year (before year 3):

\(S*(1-p) * (1-S)\)

Yet another possibility is that is lived 5 years (to year 6), went undetected after each of those years and then died before the 7th year

\([S * (1-p)]^5 *(1-S)\)

And perhaps by some alignment of the stars, despite getting captured in the first year of the study, a woodchat lived 19 years successfully avoiding recapture every year. This is a lot of coin tossing to always come up heads (1) on survival but tails (0) on detection every year.

\([S*(1-p)]^{19}\)

The probability that all the woodchats in the ‘never’ column of our m-array being like this last case would require very high survival and very low detection probability, which might be the case for some systems.

The point is, every number in the m-array is assumed to be derived from known combinations of the probabilities of survival and detection - quantities we don’t know but would like to know - our latent parameters. By specifying a model which includes equations for each cell of the m-array, we may be able to derive an estimate of the annual survival and the detection probability.

Lets examine the code from chapter 5 section 2 which is a very simple population model using the capture history data. Before we look at this code, we need a special type of m-array that separates woodchats by their age at first capture. This simply creates to m-array matrices - one built using only the capture histories from woodchats first captured as juveniles (almost 1 year old) and one using only the capture histories from woodchats first captured as adults.

marr<-marrayAge(woodchat5$ch,woodchat5$age)
# Bundle data for JAGS
jags.data <- list(marr.j=marr[,,1], marr.a=marr[,,2], n.occasions=dim(marr)[2],
    rel.j=rowSums(marr[,,1]), rel.a=rowSums(marr[,,2]), mean.f=c(2.6, 3.6), T=15)

Looking at the code from:https://github.com/mikemeredith/IPM_code/blob/main/IPM_05/IPM_05.2.R

Note the lines in the JAGS model that first call the actual m-array constructed from observations of wood-chats in the field. These are lines 54 and 55.

for (t in 1:(n.occasions-1)){ 
  marr.j[t,1:n.occasions] ~ dmulti(pr.j[t,], rel.j[t]) 
  marr.a[t,1:n.occasions] ~ dmulti(pr.a[t,], rel.a[t]) 
}

This IS THE ONLY PLACE WHERE OBSERVATIONS OF WOODCHATS ARE IN THIS MODEL. This chunk of code looks at an m-array for juveniles and adults as separate data sets (for now) as though they came from different populatoins. It simply states that the observed counts in each row an m-array are proportional to the cell probabilities of a successful recapture \((\pi)\), stored in matrices, pr.j[t,] and pr,a[t,], after release in year t. The rel.j[t] and rel.a[t] terms are simply the total number of woodchat juveniles and adults released in each year.

Where do the probabilities controlling the expected number of recaptures in each year come from?

pr.j[t,] # a vector of probabilities for woodchat juveniles released in year t, for the likelihood of being recaptured in each year of the study

pr.a[t,] # same but for adults

Both pr.j and pr.a are matrices will cells that directly correspond to the cells of the m-array. The m-array is filled with counts while the pr.j and pr.a matrices are filled with cell probabilities \((\pi)\) - the combined probabilities of survival and detection that lead to recaptures - each row holding the probabilities that a woodchat released in year t will contribute to each recapture counts held in row t of the m-array.

In lines 57 through 78, we actually construct these cell probabilities using placeholder parameters for survival and detection. For simplicity, I will first examine only the lines of this code referencing the adult probabilities. the code below is confusing because it indexes the m-array columns using a number, rather than the names of the columns, and it conveniently works out to index the columns simply using t - the index used for year. However, remember that the column numbers of the m-array actually reference recaptures in year t+1. Thes any reference to recapture indexed by t in fact refers to time t+1. To try to draw a distinction between year and recapture, I specify the recapture number as a separate placeholder.

  for (t in 1:19){ # loop over 19 times
    t_rec<-t # a nominal distinction between the years of the study, t, and the years of recapture, t_rec. The first year of recaptures is actually the second year of the study. The actual point in time of any particular t_rec is trec+1.
    q[t_rec] <- 1 - p[t_rec] # Probability of non-recapture in the following year, if still alive at time t+1 
    # Cell Probability that a woodchat released in year t shows up the next year as a recapture - these counts are the main diagonal of the m-array 
    pr.a[t,t_rec] <- sa[t] * p[t_rec].
    # Next,cell probabilities determining counts above main diagonal in the m-array these are the probabilities of next recapture starting with year t+2 until the last year of the study. All of these probabilities require the woodchat to have been missed at least once before being recaptured.
    for (j in (t_rec+1):19){ # nested j loop, for each year of release (t loop above) starting with the second year of possible recapture (e.g. year 3 after release in year 1) and going until the last possible recapture (recap session 19 in year 20).
      pr.a[t,j] <-  # getting counted as a recapture from release t in recapture session j is the product of:
        prod(sa[t:j]) * # prob of survival over the years up to recapture j starting from year of release 
        prod(q[t_rec:(j-1)]) * # prob of being missed every recapture session up until the year before recapture session j  (year j+1) AND
        p[j] # prob of being detected in recapture session j (year j+1)
    } # end j loop
    # Next, cell probabilites below main diagonal. Have to be 0 because a woodchat cannot be released in year t and recaptured in acyear before t.
    for (j in 1:(t-1)){
      pr.a[t,j] <- 0
    } #j
  } #t
  # All we need is the last column: probability of never beeing seen again
  for (t in 1:19){
    pr.a[t,20] <- 1-sum(pr.j[t,1:19]) # 1 - sum of all the other cell probabilities in the row - it's the only option left for a wood chat.
  }

This chunk of code is a model - a model of the cell probabilities that influence the counts in the m-array. Adding in the code for the cell probabilities of juveniles is very easy, the lines are identical to the adults except for the upper diagonals:

  for (t in 1:19){
    t_rec<-t 
    q[t_rec] <- 1 - p[t_rec] 
    #main diagonal
    pr.j[t,t_rec] <- sj[t] * p[t_rec]
    pr.a[t,t_rec] <- sa[t] * p[t_trec] 
    #above main diagonal
    for (j in (t_rec+1):19){ 
      pr.j[t,j] <- sj[t] * prod(sa[(t_rec+1):j]) * prod(q[t_rec:(j-1)]) * p[j] # this is the line that is different for juveniles. the only difference is the survival probabilities, the first year of survival after release is based on juvenile survival which is often lower than adult survival. After that first survival interval, survival is then driven by the adult survival probablity. 
      pr.a[t,j] <- prod(sa[t:j])* prod(q[t_rec:(j-1)])*p[j]
    } # end j loop
    for (j in 1:(t-1)){
      pr.j[t,j] <- 0
      pr.a[t,j] <- 0
    } # end j loop
  } #end t loop
  # last column
  for (t in 1:19){
    pr.j[t,20] <- 1-sum(pr.j[t,1:19])
    pr.a[t,20] <- 1-sum(pr.a[t,1:19]) 
  }

This is now a complete model for the probability that a captured and released woodchat will show up as a recapture (or a ‘never’) in the m-array we built from the data.

Lines 82-91 in 5.2 is another model, a model of population dynamics based on the estimated survival and assumed reproduction. This chunk of code is not necessary for estimating survival, and it is not referencing any observations:

# Model for population dynamics: a simple matrix population model
# Define initial stage-specific population sizes
  N[1,1] <- 1
  N[2,1] <- 1
  # Loop over time
  for (t in 1:T){
  # Population projection
  N[1,t+1] <- N[1,t] * mean.f[1] / 2 * mean.sj + N[2,t] * mean.f[2] / 2 * mean.sj
  N[2,t+1] <- (N[1,t] + N[2,t]) * mean.sa
# Annual growth rate
ann.growth.rate[t] <- (N[1,t+1] + N[2,t+1]) / (N[1,t] + N[2,t])
}
lambda <- ann.growth.rate[T]

You can see this code creates a matrix N, with two columns, one column for the number of juveniles in each year and the second column for the number of adults in each year. The number of juveniles and adults are both initialized at a population size of 1 in the first year (the number doesn’t matter because we are simply using the predicted trend to estiamate population growth rate). The number of juveniles in each subsequent year depends on the number of juveniles AND adults and their respective fecundities, mean.f[1] and mean.f[2] which are provided to the model as constants (2.6 young per female and 3.6 young per female).

However the reproductive output of each juvenile and adult is divided by 2 - so we are assuming a 50:50 sex ratio and modelling only the production of female woodchats. Finally the production of new females by both juvenile and adult woodchats is then multiplied by the estimated survival of juveniles implying all juveniles must survive one year to contribute to the count of juveniles in the next year. This is a variant of a pre-breeding survey - the counts of juveniles refer to individuals that are almost 1-year old. Shortly after being counted, these indivduals will reproduce and mature into adults. They will be considered adults in the next count (if they survive) and there offspring counted in that next count will be considered to be the offspring of juveniles, even though the parents are will now be adults.

In this variant juvenile survival refers to both the survival of hatchlings to ‘juvenilehood’ as well as the survival of juveniles to ‘adulthood’, which may not apply to some systems. It most cases it may be better to model the survival of juveniles to adulthood as part of adult survival like most prebreeding surveys.

Building m-arrays

David Christianson

2025-02-11

For survival analysis, only intervals in detection histories matter.

Detection probability and survival