Occasion covariates -RD

[sergi_20]

Hi,

I am running a robust design analysis in Rmark with 54 capture histories (12 primary occasions and secondary occasions varying from 2 to 7). The samples were collected every three months, except for two occasions, one with a six month interval and the other one with 12 month interval.

I want to test a set of occasion covariates but I am having problems to input the covariates into the design matrix. I am including 12 primary occasions covariates (effort) and I am getting a design matrix with effort for each secondary sampling occasion.

Can anyone help me to see what I am missing?

Code: Select all

#convert inp file
rd.data=convert.inp("C:/R/RD/Qi 2-3 (5th fortnight-all days)_without calves and sightings(identified dolphins).inp",use.comments=T)

#Process data specifying primary and secondary capture occasions
time.intervals=c(0,0,0,0,0,3,0,0,3,0,3,0,0,0,3,0,0,0,0,0,3,0,0,0,0,3,0,0,0,0,3,0,0,0,12,0,0,0,0,0,0,3,0,0,3,0,0,3,0,0,0,0,0)
rd.process=process.data(rd.data,begin.time=1,model="Robust",time.intervals=time.intervals)

#Create the design data
rd.ddl=make.design.data(rd.process)

df=data.frame(time=c(1:12),effort=c(2,4,5,4,5,6,17,28,89,120,91,22))
summary(rd.ddl$p)

rd.ddl$p=merge_design.covariates(rd.ddl$p,df)
summary(rd.ddl$p)

                             
> rd.ddl$p
   time par.index model.index group session Time c effort
1     1         1         188     1       1    0 0      2
2     2         2         189     1       1    1 0      4
3     3         3         190     1       1    2 0      5
4     4         4         191     1       1    3 0      4
5     5         5         192     1       1    4 0      5
6     6         6         193     1       1    5 0      6
7     1         7         194     1       4    0 0      2
8     2         8         195     1       4    1 0      4
9     3         9         196     1       4    2 0      5
10    1        10         197     1       7    0 0      2
11    2        11         198     1       7    1 0      4
12    1        12         199     1      10    0 0      2
13    2        13         200     1      10    1 0      4
14    3        14         201     1      10    2 0      5
15    4        15         202     1      10    3 0      4
16    1        16         203     1      13    0 0      2
17    2        17         204     1      13    1 0      4
18    3        18         205     1      13    2 0      5
19    4        19         206     1      13    3 0      4
20    5        20         207     1      13    4 0      5
21    6        21         208     1      13    5 0      6
22    1        22         209     1      16    0 0      2
23    2        23         210     1      16    1 0      4
24    3        24         211     1      16    2 0      5
25    4        25         212     1      16    3 0      4
26    5        26         213     1      16    4 0      5
27    1        27         214     1      19    0 0      2
28    2        28         215     1      19    1 0      4
29    3        29         216     1      19    2 0      5
30    4        30         217     1      19    3 0      4
31    5        31         218     1      19    4 0      5
32    1        32         219     1      22    0 0      2
33    2        33         220     1      22    1 0      4
34    3        34         221     1      22    2 0      5
35    4        35         222     1      22    3 0      4
36    1        36         223     1      34    0 0      2
37    2        37         224     1      34    1 0      4
38    3        38         225     1      34    2 0      5
39    4        39         226     1      34    3 0      4
40    5        40         227     1      34    4 0      5
41    6        41         228     1      34    5 0      6
42    7        42         229     1      34    6 0     17
43    1        43         230     1      37    0 0      2
44    2        44         231     1      37    1 0      4
45    3        45         232     1      37    2 0      5
46    1        46         233     1      40    0 0      2
47    2        47         234     1      40    1 0      4
48    3        48         235     1      40    2 0      5
49    1        49         236     1      43    0 0      2
50    2        50         237     1      43    1 0      4
51    3        51         238     1      43    2 0      5
52    4        52         239     1      43    3 0      4
53    5        53         240     1      43    4 0      5
54    6        54         241     1      43    5 0      6


[jlaake]

At present the merge_design.covariates doesn't work with RD to assign covariates at the session level. If it did work with session you would have to use the values 1,4,7.... that matched the session time values. It is not difficult to write code to do what you want. I just haven't gotten around to it. All it is doing is using the merge function in R but first it assigns a sequence number to the ddl dataframe for that parameter so it can sort to the original order. Look at the code for merge_design.covariates. If you can't work it out, let me know. Alternatively, you can create a vector with the effort values repeated and cbind it to the ddl dataframe.

--jeff

[sergi_20]

I wrote the following code to add the occasion covariates. I had to do it for p and also for c.

Code: Select all

rd.ddl$p$Hours[rd.ddl$p$session==1]=50.97
rd.ddl$p$Hours[rd.ddl$p$session==4]=27.01
rd.ddl$p$Hours[rd.ddl$p$session==7]=44.00
rd.ddl$p$Hours[rd.ddl$p$session==10]=28.42
rd.ddl$p$Hours[rd.ddl$p$session==13]=35.58
rd.ddl$p$Hours[rd.ddl$p$session==16]=28.15
rd.ddl$p$Hours[rd.ddl$p$session==19]=28.20
rd.ddl$p$Hours[rd.ddl$p$session==22]=36.52
rd.ddl$p$Hours[rd.ddl$p$session==34]=29.82
rd.ddl$p$Hours[rd.ddl$p$session==37]=15.58
rd.ddl$p$Hours[rd.ddl$p$session==40]=28.17
rd.ddl$p$Hours[rd.ddl$p$session==43]=45.68

rd.ddl$c$Hours[rd.ddl$c$session==1]=50.97
rd.ddl$c$Hours[rd.ddl$c$session==4]=27.01
rd.ddl$c$Hours[rd.ddl$c$session==7]=44.00
rd.ddl$c$Hours[rd.ddl$c$session==10]=28.42
rd.ddl$c$Hours[rd.ddl$c$session==13]=35.58
rd.ddl$c$Hours[rd.ddl$c$session==16]=28.15
rd.ddl$c$Hours[rd.ddl$c$session==19]=28.20
rd.ddl$c$Hours[rd.ddl$c$session==22]=36.52
rd.ddl$c$Hours[rd.ddl$c$session==34]=29.82
rd.ddl$c$Hours[rd.ddl$c$session==37]=15.58
rd.ddl$c$Hours[rd.ddl$c$session==40]=28.17
rd.ddl$c$Hours[rd.ddl$c$session==43]=45.68


And I got this:

Code: Select all

$p
   par.index model.index group time session Time c Hours Tourist Boats
1          1         188     1    1       1    0 0 50.97   20059  1224
2          2         189     1    2       1    1 0 50.97   20059  1224
3          3         190     1    3       1    2 0 50.97   20059  1224
4          4         191     1    4       1    3 0 50.97   20059  1224
5          5         192     1    5       1    4 0 50.97   20059  1224
6          6         193     1    6       1    5 0 50.97   20059  1224
7          7         194     1    1       4    0 0 27.01    7629   461
8          8         195     1    2       4    1 0 27.01    7629   461
9          9         196     1    3       4    2 0 27.01    7629   461
10        10         197     1    1       7    0 0 44.00   20332   829
11        11         198     1    2       7    1 0 44.00   20332   829
12        12         199     1    1      10    0 0 28.42   18227   822
13        13         200     1    2      10    1 0 28.42   18227   822
14        14         201     1    3      10    2 0 28.42   18227   822
15        15         202     1    4      10    3 0 28.42   18227   822
16        16         203     1    1      13    0 0 35.58   19721   759
17        17         204     1    2      13    1 0 35.58   19721   759
18        18         205     1    3      13    2 0 35.58   19721   759
19        19         206     1    4      13    3 0 35.58   19721   759
20        20         207     1    5      13    4 0 35.58   19721   759
21        21         208     1    6      13    5 0 35.58   19721   759
22        22         209     1    1      16    0 0 28.15    6873   288
23        23         210     1    2      16    1 0 28.15    6873   288
24        24         211     1    3      16    2 0 28.15    6873   288
25        25         212     1    4      16    3 0 28.15    6873   288
26        26         213     1    5      16    4 0 28.15    6873   288
27        27         214     1    1      19    0 0 28.20   15765   614
28        28         215     1    2      19    1 0 28.20   15765   614
29        29         216     1    3      19    2 0 28.20   15765   614
30        30         217     1    4      19    3 0 28.20   15765   614
31        31         218     1    5      19    4 0 28.20   15765   614
32        32         219     1    1      22    0 0 36.52    8973   370
33        33         220     1    2      22    1 0 36.52    8973   370
34        34         221     1    3      22    2 0 36.52    8973   370
35        35         222     1    4      22    3 0 36.52    8973   370
36        36         223     1    1      34    0 0 29.82   12673   547
37        37         224     1    2      34    1 0 29.82   12673   547
38        38         225     1    3      34    2 0 29.82   12673   547
39        39         226     1    4      34    3 0 29.82   12673   547
40        40         227     1    5      34    4 0 29.82   12673   547
41        41         228     1    6      34    5 0 29.82   12673   547
42        42         229     1    7      34    6 0 29.82   12673   547
43        43         230     1    1      37    0 0 15.58   10875   532
44        44         231     1    2      37    1 0 15.58   10875   532
45        45         232     1    3      37    2 0 15.58   10875   532
46        46         233     1    1      40    0 0 28.17    3694   262
47        47         234     1    2      40    1 0 28.17    3694   262
48        48         235     1    3      40    2 0 28.17    3694   262
49        49         236     1    1      43    0 0 45.68   12143   543
50        50         237     1    2      43    1 0 45.68   12143   543
51        51         238     1    3      43    2 0 45.68   12143   543
52        52         239     1    4      43    3 0 45.68   12143   543
53        53         240     1    5      43    4 0 45.68   12143   543
54        54         241     1    6      43    5 0 45.68   12143   543

$c
   par.index model.index group time session Time c Hours Tourist Boats
1          1         242     1    2       1    0 1 50.97   20059  1224
2          2         243     1    3       1    1 1 50.97   20059  1224
3          3         244     1    4       1    2 1 50.97   20059  1224
4          4         245     1    5       1    3 1 50.97   20059  1224
5          5         246     1    6       1    4 1 50.97   20059  1224
6          6         247     1    2       4    0 1 27.01    7629   461
7          7         248     1    3       4    1 1 27.01    7629   461
8          8         249     1    2       7    0 1 44.00   20332   829
9          9         250     1    2      10    0 1 28.42   18227   822
10        10         251     1    3      10    1 1 28.42   18227   822
11        11         252     1    4      10    2 1 28.42   18227   822
12        12         253     1    2      13    0 1 35.58   19721   759
13        13         254     1    3      13    1 1 35.58   19721   759
14        14         255     1    4      13    2 1 35.58   19721   759
15        15         256     1    5      13    3 1 35.58   19721   759
16        16         257     1    6      13    4 1 35.58   19721   759
17        17         258     1    2      16    0 1 28.15    6873   288
18        18         259     1    3      16    1 1 28.15    6873   288
19        19         260     1    4      16    2 1 28.15    6873   288
20        20         261     1    5      16    3 1 28.15    6873   288
21        21         262     1    2      19    0 1 28.20   15765   614
22        22         263     1    3      19    1 1 28.20   15765   614
23        23         264     1    4      19    2 1 28.20   15765   614
24        24         265     1    5      19    3 1 28.20   15765   614
25        25         266     1    2      22    0 1 36.52    8973   370
26        26         267     1    3      22    1 1 36.52    8973   370
27        27         268     1    4      22    2 1 36.52    8973   370
28        28         269     1    2      34    0 1 29.82   12673   547
29        29         270     1    3      34    1 1 29.82   12673   547
30        30         271     1    4      34    2 1 29.82   12673   547
31        31         272     1    5      34    3 1 29.82   12673   547
32        32         273     1    6      34    4 1 29.82   12673   547
33        33         274     1    7      34    5 1 29.82   12673   547
34        34         275     1    2      37    0 1 15.58   10875   532
35        35         276     1    3      37    1 1 15.58   10875   532
36        36         277     1    2      40    0 1 28.17    3694   262
37        37         278     1    3      40    1 1 28.17    3694   262
38        38         279     1    2      43    0 1 45.68   12143   543
39        39         280     1    3      43    1 1 45.68   12143   543
40        40         281     1    4      43    2 1 45.68   12143   543
41        41         282     1    5      43    3 1 45.68   12143   543
42        42         283     1    6      43    4 1 45.68   12143   543


Does this code works the same way that the merge_design.covariates? Do I also need to run the code for the S, N, etc?

So far, the models are running ok with these data.

Code: Select all

> rd.markovian_covariates.results
                                                                       model npar      AICc DeltaAICc weight  Deviance
59       S(~1)Gamma''(~time)Gamma'(~time)p(~-1 + session:time)c()N(~session)   71  786.3054  1055.254      0 1378.5068
123   S(~time)Gamma''(~time)Gamma'(~time)p(~-1 + session:time)c()N(~session)   80  792.9870  1061.936      0 1358.6123
43          S(~1)Gamma''(~time)Gamma'(~1)p(~-1 + session:time)c()N(~session)   71  807.9944  1076.943      0 1400.1957
91       S(~time)Gamma''(~1)Gamma'(~time)p(~-1 + session:time)c()N(~session)   76  822.1561  1091.105      0 1399.7603
27          S(~1)Gamma''(~1)Gamma'(~time)p(~-1 + session:time)c()N(~session)   74  826.9691  1095.918      0 1410.4615
11             S(~1)Gamma''(~1)Gamma'(~1)p(~-1 + session:time)c()N(~session)   67  852.2580  1121.207      0 1455.8468
107      S(~time)Gamma''(~time)Gamma'(~1)p(~-1 + session:time)c()N(~session)   70  862.9490  1131.898      0 1458.0210
75          S(~time)Gamma''(~1)Gamma'(~1)p(~-1 + session:time)c()N(~session)   71  863.0373  1131.986      0 1455.2387
127      S(~time)Gamma''(~time)Gamma'(~time)p(~Tourist * time)c()N(~session)   36  889.4514  1158.400      0 1573.4328
116        S(~time)Gamma''(~time)Gamma'(~time)p(~Boats * time)c()N(~session)   38  895.1463  1164.095      0 1574.3257
100           S(~time)Gamma''(~time)Gamma'(~1)p(~Boats * time)c()N(~session)   35  903.3389  1172.288      0 1589.7030
111         S(~time)Gamma''(~time)Gamma'(~1)p(~Tourist * time)c()N(~session)   34  908.8738  1177.823      0 1597.6085
63          S(~1)Gamma''(~time)Gamma'(~time)p(~Tourist * time)c()N(~session)   32  911.5907  1180.540      0 1605.0308
52            S(~1)Gamma''(~time)Gamma'(~time)p(~Boats * time)c()N(~session)   31  911.6508  1180.600      0 1607.4258
36               S(~1)Gamma''(~time)Gamma'(~1)p(~Boats * time)c()N(~session)   31  911.6508  1180.600      0 1607.4259


Thanks Jeff

[jlaake]

That works as well. Not quite the same as merge_design.covariates which uses merge to do a join in database terminology but what you did accomplishes the same thing. You only need to add the covariates to the design data which use the covariate (eg effort) in the formula. The parameters p and c can be shared and when you share them the design data for p and c are combined (rbind) and thus must have the same structure for the dataframe - which means they both must have the same covariates. Hopefully that clears things up.

--jeff