lecture 17: multistate models

class: center, middle, inverse, title-slide

.title[
# lecture 17: multistate models
]
.subtitle[
## FANR 8370
]
.author[
### <br/><br/><br/>Spring 2025
]

---

## Readings

> ### Kéry & Schaub 264-313

---
# what is a "state"

#### Alive/dead
- CJS model

--
#### Occupied/unoccupied
- Occupancy model

--
#### Geographic location (location `$1, 2, ...,S$`)

--
#### Un-infected/infected/dead

--
#### Age class

--
#### Any combination of above states*

.footnote[
*At least in theory
]

---
# two states

#### State model

`$$\Large z_{i,t} \sim Bernoulli(z_{i,t-1}\phi)$$`

---
# two states

#### Observation model

<br/>

`$$\Large y_{i,t} \sim Bernoulli(z_{i,t}p)$$`

---
# three states

#### State model

---
# three states

#### State model

---
# three states - example

#### Movement between two sites

---
# three states - example

#### Movement between two sites

- `$\large \phi_A$`: Probability alive at site A

- `$\large \phi_B$`: Probability alive at site B

- `$\large \psi_{AB}$`: Probability of moving from site A to site B

- `$\large \psi_{BA}$`: Probability of moving from site B to site A

---
# three states - example

#### State model - transition matrix

.left-column[
<br/>
<br/>
<br/>
#### True state at time *t*
]

.right-column[
<center>
#### True state at time *t + 1*
</center>
<table class="table table-striped table-hover table-condensed table-responsive" style="font-size: 14px; margin-left: auto; margin-right: auto;">
 <thead>
  <tr>
   <th style="text-align:center;">  </th>
   <th style="text-align:center;"> Site A </th>
   <th style="text-align:center;"> Site B </th>
   <th style="text-align:center;"> Dead </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:center;"> Site A </td>
   <td style="text-align:center;"> $\phi_A (1 - \psi_{AB})$ </td>
   <td style="text-align:center;"> $\phi_A \psi_{AB}$ </td>
   <td style="text-align:center;"> $1 - \phi_A$ </td>
  </tr>
  <tr>
   <td style="text-align:center;"> Site B </td>
   <td style="text-align:center;"> $\phi_B \psi_{BA}$ </td>
   <td style="text-align:center;"> $\phi_B (1 - \psi_{BA})$ </td>
   <td style="text-align:center;"> $1 - \phi_B$ </td>
  </tr>
  <tr>
   <td style="text-align:center;"> Dead </td>
   <td style="text-align:center;"> 0 </td>
   <td style="text-align:center;"> 0 </td>
   <td style="text-align:center;"> 1 </td>
  </tr>
</tbody>
</table>

]

---
# three states - example

#### Observation model - observation matrix

.left-column[
<br/>
<br/>
<br/>
#### True state at time *t*
]

.right-column[
<center>
#### Observation at time *t*
</center>
<table class="table table-striped table-hover table-condensed table-responsive" style="font-size: 14px; margin-left: auto; margin-right: auto;">
 <thead>
  <tr>
   <th style="text-align:center;">  </th>
   <th style="text-align:center;"> Site A </th>
   <th style="text-align:center;"> Site B </th>
   <th style="text-align:center;"> Not Seen </th>
  </tr>
 </thead>
<tbody>
  <tr>
   <td style="text-align:center;"> Site A </td>
   <td style="text-align:center;"> $p_A$ </td>
   <td style="text-align:center;"> 0 </td>
   <td style="text-align:center;"> $1 - p_A$ </td>
  </tr>
  <tr>
   <td style="text-align:center;"> Site B </td>
   <td style="text-align:center;"> 0 </td>
   <td style="text-align:center;"> $p_B$ </td>
   <td style="text-align:center;"> $1 - p_B$ </td>
  </tr>
  <tr>
   <td style="text-align:center;"> Dead </td>
   <td style="text-align:center;"> 0 </td>
   <td style="text-align:center;"> 0 </td>
   <td style="text-align:center;"> 1 </td>
  </tr>
</tbody>
</table>

]

---
# three states - example

#### The data

- `$1 =$` Seen alive at site A

- `$2 =$` Seen alive at site B

- `$3 =$` Not seen

``` r
head(rCH)
```

```
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    3    3    3    3    3
## [2,]    1    1    3    3    3    3
## [3,]    1    1    1    3    3    3
## [4,]    1    3    1    1    3    3
## [5,]    1    1    2    3    3    3
## [6,]    1    1    3    3    3    3
```

---
# three states - example

#### The data

- `$f =$` First capture occasion

``` r
tail(rCH)
```

```
##        [,1] [,2] [,3] [,4] [,5] [,6]
## [795,]    3    3    3    3    2    1
## [796,]    3    3    3    3    2    3
## [797,]    3    3    3    3    2    2
## [798,]    3    3    3    3    2    2
## [799,]    3    3    3    3    2    2
## [800,]    3    3    3    3    2    3
```

``` r
tail(f)
```

```
## [1] 5 5 5 5 5 5
```

---
# three states - example

#### The likelihood
- The Bernoulli distribution describes the probability of one (of two) outcomes

- The categorical distribution describes the probability of one (of `$N$`) outcomes
    + `$\Omega$` and `$\Theta$` are matrices defining the state/observation probabilities

--
`$$\large z_{i, t+1}|z_{i,t} \sim categorical(\Omega_{z_{i,t},1:S})$$`

`$$\large y_{i, t}|z_{i,t} \sim categorical(\Theta_{z_{i,t},1:S})$$`

---
# analysis in nimble - priors

``` r
ms <- nimbleCode({
  # Priors 
   phiA ~ dbeta(1, 1)
   phiB ~ dbeta(1, 1)
   psiAB ~ dbeta(1, 1)
   psiBA ~ dbeta(1, 1)
   pA ~ dbeta(1, 1)
   pB ~ dbeta(1, 1)
```

---
# analysis in nimble - transition matrix

``` r
  # Define state-transition and observation matrices
   # Define probabilities of state S(t+1) given S(t)
      ps[1,1] <- phiA * (1-psiAB)
      ps[1,2] <- phiA * psiAB
      ps[1,3] <- 1-phiA
      ps[2,1] <- phiB * psiBA
      ps[2,2] <- phiB * (1-psiBA)
      ps[2,3] <- 1-phiB
      ps[3,1] <- 0
      ps[3,2] <- 0
      ps[3,3] <- 1
```

---
# analysis in nimble - observation matrix

``` r
      # Define probabilities of O(t) given S(t)
      po[1,1] <- pA
      po[1,2] <- 0
      po[1,3] <- 1-pA
      po[2,1] <- 0
      po[2,2] <- pB
      po[2,3] <- 1-pB
      po[3,1] <- 0
      po[3,2] <- 0
      po[3,3] <- 1
```

---
# analysis in nimble - likelihood

``` r
# Likelihood 
for (i in 1:nInd){
   # Define latent state at first capture
   z[i,f[i]] <- y[i,f[i]]
   for (t in (f[i]+1):nOcc){
      # State process: draw S(t) given S(t-1)
      z[i,t] ~ dcat(ps[z[i,t-1],1:nStates])
      # Observation process: draw O(t) given S(t)
      y[i,t] ~ dcat(po[z[i,t],1:nStates])
      } #t
   } #i
})
```

---
# Multistate models

#### Very flexible
- Most capture-recapture models can be formulated as multi-state models (we will explore this more in the next lecture and lab)

- Can accommodate many states

- Can include co-variates and individual/temporal variation

#### But complexity comes at a cost
- Data hungry

- Often not fully identifiable

- Computationally intensive

---
# Example - Montana Toads (Our Favorite!)

Consider a data set we've used before - toads in Montana being tested for chytrid.

``` r
data('WyomingFrogs')
head(frog_caps, n = 3)
```

```
## # A tibble: 3 × 6
##   Ind.ID  Survey.Date Project Species  Sex   Bd.presence
##   <chr>   <chr>       <chr>   <chr>    <chr> <chr>      
## 1 3265127 5/12/16     WY      A.boreas male  <NA>       
## 2 3265127 5/17/15     WY      A.boreas male  negative   
## 3 3265127 5/19/14     WY      A.boreas male  negative
```

If we're interested both in survival AND infection probability, a multistate model is one option for analyzing the data.

---
# Example - Montana Toads (Our Favorite!)

Imagine we want to know both survival and infection probability by sex. This means we have 5 total states: Neg M, Neg F, Pos M, Pos F, and Dead. First let's just look at the true state probabilities. We will define the states in the order listed above.

``` r
   # Define probabilities of state S(t+1) given S(t)
      ps[1,1] <- phiNM * (1-psiNM) #prob Neg M to Neg M
      ps[1,2] <- 0 #prob Neg M to Neg F
      ps[1,3] <- phiNM * psiNM #prob Neg M to Pos M
      ps[1,4] <- 0 #prob Neg M to Pos F
      ps[1,5] <- 1-phiNM #prob Neg M to death
      ps[2,1] <- 0 #prob Neg F to Neg M
      ps[2,2] <- phiNF * (1-psiNF) #prob Neg F to Neg F
      ps[2,3] <- 0 #prob Neg F to Pos M
      ps[2,4] <- phiNF * psiNF #prob Neg F to Pos F
      ps[2,5] <- 1-phiNF #prob Neg F to death
      ...
      etc. 
```

---
# Example - Montana Toads (Our Favorite!)

The next part of the model would need to describe our observation process. Even if we test every animal we capture, there's always a possibility of false negatives or false positives (tests can be wrong!). We can also just miss the animal. Let's assume we're very good at sexing animals and we capture males and females with equal probability.

.smaller[

``` r
 # Define probabilities of O(t) given S(t)
# Account for both capture probability and testing error
      po[1,1] <- pCap * (1-pFalseP) #observed as Neg M given it was a Neg M
      po[1,2] <- 0 #observed as Neg M given it was a Neg F. 
      po[1,3] <- pCap * pFalseP #observed as Neg M given it was actually positive
      po[1,4] <- 0 #observed as Neg M given it was a Pos F. 
      po[1,5] <- 0 #observed as Neg M given it was dead
      po[2,1] <- 0 #observed Neg F given it was a Neg M
      po[2,2] <- pCap * (1-pFalseP) #observed Neg F given it was a Neg F
      ...
      ...
      po[5,1] <- 1-pCap #not seen given it was Neg M 
      po[5,2] <- 1-pCap #not seen given it was Neg F
      po[5,3] <- 1-pCap #not seen given it was Pos M
      po[5,4] <- 1-pCap #not seen given it was Pos F
      po[5,5] <- 1 #not seen given it was dead
      ... 
      etc
```
]