From Lincoln-Petersen to the Cormack-Jolly-Seber Model

Author

Deon Roos

Published

July 24, 2026

A local hero (two of them, actually)

Before we get into the model, I want to take a moment to point something out that I think is genuinely worth appreciating.

The Cormack-Jolly-Seber model, or CJS model, is one of the most widely used statistical models in ecology. It underpins decades of research into animal survival, population dynamics, and conservation. If you have ever read a paper reporting survival estimates for a wild animal population, there is a very good chance CJS was involved somewhere.

Richard Cormack and George Jolly both developed this model in 1964. Both were working at the University of Aberdeen, in the ARC Unit of Statistics. They worked in the same corridor. They had coffee together most mornings. They even played chess together regularly (specifically a version of chess called kriegspiel which is kind of like battleships + chess).

And neither had any idea the other was working on the same problem.

Cormack later recalled that it was “completely unknown to the two of us that we were working in the same area.” They were sitting across a chess board from each other and somehow the topic of mark-recapture never came up. It was only when Cormack submitted his paper to Biometrika in 1964 that the overlap became apparent, through the referee’s comments rather than through any conversation between the two of them.

It’s also worth noting that Cormack went to Cambridge at 17, did his undergraduate in two years, and became a lecturer in Aberdeen at 21…

George Seber developed his version independently in New Zealand at around the same time, which makes three separate derivations of essentially the same model, across two continents, by people who were unaware of each other’s work. Hence all three names on the tin.

So the next time someone asks you what Aberdeen is known for, you have a genuinely good answer. Two of them, actually, sitting in the same corridor, drinking coffee, playing chess, and accidentally doing the same groundbreaking statistics.

What LP could not do

On the last page we saw that Lincoln-Petersen hands you a snapshot of population size \(N\) from two sampling occasions, as long as you’re willing to assume the population is closed. No births, no deaths, nobody moving in or out.

That closure assumption is doing a lot of heavy lifting. It holds up fine over a short window, then falls apart the moment your study runs over weeks or months. And if the thing you actually care about is survival, then closure is the wrong assumption full stop. You want animals to be able to die between occasions, because dying is the exact thing you’re trying to measure.

CJS (Cormack-Jolly-Seber) throws the closure assumption straight out. Animals can die between occasions. The population is open. And what you get back for giving up closure is a genuine estimate of apparent survival probability \(\phi\).

The catch is that you can no longer get \(N\) directly. CJS conditions on first capture, so it only ever looks at an animal’s history after it’s been caught at least once. Animals you never caught don’t show up in the model at all. That neatly ducks the problem of counting the never-seen, but it also means \(N\) is off the table.

Don’t worry, we get \(N\) back eventually, and clawing it back is one of the big selling points of the robust design. For now though, let’s build CJS up slowly.

The capture history

The core data structure in CJS is the capture history, and it’s genuinely new compared to anything you worked with in BI3010, so it’s worth slowing down for a second.

Each tagged animal gets its own row. Each column is a sampling occasion. A 1 means you detected the animal that occasion, a 0 means you didn’t. That’s the whole format.

Here’s a small example with eight animals surveyed over four occasions:

Code 
library(ggplot2)
library(dplyr)
library(tidyr)

ch <- data.frame(
  animal = paste0("Animal ", 1:8),
  occ1 = c(1, 1, 1, 0, 1, 1, 0, 1),
  occ2 = c(1, 0, 1, 1, 0, 1, 1, 0),
  occ3 = c(0, 1, 1, 1, 1, 0, 1, 1),
  occ4 = c(1, 1, 0, 1, 1, 1, 1, 0)
)

ch
    animal occ1 occ2 occ3 occ4
1 Animal 1    1    1    0    1
2 Animal 2    1    0    1    1
3 Animal 3    1    1    1    0
4 Animal 4    0    1    1    1
5 Animal 5    1    0    1    1
6 Animal 6    1    1    0    1
7 Animal 7    0    1    1    1
8 Animal 8    1    0    1    0

Read each row as a story. Animal 1 was detected on occasions 1, 2, and 4, but not 3. Animal 2 was detected on 1, 3, and 4, but not 2. And so on.

Now, for every 0 in that table, ask yourself the question from the last page: is this animal alive and just missed, or is it dead? A single zero can’t tell you. But the pattern of detections across all the occasions does carry that information, and digging it out is exactly what CJS is for.

Code 
ch_long <- ch |>
  pivot_longer(cols = starts_with("occ"),
               names_to = "occasion",
               values_to = "detected") |>
  mutate(occasion = as.integer(gsub("occ", "", occasion)),
         detected_label = if_else(detected == 1, "Detected", "Not detected"))

ggplot(ch_long, aes(x = occasion, y = animal, fill = detected_label)) +
  geom_tile(colour = "white", linewidth = 1.5) +
  scale_fill_manual(values = c("Detected" = "#d19527",
                                "Not detected" = "#FF5733")) +
  scale_x_continuous(breaks = 1:4, labels = paste("Occasion", 1:4)) +
  labs(x = NULL, y = NULL, fill = NULL) +
  theme_dark_site() +
  theme(legend.position = "bottom",
        panel.grid = element_blank())

The gold tiles are detections. The red tiles are non-detections. For every red tile, there are two possible stories. The CJS model is a machine for working out which story is more likely, given the full pattern of the data.

Have a look at the figure above. Look at Animal 1 at the bottom. We never caught animal 1 on occasion 3 but we did capture it on occasion 4. We’ve said that non-detection can mean two things; 1) alive but not seen or 2) dead. Was animal 1 dead at occasion 3? Hopefully you conclude it was not, but how do you know? That’s the logic that CJS uses.

The two processes

Same as the last page, CJS models two separate processes out in the open.

The first is survival. Between any two consecutive occasions, an animal either survives with probability \(\phi\) (say \(0.9\), or 90%) or dies with probability \(1 - \phi\) (so \(1-0.9=0.1\), or 10%). Once it’s dead, it stays dead. It can’t come back. Zombies don’t exist in this house.

The second is detection. On each occasion, an animal that’s alive gets detected with probability \(p\), or slips past you with probability \(1 - p\).

To actually get a detection in your data, both processes have to go your way: the animal had to survive to that occasion, and you had to spot it. A non-detection is either one of them failing.

It’s the same logic as the last page, just written down now as a model with two named parameters.

Building the likelihood from a capture history

This is where the GLM framing from the first page comes back to bite, in a good way. CJS estimates \(\phi\) and \(p\) by writing down the probability of seeing each animal’s capture history, then hunting for the values of \(\phi\) and \(p\) that make the histories you actually observed as likely as possible. That’s maximum likelihood estimation, which you met in BI3010 (and probably forgot - fair enough).

Let’s walk through one animal’s history to make it concrete.

Take an animal with the history 1 0 1 1. Caught on occasion 1 (that’s how it got its tag), missed on occasion 2, then caught again on 3 and 4.

So what’s the probability of seeing that history? Go occasion by occasion, starting after first capture:

  • Occasion 1 to 2: the animal survived (\(\phi\)) but was not detected (\(1 - p\)). Probability: \(\phi \times (1 - p)\)

  • Occasion 2 to 3: the animal survived (\(\phi\)) and was detected (\(p\)). Probability: \(\phi \times p\)

  • Occasion 3 to 4: the animal survived (\(\phi\)) and was detected (\(p\)). Probability: \(\phi \times p\)

Putting it together, the probability of the full history 1 0 1 1 (after first capture) is:

\[P(\text{1 0 1 1}) = [\phi(1-p)] \times [\phi \cdot p] \times [\phi \cdot p]\]

\[= \phi^3 \times p^2 \times (1-p)\]

Now build your own capture history below and watch the likelihood surface react. For every combination of \(\phi\) and \(p\), the surface shows how probable your chosen history would be if those were the true values. The white dot is the maximum likelihood estimate, the single combination that makes your history most probable of all.

The brighter the patch, the better that combination of \(\phi\) and \(p\) explains the history you built. Maximum likelihood estimation goes looking for that peak across every animal in the dataset at once. Notice how some histories give you a tight, sharp peak; those are the informative ones. Others leave the surface almost flat, which means loads of different \(\phi\) and \(p\) combinations explain the data about equally well, and the model just can’t pin the estimates down with any confidence.

But what about the final non-detection?

There’s a subtlety I’ve been quietly skating over. What if an animal’s last detection isn’t on the last occasion? Take an animal with the history 1 1 0 0.

After that last detection on occasion 2, two different things could be going on. It might have died somewhere around occasion 3. Or it might have sailed right through to the end of the study and simply never been caught again.

The probability of never seeing an animal again from occasion \(t\) onward, usually written \(\chi_t\) (chi, said “ky”), is what handles this. In our four-occasion example \(\chi_4 = 1\), because once you’re past the final occasion there are no more chances to see anything. Then work backwards:

\[\chi_3 = (1 - \phi) + \phi(1-p)\chi_4\]

In words: the chance of never seeing the animal again from occasion 3 is either that it died between 3 and 4 (\((1 - \phi)\)), or that it survived but you missed it on occasion 4 and then never saw it again (\(\phi(1-p)\chi_4\)).

That recursion keeps rolling back through all the occasions. CJS does this for you automatically, but it’s worth knowing it’s under the hood, because it’s what lets the model be honest about the animals that vanish from the data after their first capture.

The model equations

Putting it all together, the CJS model can be written as two linked processes:

\[z_{i,t} \sim Bernoulli(\phi_{t-1} \times z_{i,t-1})\]

\[y_{i,t} \sim Bernoulli(p_t \times z_{i,t})\]

where:

  • \(z_{i,t}\) is the true alive/dead state of individual \(i\) at occasion \(t\). This is a latent variable: we never observe it directly, we only observe detections.

  • \(\phi_{t-1}\) is the probability of surviving from occasion \(t-1\) to \(t\)

  • \(y_{i,t}\) is the observed detection of individual \(i\) at occasion \(t\) (the 0 or 1 in the capture history)

  • \(p_t\) is the probability of detecting individual \(i\) at occasion \(t\), given it is alive

The first equation says the true state at time \(t\) is a coin flip (a Bernoulli trial) that hangs on whether the animal was alive at \(t-1\) and on the survival probability. If it was already dead (\(z_{i,t-1} = 0\)) it stays dead, because \(\phi \times 0 = 0\). If it was alive, it makes it to \(t\) with probability \(\phi\).

The second equation says what you actually see is another coin flip, this one hanging on whether the animal is alive at \(t\) and on the detection probability. You can’t detect a dead animal, since \(p \times 0 = 0\). You detect a live one with probability \(p\).

If those two equations feel like déjà vu, good, they should. It’s the exact thing we said in plain words a moment ago, that seeing an animal needs it to have both survived and been detected. All we’ve done is stop waving our hands and write it down properly.

Fitting a CJS model in R

Let’s simulate some capture history data and fit a CJS model to it, so you can see the whole thing actually work.

Code 
set.seed(42)

n_animals <- 100
n_occasions <- 5
phi_true <- 0.80
p_true   <- 0.45

# Simulate true survival states
alive <- matrix(0, nrow = n_animals, ncol = n_occasions)
alive[, 1] <- 1  # All animals alive at first capture

for (t in 2:n_occasions) {
  alive[, t] <- rbinom(n_animals, 1, alive[, t-1] * phi_true)
}

# Simulate detections
detected <- matrix(0, nrow = n_animals, ncol = n_occasions)
for (t in 1:n_occasions) {
  detected[, t] <- rbinom(n_animals, 1, alive[, t] * p_true)
}

# Force first occasion to 1 (all animals were caught at least once to be tagged)
detected[, 1] <- 1

# Show first 10 capture histories
head(detected, 10)
      [,1] [,2] [,3] [,4] [,5]
 [1,]    1    0    0    0    0
 [2,]    1    0    0    0    0
 [3,]    1    1    0    1    0
 [4,]    1    0    0    0    0
 [5,]    1    0    0    0    1
 [6,]    1    1    1    1    1
 [7,]    1    0    0    0    0
 [8,]    1    1    0    0    0
 [9,]    1    0    0    0    0
[10,]    1    1    0    0    0

Now we can fit the CJS model. We will use the marked package, which provides a clean interface for CJS models in R.

Code 
library(marked)

# Format as data frame with capture history string
ch_strings <- apply(detected, 1, paste, collapse = "")
cjs_data <- data.frame(ch = ch_strings)

# Fit a simple CJS model with constant phi and p
cjs_fit <- crm(cjs_data,
               model = "CJS",
               model.parameters = list(
                 Phi = list(formula = ~ 1),
                 p   = list(formula = ~ 1)
               ),
               hessian = TRUE)

 Number of evaluations:  100  -2lnl: 404.1465438
Code 
cjs_fit

crm Model Summary

Npar :  2
-2lnL:  404.1465
AIC  :  408.1465

Beta
                  Estimate        se        lcl          ucl
Phi.(Intercept)  1.2898696 0.2526289  0.7947169  1.785022249
p.(Intercept)   -0.4241131 0.2133404 -0.8422602 -0.005965891
Code 
# Back-transform estimates from logit scale to probability scale
results <- data.frame(
  parameter = c("Survival (phi)", "Detection (p)"),
  true_value = c(phi_true, p_true),
  estimate = plogis(c(cjs_fit$results$beta$Phi,
                      cjs_fit$results$beta$p))
)

results
       parameter true_value  estimate
1 Survival (phi)       0.80 0.7841251
2  Detection (p)       0.45 0.3955330

Those estimates should come out reassuringly close to the true values we baked into the simulation. They won’t be identical, because we’ve only got a finite sample to work from, but they should sit in the right ballpark. All of which is just to show you that CJS actually does what it says on the tin.

What CJS gives us, and what it still cannot

CJS is a genuine step forward from Lincoln-Petersen. We now have:

  • An estimate of apparent survival \(\phi\) that is not contaminated by detection probability
  • An explicit model for both the biological process (survival) and the observation process (detection)
  • The ability to use more than two sampling occasions, which gives us more information and more precise estimates
  • The ability to model \(\phi\) and \(p\) as functions of covariates, just like in a standard GLM

But CJS still leaves a hole. Because it conditions on first capture and only ever deals with marked animals, it can’t estimate \(N\). We gave that up the moment we walked away from closure.

There’s also a wrinkle in what \(\phi\) really means. CJS is built from animals you caught inside your study area and then either recaptured or didn’t. An animal that permanently clears off looks exactly like a dead one in the capture history. So \(\phi\) is apparent survival: the probability of surviving and sticking around in your study area. True survival is higher than that, by some unknown amount that depends on how much permanent emigration is going on.

So here’s where we’ve landed. CJS gives us survival but can’t count the population. Lincoln-Petersen could count the population, but only by freezing it shut, and it told us nothing about survival. What we actually want is both numbers at once, pulled from a single dataset, and so far neither model can manage that.

That’s exactly what the robust design pulls off. And the trick it uses is being clever about the timing of your sampling, nesting closed sampling occasions inside an open population framework.

That’s the next page.