Burnham model in RMARK with dot notation in capture history

[luke.eberhart]

Hi team,

I'm studying sex-specific survival in snowy plover chicks and I'm having a few problems in RMARK while using Burnham's Joint Live and Dead Encounters model with encounter histories that include "dot notation".

To introduce you to my study system, plover chicks are brooded for 27 days since hatch, after which they fledge. Thus, my dataset has 27 encounter occasions that represent each day of a chick's life. Each chick has one of three fates: 1) dead (i.e., known death before fledging), 2) fledged (i.e., made it past 27 days), or 3) unknown. Because of the "joint" nature of these fates in my data, I figured that Burnham's Joint Live and Dead Encounters model would be most appropriate since it would utilize as much fate information as possible to produce the most precise parameter estimates. As such, I coded my data according to the standard LD format. With 27 encounter occasions, there are thus a total of 54 characters in the LD capture history string of each individual.

To give you a taste of my data, here are the first 20 rows (i.e., individuals) of the LD dataframe in R (Note: the first column is the LD capture history string for each individual, and the second column is the sex of each individual):

Code: Select all

> plover[c(1:20),]
                                                    ch Sex
10..001010000010..000010101000100000101000..1000001010   M
10..001010000010..000010101000100000101000..1000001010   M
10..001010000010..000010101000100000101000..1000001010   F
..10101000..10..10101010101010101010101000101000000010   F
..10101000..10..10101010101010101010101000101000000010   F
1010..10..10001010101010101000100000001000000010000010   M
1010..10..10101010100000101000101000100000101010001010   F
1010..01..00000000000000000000000000000000000000000000   M
0010..10..10001010101010101000100000001000000010000010   M
100010101000..............................00..00..0001   M
001010101000..10..10100010100010........01....00000000   F
001010101000..10..10100010100001..............00000000   F
001001000000..00..00000000000000..............00000000   M
100010100100..0000000000..00..0000..0000....0000000000   M
100000100000..0001000000..00..0000..0000....0000000000   F
100010100100..0000000000..00..0000..0000....0000000000   M
1000101010010000000000000000000000000000..000000..00..   M
1000101010010000000000000000000000000000..000000..00..   M
101000100010010000000000..0000000000000000000000000000   F
101000101010100010101000..1010000010000000100010001010   F

At first glance, you may be wondering "what the heck are those dots in the capture history?" The simple answer is that a ".." codes for a day (in LD format) in which there was no survey conducted, and thus it was not possible to resight (or sight) a chick. I coded the capture history this way to acknowledge the fact that days in which a survey was not conducted vs. days in which a survey was done but the chick was not detected are not the same. If these two cases were coded the same in the capture history (i.e., with a "0"), then the resight history would be negatively biased and the precision of the survival estimate would be reduced. This is simply my intuition, but I welcome your thoughts and expertise about this issue.

After reading the data into R, I then followed the standard RMARK pipeline by first processing the data via specifying the "Burnham" model and defining the "Sex" variable as a group:

Code: Select all

plover.proc <- process.data(plover, model = "Burnham", groups = c("Sex"))

Then I created the design matrix:

Code: Select all

plover.ddl <- make.design.data(plover.proc)

Then I tested a simple model (S(~Sex)p(~1)r(~1)F(~1)) to see if the data were being properly read:

Code: Select all

S.sex <- list(formula=~Sex) # S(sex)
p.dot <- list(formula=~1) # p(.)
F.dot <- list(formula=~1) # F(.)
r.dot <- list(formula=~1) # r(.)
Simple.model <- mark(data = plover.proc,
                  ddl = plover.ddl,
                  model.parameters = list(S = S.sex, p = p.dot, F = F.dot, r = r.dot),
                  invisible = FALSE, model = "Burnham")

However, when I ran this simple model, I got the following error message:

Code: Select all

STOP ERROR
MARK did not run properly.  If error message was not shown, re-run MARK with invisible=FALSE

Model output is not available

Error in wtable[i, cols] = real[pim[i, cols]] :
  number of items to replace is not a multiple of replacement length

Does anyone have any ideas what this error message could be referring to or how to solve this problem? My gut feeling is that RMARK doesn't like that I've used the dot notation in my capture history. However, I tried these same data (but in LL data format) in a classic CJS model in RMARK using very similar code and it worked just fine. Could it be that the Burnham model doesn't recognise dots in the capture history?

Your help is greatly appreciated!

Cheers,
Luke

[jlaake]

Look at the most recent .out file in the directory where you ran R with RMark. There should be an error message in the out file.

[luke.eberhart]

Hi Jeff,

Many thanks for your quick reply. You were totally right, after checking the .OUT file, the error message says:

Code: Select all

The character '.' was found in the encounter history.

 ERROR -- Encounter history must consist of only these characters: '120'

Any thoughts on how to proceed? Am I correct to think that MARK is not recognizing the "." in a Burnham model (i.e., as opposed to RMARK not recognizing the ".")?

I appreciate any further suggestions!

Cheers,
Luke

[jlaake]

Yes. Any message in the .out file is from MARK. I must be letting it slip through RMark. I'm surprised that it says 120 are allowed characters. 2 is for resighting which is only in the Barker model.

Your option is to replace the .'s with 0's and fix the parameters - p to 0 and r to 0. If it is a gap in all of the capture histories you can skip that occasion and set time.intervals to >1 for that span between occasions but it looks like they vary across individuals. Can you group the sets of chicks such that a set of capture histories has the same missing values? If so then you can fix parameters by group. Only way I can see to do it. In my R package marked I allow time.intervals and fixed parameters to vary by individual animal but I don't have the Burnham model programmed into it. This comes up a lot and Gary may want to consider allowing time.intervals to vary across individuals or groups at least. Or allow . notation for this model.

See http://www.phidot.org/forum/viewtopic.php?f=21&t=3035 starting at the 4th posting.

regards --jeff

[luke.eberhart]

Hi Jeff,

This is helpful, thanks. Regarding your comment:

Your option is to replace the .'s with 0's and fix the parameters - p to 0 and r to 0.

Yes, this does appear to be the only option. However (as you noticed) the location of missing occasions in the capture histories vary across individuals. This is because the structure of each individual's capture history is determined by three parameters I hadn't mentioned until now: 1) the chick's hatch date (i.e., the first capture occasion is the first day of a chick's life), 2) the site where the chick is being brooded (i.e., not every site was surveyed daily and thus the missing occasions), and 3) the year (i.e., the survey schedule varies by year). To illustrate the importance of these three parameters, check out the first 6 individuals of the dataframe, this time with three extra columns (the hatch date (Julian), the brooding site, and year):

Code: Select all

                                                    ch  Sex  Hatch_date Site  Year
10..001010000010..000010101000100000101000..1000001010   M       12       C   2009
10..001010000010..000010101000100000101000..1000001010   M       12       C   2009
10..001010000010..000010101000100000101000..1000001010   F       12       C   2009
..10101000..10..10101010101010101010101000101000000010   F       13       A   2009
..10101000..10..10101010101010101010101000101000000010   F       13       A   2009
1010..10..10001010101010101000100000001000000010000010   M       14       A   2009

You can clearly see that the first three individuals all hatched on the same date and year (12, 2009) and at the same site (C), thus their missing occasions all line up nicely (i.e., ".." appear on the 2nd, 9th, and 22nd day of each of these three chick's lives). Individuals 4 and 5 follow this same line of logic, as they both hatched on the same date and year (13, 2009) at the same site (A). Individual 6 is interesting because it hatched on site A one day later than individuals 4 and 5, thus the locations of the ".." in it's capture history are shifted one day to the left. Does that make sense?

The posting you referred me to demonstrated how to "pad" a capture history with "0"s and then fix parameters by group to get around similar issues that I'm facing. I think this method would work if the structure of my data didn't have so many hierarchical groupings (i.e., hatch date, then site, then year). What I mean to say is that based on these three group levels, there are a total of 124 "survey groups" that I would need to fix the p and r parameters to 0:

Code: Select all

plover$Survey_Group <- as.factor(paste(plover$Year, plover$Hatch_date, plover$Site, sep = "_"))
length(levels(plover$Survey_Group))
[1] 124


Maybe fixing the the parameters for 124 groups is not as challenging as I think it is, but I anticipate that it would be a big hassle to wrangle the data occasion-by-occasion and group-by-group.

So, given this extra information about the structure of my data, do you still think there is a way to work around the problem that "."s are not recognized in the Burnham model? Any more suggestions would be appreciated. Or perhaps you could gently "nudge" Gary to write the FORTRAN code that would enable MARK to recognize "." in the Burnham model.

Cheers,
Luke

[jlaake]

I think you can trim the number of groups to year and site. Leading 0's to the first 1 are ignored because the Burnham model is like CJS being conditioned on the first release (first 1). Thus you can structure the capture histories starting with the first day a chick was born that year. Chicks not born until later are given 00 for each occasion until they are hatched. That should allow you to have year by site groups.

Now because you aren't searching the site each day, you have to be careful with how you model the parameters. Clearly they can't be estimated for each occasion. Otherwise this should work but I haven't thought it through completely so you may run into issues. Have you considered using the nest survival model?

--jeff

[luke.eberhart]

Hi Jeff,

In regards to:

...you can structure the capture histories starting with the first day a chick was born that year. Chicks not born until later are given 00 for each occasion until they are hatched. That should allow you to have year by site groups.

This could be a good solution. However, to code the capture histories correctly this way, wouldn't I need to introduce a new variable describing the age at which the chick was first marked? For example, in many cases the first "1" of a capture history wasn't on the day of hatch, but was when the chick was two or three days old. One of the main parameters I'm interested in is how daily survival varies across age (i.e., over the 27 days from hatch to fledge) and so it would be important to code this correctly.

Okay, so would this be the correct workflow for my next stage of data-wrangling:

1) make a new variable describing the age (in days) that the chick was first encountered (i.e., as mentioned in the previous paragraph)

2) group individuals according to site and year

3) insert the correct numbers of "0"s at the beginning and end of each capture history

4) use your code found here to fix (by group!) the p and r to 0 for days that had no surveys (i.e., those days with the ".." in my original dataframe.)

5) simply replace all "." characters in all capture history strings with "0" (Am I correct to assume that the missing data info is now redundant since it is embedded in the time-and-group-varying fixed parameter estimates defined in step 4?)

6) try rerunning the Burnham model, this time using "age" as the covariate as opposed to "time" (which is what I was using in the first attempt since all capture histories were 27 days long ranging from hatch to fledge).

Have you considered using the nest survival model?

Yes, I considered using a classic nest survival analysis (which has been done before for analysis of chick survival (e.g., Colwell et al. 2007 Auk)), however by using a nest survival approach, one must be super confident with resighting all living chicks on each survey (the main difference between a nest and a chick being that former is a stationary object and the later is a critter on the move). In my study system, I'm not confident with this assumption since plover chicks are mighty cryptic, hence why I choose to use a Burnham model. I'm not set in stone on using a Burnham model however, and so I welcome any further suggestions or alternative models that may fit my research aims (being: to assess sex- and age-specific chick survival using mark-recapture methods)

Thanks again Jeff, I appreciate your time.

Cheers,
Luke

[jlaake]

Understood with regard to the nest survival. I wasn't certain of the situation with your chicks.

I think you are going to what to use a time-varying individual covariate for age rather than using initial.age as another group variable although that is possible as well. Here is an example

10100001 0 1 2 3
00100000 2 3 4 5

The first chick was born on the first day and then aged one day for each occasion. The second chick was 3 days old the first day it was seen so it was 2 on the first occasion. Those values prior to the first one aren't important however because they aren't used in the likelihood. If your begin.time is 1 then your age variable should be named age1,age2,age3,age4,... for S. If age is also going to affect p then you need to have the last value (age4 in this example) because age1,age2,age3 will be used for S and age2,age3,age4 will be used for p because S is an interval parameter and p is an occasion parameter. If you have a different begin.time value then the variable names need to be changed to match the times in the design data.

Where you used to have "." you'll have "0" and you'll set p=0 and r=0. What I'm not certain about is recovery. A chick that dies between occasions when you don't survey could be recovered the next time you survey. Not sure what this will do.

I'm not sure I follow your #3. Are the 0s at the end of the capture history because there are different numbers of survey days across years? The 0s at the beginning should just be replacing where you previously had "." but maybe I'm not understanding something. I didn't understand your comment in 5 about missing data info. You'll always know the chicks age even if it was not seen or the site was not surveyed.

regards --jeff