Applications of hidden hybrid Markov/semi-Markov models: from stopover duration to breeding success dynamics

Usually in capture–recapture, a model parameter is time or time since first capture dependent. However, the case where the probability of staying in one state depends on the time spent in that particular state is not rare. Hidden Markov models are not appropriate to manage these situations. A more convenient approach would be to consider models that incorporate semi-Markovian states which explicitly define the waiting time distribution and have been used in previous biologic studies as a convenient framework for modeling the time spent in a given physiological state. Here, we propose hidden Markovian models that combine several nonhomogeneous Markovian states with one semi-Markovian state and which (i) are well adapted to imperfect and variable detection and (ii) allow us to consider time, time since first capture, and time spent in one state effects. Implementation details depending on the number of semi-Markovian states are discussed. From a user’s perspective, the present approach enhances the toolbox for analyzing capture–recapture data. We then show the potential of this framework by means of two ecological examples: (i) stopover duration and (ii) breeding success dynamics.


Introduction
In the study of population dynamics, the generation (or life) time of species has long been used to classify life-history trait variations. Species may be placed along the r-and K-strategy gradient (Pianka 1970) depending on their reproductive and survival rates (Saether et al. 2004). High reproductive species with short life expectancy are at one end of the continuum (fast), and species that often produce a single offspring but are long lived are at the other (slow). This classification has, for instance, been used to predict how demographic stochasticity affects population dynamics (Saether et al. 2004) or to detect cost of reproduction (Hamel et al. 2010). The time spent by an individual in a particular biologic state (e.g., alive or breeding) is of critical importance for evolutionary ecologists. Beyond the time spent alive or breeding, remaining lifetime is spent in different phases (e.g., dormancy in plants) or in specific behavior (courtship display, resting, hibernating, molting in animals). In general, the probability of staying, and therefore exiting, in a given state depends on the time elapsed in this state. It may present little individual variation due to certain physiological constraints, or because it has been shaped by evolution (e.g., incubation or gestation). In contrast, the time spent in other biologic states may be very variable because of either intrinsic individual variations, or heterogeneous response to density dependence or environmental effects [e.g., courtship, refueling, molting; Dietz et al. (2013)].
Empirical knowledge of some state durations may help in estimating the causes of life-history trait variations. For example in Schmaltz et al. (2011), consideration of incubation duration led to an understanding of variations in greater flamingo breeding performance. After hatching, variations in adult provision rates also affect the duration of fasting periods in the juveniles of several bird species and may help diagnose sex-specific costs of parental care (e.g., Booth et al. 2000;Amat et al. 2007). Tracking temporary and definitive changes in biologic stages is of critical importance in several species. Some species will respond to seasonal variations in environmental conditions by changing physiological state, for example, hibernation (Bryant and Page 2005) or by migrating, while others will permanently change their life-history stage: caterpillars becoming butterflies, tadpoles becoming frogs, etc. Such transitional or definitive stages constitute critical life phases, and their duration is central to the understanding of population dynamics.
A clear understanding of variation in life-history traits has numerous applications. Estimating of the stopover duration in migrating birds or time to extinction in population viability analyses is critical to conservation biologists. Similarly, an understanding of the dynamics of epidemic spread is also dependent on in-depth knowledge of the time spent by individuals in the infected state (Shaman and Karspeck 2012).
The state of an individual is determined from different clues (visual appearance, behavior, measurements, etc.). But, it may be observed incorrectly or remain unobserved. Even in the simplest situation, state duration may be difficult to estimate. For example, knowledge of life duration is conditional on knowledge of times of birth and death. In the wild, time of death is rarely known, but classical CR models can deal with this (Lebreton et al. 1992;Choquet et al. 2011) by estimating a capture probability. Thus, we can deal with the lack of information after the last capture of an individual by stating that after that date, an individual may be dead or missed at capture. But sometimes the age of an individual is imperfectly known, and time of birth should also be modeled (Colchero and Clark 2012;Matechou et al. 2013).
A major drawback of hidden Markov models is the inflexible description of the time spent in a given state, as sojourn time (state occupancy) distributions are implicitly geometric. To overcome this limitation, a semi-Markovian framework may be considered where parametric sojourn time distributions are incorporated in the model, or where states are replaced by series-parallel networks of states with a common observation distribution; see Gu edon (2005), Langrock and Zucchini (2011) and references therein. The sojourn time distributions of the macro-states defined in this way are built from the implicit geometric sojourn time distributions of the elementary Markovian states. These geometric distributions are combined either by convolution for states in series or by mixture for (series of) states in parallel. Gu edon (2005) showed that hidden Markov models with macro-states are not a valid alternative to hidden semi-Markov models because of higher algorithm space complexity and strong constraints in the definition of parametric macro-state sojourn time distributions. It should be noted that Markovian and semi-Markovian states can easily be combined in hidden hybrid models. In particular, Choquet et al. (2013b), unified the models proposed by Pledger et al. (2009), Pradel (2009 and Fewster and Patenaude (2009) and proposed a new model combining nonhomogeneous Markovian states with a semi-Markovian state and Markov chain observation models to represent trap effects.
In this article, we consider hidden Markovian models that combine several nonhomogeneous Markovian states with one semi-Markovian state. We show how the length of migratory birds stop overs can be easily written using this extended framework. In the second application, we incorporate accurately known duration of incubation to study the breeding dynamics of a colonial waterbird. The approach proposed is unusual in CR models. Finally, possible generalizations and applications of this work are discussed.

Hidden Hybrid Markov/Semi-Markov Models
Assume that we have T capture occasions and N individuals captured at least once. Let the encounter history for

Model definition
Let S t be a hybrid Markov/semi-Markov model with finite state space 1,…,J; see Kulkarni (1995) for a general description of Markov and semi-Markov models. This J-state hybrid Markov/semi-Markov model is defined by the following parameters: • initial probabilities p j = P (S 1 = j) with ∑ j p j = 1, • transition probabilities semi-Markovian state j: for each k 6 ¼ j; / jk ¼ P ðS tþ1 ¼ kjS tþ1 6 ¼ j; S t ¼ jÞ with ∑ k6 ¼j / jk = 1 and / jj = 0, -Markovian state j: An explicit sojourn time distribution is attached to each semi-Markovian state where M j denotes the upper bound of the time spent in state j. Hence, we assume that sojourn time distributions are concentrated on finite sets of time points. The output process {O t }, here the encounter history h i for individual i, is related to the hybrid Markov/semi-Markov chain {S t } by the observation (or emission) distributions ) ′ denote the column vector of observation probabilities of history h and p denote the vector of initial probabilities. Let U denote the matrix of transition probabilities. In the following, we assume that some observation distributions are a function of the index parameter t introducing some nonhomogeneity into the model.

Estimation problem with one semi-Markovian state
By convention, the semi-Markovian state number is SM, the Markovian transient states that precede this semi-Markovian state are numbered from 1 to SMÀ1 and the successive Markovian states that follow the semi-Markovian state are numbered from SM+1 to J. We consider a model that is conditional on the first occasion. The probability of an encounter history h = (o 1 ,…,o K ) is where diag (x) is a diagonal matrix with ith diagonal element x i and 1 i:j is a J column vector of one's between i and j and 0 elsewhere. The first line of equation (1) corresponds to the probability of no occurrence of the semi-Markovian state. Rows 2 to 4 correspond to the sum of probabilities of occurrence of the semi-Markovian state between the time of arrival ar and the time of departure dep. The likelihood is conditional on the set of individuals seen at least once in the experiment. With the usual assumption of independence among individuals, the likelihood L is the product of the N individual contributions where ø is a history with no capture. In Choquet et al.
(2013b), developed specifically for stopover duration, individuals enter the semi-Markovian state before the first observation. Equation (1) extends this work by allowing individuals to enter or not enter the semi-Markovian state at any occasion. We can also easily consider conditioning on the first capture occasion. To do this, we need to start equation (1) from the first event and define the likelihood

Applications
We assume that the nonobservable Markov/semi-Markov model comprises several transient states and a final absorbing state, with all transient states being nonhomogeneous Markovian except one that is semi-Markovian. For the stopover duration application, the nonobservable Markov/semi-Markov model comprises three states: an initial transient Markovian state followed by a transient semi-Markovian state, and a final absorbing state. A state is said to be transient if the probability of definitively leaving this state without the possibility of returning to it is >0. The states are thus ordered, and each state can be visited exactly once ("left-right" model); see Fig. 1.
For the breeding attendance application, the nonobservable Markov/semi-Markov model comprises six states: two Markovian states forming an initial transient class followed by a transient semi-Markovian state, then two successive transient Markovian states, and a final absorbing state. The initial transient class and the 4 following states are ordered; the states are thus partially ordered (no ordering between states 1 and 2 in the initial transient class), and only the first two states can be visited more than once; see Fig. 2.
The specific assumptions of the two applications are as follows: • We assume that individuals enter in the experiment after the first sampling occasion. We assume that the last sampling occasion individual state is not the semi-Markovian state, that is, sampling covers all the time when individuals are present in the semi-Markovian state of interest (i.e., no right censoring of the time spent in this state).
• We assume that individuals are in a superpopulation and enter the experiment where captures are possible. Likelihood (2) is conditional on an animal being seen somewhere in the experiment, as in Link and Barker (2005); Pradel (1996).
(1) We also need to know for how long an individual has been present. One convenient method for representing stopover duration is to introduce semi-Markovian states where the time spent in these states is explicitly modeled by appropriate sojourn time distributions.
To do this, we consider a three-state model made up of an initial "not yet arrived" transient state which is the state before arrival at the site, an intermediate "arrived" transient state which is the state of arrival at the site, and a final "departed" absorbing state which is the state of leaving the site. The intermediate "arrived" state is semi-Markovian, while the other two "not yet arrived" and "departed" states are Markovian. We assume that "not yet arrived" state is the only possible initial state, that is, p 1 =1. For this initial transient Markovian state, we assume that only the self-transition and the transition to the semi-Markovian "arrived" state are possible with a transition distribution that depends on the time index t. For this nonhomogeneous Markovian state, full timedependent models are used for the transition distribution. Hence, the initial probability vector and the transition probability matrix are given by with / ðtÞ 12 the probability of an individual not present at the site actually arriving between occasions t and t+1. We assume that / ðTÀ1Þ 12 = 1, that is, all individuals are present at the site before occasion T. A parametric or nonparametric sojourn time distribution is attached to the semi-Markovian "arrived" state. In Choquet et al. (2013b), several nongeometric distributions were used and compared.
The set of observations is O t = {0 for "not captured", 1 for "captured"}. Concerning the observation models, we assume that only the "not captured" output can be observed in the initial "not yet arrived" state. The corresponding observation distribution is thus degenerate. For the intermediate "arrived" state, the observation distribution depends on the index parameter t. The observation probability matrix is thus given by with b ðtÞ j ðyÞ ¼ ðB t ð Þ Þ j;yþ1 the probability of an individual being in state j at time t to be observed as y. This model, represented in Fig. 1, can be easily generalized to deal with different sources of heterogeneity. For example, the output process can be adapted to the case where trap effects are present (Choquet et al. 2013b). Neglecting trap effects, sometimes due to the food used to catch the birds and stress caused by handling, can lead to marked bias in estimating stopover duration.

Incorporating incubation period duration
Since 1974, greater flamingos have bred on an artificial island located in the commercial saltpans of Salinde-Giraud (Camargue, southern France). Flamingos generally start incubating in April. At age of 10 days, the chicks gather in small cr eches on the breeding island before moving into the water to join a larger cr eche which may contain up to several thousand individuals (Johnson and Cezilly 2007). Since 1977, 12% of the chicks fledged in the Camargue have been ringed each year with PVC rings engraved with alphanumeric codes that can be read from a distance of up to 300 m. We used the subsequent resightings of these ringed birds as breeding adults at the colony in 1991. In that year, an average of 10 hours per day was spent in the towel observing the colony from April 2 (first egg observed) to July 16, and 2 hours per day were spent at the cr eche, from July 16 to September 8 (end of observation period). Ringed birds were resighted by means of a telescope located in a hide 70 m from the colony. Only flamingos observed at least once as breeders were considered. We divided the breeding season into 18 intervals, each of 10 days. The set of observations is O t = {0 for "not captured", 1 for "seen incubating an egg", 2 for "feeding a chick on the breeding island", 3 for "feeding a chick in the cr eche"}.
We considered a six-state model made up of a transient class with two states "elsewhere/preincubating", and "incubation with failure" , three successive transient states "incubation with success", "feeding on the island" and "feeding in the cr eche", and a final absorbing state "departed after incubation with success". The intermediate transient state "incubation with success" is semi-Markovian, while the other states are Markovian; see Fig. 2 for a representation of this model. We used the fact that the length of the incubation period is fixed and is known to last 29 days. The sojourn time distribution of incubation becomes d SM (3) = 1 with SM= "incubation with success". An individual that has not succeeded in the incubation may move back to the initial state, "elsewhere/preincubating".
We assume here that "elsewhere/preincubating" is the only possible initial state, that is, p 1 = 1. Hence, the initial probability vector and the transition probability matrix are given by: The observation probability matrix is: Because the reproduction pattern changes with time, we need specific transition probabilities for the beginning, the middle, and the end of the season. Furthermore, Schmaltz et al. (2011) demonstrated that rainfall has an effect on reproduction at intervals 7 and 8 (i.e., between occasions 7 and 9). Thus, for states 1 to 3, we consider successive time periods corresponding to groupings of intervals 1 to 3, 4 to 6, 7, 8, 9 to 10, 11 to 15, and 16 to 18 for which transitions are set constant. Furthermore, we assume that / ðt¼15:18Þ 13 ¼ 0 and / ðt¼15:18Þ 23 ¼ 0. For states 4 to 6, we consider three successive time periods, 4 to 6, 7 to 8, and 9 to 18.
We consider also specific capture probabilities for successive time periods, occasions 4 to 7, 8, 9 and 10 to 19 for which capture probabilities are constant. Furthermore, as there is no capture for the first period, b t¼1;...;3 i ¼ 0. Models were built using the E-SURGE program (Choquet et al. 2009), and we computed maximum likelihood estimates (MLE) for each model using a quasi-Newton algorithm (Dennis and Schnabel 1983) with multiple starting points to avoid spurious local minima. MLE and confidence intervals are given in Tables 1 and 2. We did not attempt to find the most parsimonious model here using Akaike information criteria (Akaike 1987) as this was not the goal of this study. Gu edon (1999) proposed a validation methodology relying on the fit of different types of characteristic distributions computed from model parameters to their empirical equivalents extracted from data. In particular, we used the fit of output distributions conditional on the observed sequence at each capture occasion; see Gu edon (2005) for the recursive algorithms for computing state and output distributions conditional on the observed sequence for each successive time t in the case of hidden hybrid Markov/semi-Markov models. We compared the proportion of estimated and observed individuals for each possible output. Figure 3 shows that the model successfully reproduces the general pattern of the observations. We show in Fig. 4 the state distributions at each capture occasion in order to illustrate the dynamics of incubation. The general pattern is shown in Fig. 4 corresponds to expectations regarding the successive steps of greater flamingo breeding dynamics. Birds incubated (successfully or not) from the start of the breeding period to occasion 14 at which point the island was deserted. At this occasion, no birds started to incubate as they would not have had enough time to complete the full breeding cycle before food and climatic conditions deteriorate in the fall. Hatching was first observed at occasion 6 (Fig. 3). From then on, the proportion of birds feeding a chick on the island then in the cr eche increased gradually until the end of the observation period (Fig. 4). Regarding incubation, occasions 7 and 8 appeared to be pivotal given that the proportion of incubating birds that managed to hatch their egg (incubation with success) increased slowly until occasion 7 then peaked at occasions 8 and 9. In contrast after occasions 7 and 8, the proportion of unsuccessful incubating birds (incubation with failure) dropped almost to zero. Heavy rains at occasions 7 and 8 [corresponding to occasions 4 and 5 in Schmaltz et al. (2011)] were shown to cause substantial nest desertion (flooding nests and eggs) and were likely the cause of the failure of early incubating birds. The peak seen in the proportion of successful incubating birds suggests that a second wave of incubating birds settled after the rainy period, with most of these birds being successful. As birds with failed incubation may make a second attempt, the second wave may have been partially made up of renesters. We also computed E (N 21 |o 1 o 2 …o T ) the mean number of times that state 2: "incubation with failure" was followed by state 1: "elsewhere/preincubating". The estimated value for E(N 21 | o 1 o 2 …o T ) was 1.1 [0;2.2] showing that an individual left the incubation area once on average. This result strengthens the proposed model which does not differentiate between individuals arriving for the first time or not in the incubation area.

Extension to Models Incorporating Several Semi-Markovian States
The models in section 2 can be extended by incorporating several semi-Markovian states. In the following, we give several potential applications with more than one semi-Markovian state. If more than one semi-Markov state is incorporated, the approach used for building the Outputs 0: Not captured; 1: Seen incubating an egg; 2: Feeding a chick on the breeding island; 3: Feeding a chick in the crèche.
Non-observable states 1: Elsewhere/pre-incubating; 2: Incubation with failure; 3: Incubation with success; 4: Feeding on the breeding island; 5: Feeding in the crèche; 6: Departed after incubation with success. likelihood 1 becomes intractable. Gu edon (2005) developed a specific expectation-maximization (EM) algorithm for models with semi-Markovian states. In the following, we give several potential applications with more than one transient semi-Markovian state. We distinguish two cases: • semi-Markovian state visited at most once (this semi-Markovian state constitutes a transient class); • semi-Markovian state visited more than once (this semi-Markovian state belongs to a transient class made up of more than one state).
Potential applications with semi-Markovian states visited at most once Lorentzen et al. (2012) estimated hatching success from photo series data in a cliff-nesting seabird. The model comprised two states: "egg" and "chick" and made the assumption that the incubation process started the first time an egg was seen. This is obviously prone to errors because an egg viewed for the first time may have been missed before. In this particular case, we could use one Markovian state "no egg" and two semi-Markovian states: "egg" and "chick". Any information on the time spent in each "egg" and/or "chick" state could also be incorporated.
As heterogeneity is very common in ecological data, it might be relevant to consider two or more classes, each with a specific duration. In this case, the generalization of the likelihood (1) is quite straightforward.
Potential applications with semi-Markovian state visited more than once A hidden hybrid Markov/semi-Markov model has been used by Gu edon (2005) to analyze branching and flowering patterns in plants. In this case, the states corresponded to the successive homogeneous branching and axillary flowering zones along a shoot, and the sojourn time distributions to the lengths of each successive zone.
In plant dynamics, the probability of finding species i in a given quadrat at time t may depend not only on time tÀ1 but also on times before (Van Hulst 1979). We can wonder whether knowledge about the past is required and how much time the species i will remain in that quadrat. This question may be better considered using semi-Markovian states. Biologic cycles are made up of different stages which may have different durations. This is the case for ovarian cycle in cattle . We should also mention here that measurements were used in  as observations given that progesterone concentrations were the measurements from which the state of an individual was deduced. This may also be the case in eco-epidemiology (Choquet et al. 2013b) where the estimated duration of an illness may be time since being ill dependent. As pointed out by Zipkin et al. (2010), care should be taken regarding the time scale over which the disease occurs. As an important remark, we note that we might condition on the first capture for eco-epidemiology CR data. As an important remark, we note that we might condition on the first capture for eco-epidemiology CR data.

Conclusion
In this article, we have developed a new framework for capture-recapture data that can be used to explicitly model the time spent in one state. This formulation is well suited to dealing with imperfect and variable detection. These new models provide a generalization of hidden Markov models which are increasingly used to analyze capture-recapture data. We believe that this framework has many potential applications (see the previous section).
However, work remains to be performed on this kind of model in order to render it more accessible to ecologists. Concerning model estimation, Gu edon (2003Gu edon ( , 2005 proposed an efficient algorithm for hidden hybrid Markov/semi-Markov models potentially incorporating several successive semi-Markovian states. Chaubert-Pereira et al. (2010) extended these algorithms to models incorporating fixed and random effects in the observation process. O'Connell and Hojsgaard (2011) developed an R-package for hidden semi-Markov models. And we implemented the likelihood (1) in E-SURGE (Choquet et al. 2009). However, in a broader framework, we need not only to  adapt these algorithms to the specificity of CR data, but also design a software package.