skagit: An example of the Multistrata (multi-state) model in which...

In RMark: R Code for Mark Analysis

An example of the Multistrata (multi-state) model in which states are routes taken by migrating fish.

Description

An example of the Multistrata (multi-state) model in which states are routes taken by migrating fish.

Author(s)

Megan Moore <megan.moore at noaa.gov>

Examples

# There are just two states which correspond to route A and route B. There are 6 occasions
# which are the locations rather than times. After release at 1=A there is no movement 
# between states for the first segment, fish are migrating downriver together and all pass 2A. 
# Then after occasion 2, migrants go down the North Fork (3A) or the South Fork (3B), 
# which both empty into Skagit Bay. Once in saltwater, they can go north to Deception Pass (4A)
# or South to a receiver array exiting South Skagit Bay (4B). Fish in route A can then only go
# to the Strait of Juan de Fuca, while fish in route B must pass by Admiralty Inlet (5B). 
# Then both routes end with the array at the Strait of Juan de Fuca.
#
#       1A
#        |
#        2A
#      /     \
#    3A        3B
#   /  \      /  \ 
# 4A   4B    4A  4B       
#  |     \    /   |
#   5A    5B  5A   5B
#      \   \   /    /
#            6
# 
# from 3A and 3B they can branch to either 4A or 4B; branches merge at 6  
# 5A does not exist so p=0; only survival from 4A to 6 can be 
# estimated which is done by setting survival from 4A to 5A to 1 and
# estimating survival from 5A to 6 which is then total survival from 4A to 6.

pathtodata=paste(path.package("RMark"),"extdata",sep="/")
skagit=import.chdata(paste(pathtodata,"skagit.txt",sep="/"),field.types=c("f"),header=TRUE)
skagit.processed=process.data(skagit,model="Multistrata",groups=c("tag"))
skagit.ddl=make.design.data(skagit.processed)
#
# p
#
# Can't be seen at 5A or 2B,6B (the latter 2 don't exist)
skagit.ddl$p$fix=ifelse((skagit.ddl$p$stratum=="A"&skagit.ddl$p$time==5) | 
 (skagit.ddl$p$stratum=="B"&skagit.ddl$p$time%in%c(2,6)),0,NA)
# Estimated externally from current data to allow estimation of survival at last interval
skagit.ddl$p$fix[skagit.ddl$p$tag=="v7"&skagit.ddl$p$time==6&skagit.ddl$p$stratum=="A"]=0.687
skagit.ddl$p$fix[skagit.ddl$p$tag=="v9"&skagit.ddl$p$time==6&skagit.ddl$p$stratum=="A"]=0.975
#
# Psi
#
# only 3 possible transitions are A to B at time interval 2 to 3 and 
# for time interval 3 to 4 from A to B and from B to A
# rest are fixed values
skagit.ddl$Psi$fix=NA
# stay in A for intervals 1-2, 4-5 and 5-6
skagit.ddl$Psi$fix[skagit.ddl$Psi$stratum=="A"&
  skagit.ddl$Psi$tostratum=="B"&skagit.ddl$Psi$time%in%c(1,4,5)]=0
# stay in B for interval 4-5
skagit.ddl$Psi$fix[skagit.ddl$Psi$stratum=="B"&skagit.ddl$Psi$tostratum=="A"
  &skagit.ddl$Psi$time==4]=0
# leave B to go to A for interval 5-6
skagit.ddl$Psi$fix[skagit.ddl$Psi$stratum=="B"&skagit.ddl$Psi$tostratum=="A"&
 skagit.ddl$Psi$time==5]=1
# "stay" in B for interval 1-2 and 2-3 because none will be in B
skagit.ddl$Psi$fix[skagit.ddl$Psi$stratum=="B"&skagit.ddl$Psi$tostratum=="A"&
 skagit.ddl$Psi$time%in%1:2]=0
# 
# S
#
# None in B, so fixing S to 1
skagit.ddl$S$fix=ifelse(skagit.ddl$S$stratum=="B"&skagit.ddl$S$time%in%c(1,2),1,NA)
skagit.ddl$S$fix[skagit.ddl$S$stratum=="A"&skagit.ddl$S$time==4]=1
# fit model
p.timexstratum.tag=list(formula=~time:stratum+tag,remove.intercept=TRUE)
Psi.sxtime=list(formula=~-1+stratum:time)
S.stratumxtime=list(formula=~-1+stratum:time)
#
S.timexstratum.p.timexstratum.Psi.sxtime=mark(skagit.processed,skagit.ddl,
 model.parameters=list(S=S.stratumxtime,p= p.timexstratum.tag,Psi=Psi.sxtime),delete=TRUE)
# calculation of cummulative survival for entire route
Sest=plogis(coef(S.timexstratum.p.timexstratum.Psi.sxtime)$estimate)
# A
prod(Sest[c(1:3,6)])
#[1] 0.1644
# B
prod(Sest[c(1,2,4,5,7)])
#[1] 0.1154

RMark documentation built on Aug. 14, 2022, 1:05 a.m.