time-varying individual covariates with Burnham model

posts related to the RMark library, which may not be of general interest to users of 'classic' MARK

time-varying individual covariates with Burnham model

Postby jbauder » Sat Oct 28, 2017 12:53 pm

Hello,

I have a mix of PIT-tagged and telemetered individuals with 22 capture occasions. While all individuals were PIT-tagged for the duration of the study, some were also telemetered at varying times for varying durations. I am using Burnham’s joint recapture/recovery model and would like to fix p=1 for the telemetered individuals during the intervals in which they were telemetered. It seems like a time-varying individual covariate (say “cov”) would be the way to accomplish this, where “cov” was 1 if the individual had a functioning transmitter and 0 otherwise. For starters, I tried fitting a model where p was modeled as a function of “Telem” [i.e., S(.)p(cov)r(.)F(.)].

However, when I fit this model I had 21 (k-1) real parameter estimates for p. This wasn’t what I expected so I wanted to see if my results are indeed legitimate or I have set my data up incorrectly.

My data consists of a 23-column data frame, the first column containing the capture histories (44 values in LD format for 22 occasions) and 22 columns denoting the “cov” status (0/1) for each individual on that occasion. The 22 “cov” columns are labeled “cov1”, “cov2”, …, “cov22”.


Code: Select all
> str(CH)
'data.frame':   99 obs. of  23 variables:
 $ ch   : chr  "10000000000000000000000000000000000000000000" "10000000000000000000000000000000000000000000" "11000000000000000000000000000000000000000000" "10000000000000000000000000000000000000000000" ...
 $ cov1 : num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov2 : num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov3 : num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov4 : num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov5 : num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov6 : num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov7 : num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov8 : num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov9 : num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov10: num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov11: num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov12: num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov13: num  0 0 0 0 1 0 0 0 0 0 ...
 $ cov14: num  0 0 0 0 1 0 0 0 0 0 ...
 $ cov15: num  0 0 0 0 1 0 0 0 0 0 ...
 $ cov16: num  0 0 0 0 1 0 0 0 0 0 ...
 $ cov17: num  0 0 0 0 1 0 0 0 0 0 ...
 $ cov18: num  0 0 0 0 1 0 0 0 0 0 ...
 $ cov19: num  0 0 0 0 0 1 0 0 0 0 ...
 $ cov20: num  0 0 0 0 0 1 0 0 0 0 ...
 $ cov21: num  0 0 0 0 0 0 0 0 0 0 ...
 $ cov22: num  0 0 0 0 0 0 0 0 0 0 ...

t1_process <- process.data(data = CH,model="Burnham")
t1.ddl <- make.design.data(t1_process)
> head(t1_process$data)
                                            ch cov1 cov2 cov3 cov4 cov5 cov6 cov7 cov8 cov9 cov10 cov11 cov12
1 10000000000000000000000000000000000000000000    0    0    0    0    0    0    0    0    0     0     0     0
2 10000000000000000000000000000000000000000000    0    0    0    0    0    0    0    0    0     0     0     0
3 11000000000000000000000000000000000000000000    0    0    0    0    0    0    0    0    0     0     0     0
4 10000000000000000000000000000000000000000000    0    0    0    0    0    0    0    0    0     0     0     0
5 10000010101010100000100010101010101000000000    0    0    0    0    0    0    0    0    0     0     0     0
6 10101010101010000000000000000000000010100100    0    0    0    0    0    0    0    0    0     0     0     0
  cov13 cov14 cov15 cov16 cov17 cov18 cov19 cov20 cov21 cov22 freq
1     0     0     0     0     0     0     0     0     0     0    1
2     0     0     0     0     0     0     0     0     0     0    1
3     0     0     0     0     0     0     0     0     0     0    1
4     0     0     0     0     0     0     0     0     0     0    1
5     1     1     1     1     1     1     0     0     0     0    1
6     0     0     0     0     0     0     1     1     0     0    1

t1 <- mark(data = t1_process,
           ddl = t1.ddl,
           model.parameters = list(S = list(formula=~1), p = list(formula=~cov),
                                   F = list(formula=~1,fixed=1), r = list(formula=~1)),
           model = "Burnham")
> t1$results$beta
               estimate           se           lcl          ucl
S:(Intercept)  1.830187    0.1464575     1.5431304    2.1172437
p:(Intercept) -0.438305    0.1407934    -0.7142601   -0.1623500
p:cov         20.501717 2454.5176000 -4790.3529000 4831.3564000
r:(Intercept) -1.023930    0.2801911    -1.5731046   -0.4747556
> t1$results$real
                 estimate         se           lcl       ucl fixed                 note
S g1 c1 a0 t1   0.8617840  0.0174449  8.239193e-01 0.8925679                           
p g1 c1 a1 t2   0.4424519  6.1162633  6.236393e-22 1.0000000                           
p g1 c1 a2 t3   0.5456126 18.4401140  5.839770e-64 1.0000000                           
p g1 c1 a3 t4   0.3921449  0.0335605  3.286582e-01 0.4595014                           
p g1 c1 a4 t5   0.4424519  6.1162633  6.236393e-22 1.0000000                           
p g1 c1 a5 t6   0.3921449  0.0335605  3.286582e-01 0.4595014                           
p g1 c1 a6 t7   0.4424519  6.1162633  6.236393e-22 1.0000000                           
p g1 c1 a7 t8   0.4424519  6.1162633  6.236393e-22 1.0000000                           
p g1 c1 a8 t9   0.4424519  6.1162633  6.236393e-22 1.0000000                           
p g1 c1 a9 t10  0.4424519  6.1162633  6.236393e-22 1.0000000                           
p g1 c1 a10 t11 0.4424519  6.1162633  6.236393e-22 1.0000000                           
p g1 c1 a11 t12 0.4424519  6.1162633  6.236393e-22 1.0000000                           
p g1 c1 a12 t13 0.6908772 31.7697870 5.288171e-127 1.0000000                           
p g1 c1 a13 t14 0.6908772 31.7697870 5.288171e-127 1.0000000                           
p g1 c1 a14 t15 0.7332767 33.9436050 5.116009e-148 1.0000000                           
p g1 c1 a15 t16 0.6908772 31.7697870 5.288171e-127 1.0000000                           
p g1 c1 a16 t17 0.7717817 34.9354160 4.949429e-169 1.0000000                           
p g1 c1 a17 t18 0.8061972 34.8638080 4.788258e-190 1.0000000                           
p g1 c1 a18 t19 0.9561522 17.6707520 2.150015e-303 1.0000000                           
p g1 c1 a19 t20 0.9561522 17.6707520 2.150015e-303 1.0000000                           
p g1 c1 a20 t21 0.9561522 17.6707520 2.150015e-303 1.0000000                           
p g1 c1 a21 t22 0.9640592 15.4630380 2.644714e-303 1.0000000                           
r g1 c1 a0 t1   0.2642626  0.0544769  1.717743e-01 0.3834913                           
F g1 c1 a0 t1   1.0000000  0.0000000  1.000000e+00 1.0000000 Fixed

> t1$design.matrix
                S:(Intercept) p:(Intercept) p:cov   r:(Intercept)
S g1 c1 a0 t1   "1"           "0"           "0"     "0"         
p g1 c1 a1 t2   "0"           "1"           "cov2"  "0"         
p g1 c1 a2 t3   "0"           "1"           "cov3"  "0"         
p g1 c1 a3 t4   "0"           "1"           "cov4"  "0"         
p g1 c1 a4 t5   "0"           "1"           "cov5"  "0"         
p g1 c1 a5 t6   "0"           "1"           "cov6"  "0"         
p g1 c1 a6 t7   "0"           "1"           "cov7"  "0"         
p g1 c1 a7 t8   "0"           "1"           "cov8"  "0"         
p g1 c1 a8 t9   "0"           "1"           "cov9"  "0"         
p g1 c1 a9 t10  "0"           "1"           "cov10" "0"         
p g1 c1 a10 t11 "0"           "1"           "cov11" "0"         
p g1 c1 a11 t12 "0"           "1"           "cov12" "0"         
p g1 c1 a12 t13 "0"           "1"           "cov13" "0"         
p g1 c1 a13 t14 "0"           "1"           "cov14" "0"         
p g1 c1 a14 t15 "0"           "1"           "cov15" "0"         
p g1 c1 a15 t16 "0"           "1"           "cov16" "0"         
p g1 c1 a16 t17 "0"           "1"           "cov17" "0"         
p g1 c1 a17 t18 "0"           "1"           "cov18" "0"         
p g1 c1 a18 t19 "0"           "1"           "cov19" "0"         
p g1 c1 a19 t20 "0"           "1"           "cov20" "0"         
p g1 c1 a20 t21 "0"           "1"           "cov21" "0"         
p g1 c1 a21 t22 "0"           "1"           "cov22" "0"         
r g1 c1 a0 t1   "0"           "0"           "0"     "1"         
F g1 c1 a0 t1   "0"           "0"           "0"     "0"


It seems (and I would like) to have two real parameter estimates, one when “cov”=0 and the other when “cov”=1 but I haven’t been able to do this, even though I feel like this must be a simple fix. Any guidance would be much appreciated!

Also, if I do get this to work correctly, would the following syntax be used to fix p=1 when “cov=1”? It threw an error message when I tried it with my current set up, again suggesting something is wrong with my data formatting(?).

Code: Select all
> t1 <- mark(data = t1_process,
+            ddl = t1.ddl,
+            model.parameters = list(S = list(formula=~1),
+                                    p = list(formula=~cov,fixed=list(cov=1,value=1)),
+                                    F = list(formula=~1,fixed=1), r = list(formula=~1)),
+            model = "Burnham")
Error in make.mark.model(data.proc, title = title, parameters = model.parameters,  :
 
Unrecognized structure for fixed parameters = 1
Unrecognized structure for fixed parameters = 1
Error in mark(data = t1_process, ddl = t1.ddl, model.parameters = list(S = list(formula = ~1),  :
  Misspecification of model or internal error in code


Thanks so much,
Javan
jbauder
 
Posts: 56
Joined: Wed May 25, 2011 12:01 pm

Re: time-varying individual covariates with Burnham model

Postby jlaake » Sat Oct 28, 2017 6:31 pm

I you look at the betas you'll see that you are getting what you want. You are misinterpretingoing the real parameter estimates that are being reported. When MARK has an individual covariate it reports the real parameter for the mean value of the covariance which in this case is the proportion that are telemetered which varies across time.
jlaake
 
Posts: 1417
Joined: Fri May 12, 2006 12:50 pm
Location: Escondido, CA

Re: time-varying individual covariates with Burnham model

Postby jbauder » Sat Oct 28, 2017 6:47 pm

jlaake wrote:I you look at the betas you'll see that you are getting what you want. You are misinterpretingoing the real parameter estimates that are being reported. When MARK has an individual covariate it reports the real parameter for the mean value of the covariance which in this case is the proportion that are telemetered which varies across time.


Hi Jeff,

Thanks for your quick response (glad to hear I was setting up my data correctly!). I think I understand how these real parameter estimates are obtained but am still a little confused on how they are interpreted. If the p=0.39 for cov6 (occasion 6), is that value still interpreted as the probability of recapture for non-telemetered individuals? Or is it something else? Would the same interpretation hold if I used this covariate for modeling survival?

How would you recommend interpreting my very large SE's/CI? For example, on occasion 1, p=0.44 (which is reasonable) but SE=6.16. Could this be because I have so few telemetered snakes on that occasion?

If I could tag on one more question, how should I fix the recapture probability for telemetered snakes to 1 (since we were always able to locate a telemetered snake on a given occasion)?

Thanks again so much for your reply!

Javan
jbauder
 
Posts: 56
Joined: Wed May 25, 2011 12:01 pm

Re: time-varying individual covariates with Burnham model

Postby jlaake » Mon Oct 30, 2017 2:37 pm

The large std errors are due to the beta for telemetered animals being at a boundary. Don't even try to interpret the real problem estimates. They are nonsense. p=1 for telemetered animals and for non telemetered p is inverse logit of the beta intercept for p. That is interpretable and will have a reasonable std error. You should read the material on individual covariates in MARK book and in RMark workshop notes.
jlaake
 
Posts: 1417
Joined: Fri May 12, 2006 12:50 pm
Location: Escondido, CA

Re: time-varying individual covariates with Burnham model

Postby jbauder » Mon Oct 30, 2017 4:57 pm

I see that now! Sorry for not catching that inverse logit trick earlier. Not sure what I was missing before, but thank you very much for helping to clarify it for me.
jbauder
 
Posts: 56
Joined: Wed May 25, 2011 12:01 pm


Return to RMark

Who is online

Users browsing this forum: No registered users and 13 guests