Introduction

This document demonstrates the usage of package functions in the analysis of animal tracking data. As a complete analysis consists of a series of sequential steps, it is more informative to show how the package’s functions fit into this workflow than giving separate code examples in function help files.

First we load the libraries and prepare the data.

Note When importing your own data, ctmm::as.telemetry by default will return single telemetry object with single animal, and a list of telemetry objects with multiple animals. To make code consistent we always work with a list of telemetry objects. You can use drop = FALSE parameter to make sure it’s a proper list.

library(ctmm)
library(ctmmweb)
library(magrittr)

data(buffalo)

# for importing data, use drop = FALSE to make sure it's a proper list.
# data <- as_telemetry(data_file, drop = FALSE)

# take a 100 point sample from each animal to speed up model fitting etc
data_sample <- pick(buffalo, 100)

Basic data structure

To plot the locations of multiple animals simultaneously with ggplot2, we need to merge all location data into a single data.frame. merge_tele will merge ctmm telemetry objects/lists into a list of two data.tables: data for animal location, info for some summaries of the data.

data.table has much better performance than data.frame though uses a different syntax. You can still use data.frame syntax on data.table if you are not familiar with it.

# Error will occur if single telemetry object was provided. See ?as_tele_list and ?report for details.
collected_data <- collect(buffalo)
# a list of locations data.table/data.frame and information table
loc_data <- collected_data$data
info <- collected_data$info

In loc_data:

  • identity are animal names
  • id are animal names as a factor. We will need this when we want to maintain the level information.
knitr::kable(head(loc_data))
timestamp longitude latitude t x y identity row_name id row_no
2005-07-14 05:35:00 31.88776 -24.96738 1121319300 35215.76 836.7338 Cilla 4109 Cilla 1
2005-07-14 07:35:00 31.85942 -24.94288 1121326500 32127.59 -1629.8224 Cilla 4110 Cilla 2
2005-07-14 08:34:00 31.85512 -24.94914 1121330040 32759.68 -2153.6281 Cilla 4111 Cilla 3
2005-07-14 09:35:00 31.85694 -24.95424 1121333700 33347.02 -2048.1883 Cilla 4112 Cilla 4
2005-07-14 10:35:00 31.85977 -24.94735 1121337300 32625.39 -1661.8450 Cilla 4113 Cilla 5
2005-07-14 11:34:00 31.84466 -24.93435 1121340840 30986.23 -2978.3344 Cilla 4114 Cilla 6

You can calculate distance and speed outliers in the data via:

# To save memory, data.table modify reference by default. The incoming data.table was modified and distance_center column is added after calculation. Note the telemetry list is needed for error information
assign_distance(loc_data, buffalo)
knitr::kable(head(loc_data))
timestamp longitude latitude t x y identity row_name id row_no inc_x inc_y inc_t group_index median_x median_y distance_center_raw error distance_center
2005-07-14 05:35:00 31.88776 -24.96738 1121319300 35215.76 836.7338 Cilla 4109 Cilla 1 -3088.16863 -2466.55624 7200 Cilla_0 41637.19 171.4182 6455.809 50 6455.805
2005-07-14 07:35:00 31.85942 -24.94288 1121326500 32127.59 -1629.8224 Cilla 4110 Cilla 2 632.09224 -523.80571 3540 Cilla_0 41637.19 171.4182 9678.689 50 9678.686
2005-07-14 08:34:00 31.85512 -24.94914 1121330040 32759.68 -2153.6281 Cilla 4111 Cilla 3 587.33596 105.43981 3660 Cilla_0 41637.19 171.4182 9176.930 50 9176.927
2005-07-14 09:35:00 31.85694 -24.95424 1121333700 33347.02 -2048.1883 Cilla 4112 Cilla 4 -721.62814 386.34329 3600 Cilla_0 41637.19 171.4182 8582.171 50 8582.168
2005-07-14 10:35:00 31.85977 -24.94735 1121337300 32625.39 -1661.8450 Cilla 4113 Cilla 5 -1639.16106 -1316.48941 3540 Cilla_0 41637.19 171.4182 9196.382 50 9196.380
2005-07-14 11:34:00 31.84466 -24.93435 1121340840 30986.23 -2978.3344 Cilla 4114 Cilla 6 -46.20168 -50.72027 3660 Cilla_0 41637.19 171.4182 11106.934 50 11106.932
# always calculate distance first then calculate speed outlier with distance columns added
assign_speed(loc_data, buffalo)
knitr::kable(head(loc_data))
timestamp longitude latitude t x y identity row_name id row_no inc_x inc_y inc_t group_index median_x median_y distance_center_raw error distance_center assigned_speed
2005-07-14 05:35:00 31.88776 -24.96738 1121319300 35215.76 836.7338 Cilla 4109 Cilla 1 -3088.16863 -2466.55624 7200 Cilla_0 41637.19 171.4182 6455.809 50 6455.805 0.5489290
2005-07-14 07:35:00 31.85942 -24.94288 1121326500 32127.59 -1629.8224 Cilla 4110 Cilla 2 632.09224 -523.80571 3540 Cilla_0 41637.19 171.4182 9678.689 50 9678.686 0.2318817
2005-07-14 08:34:00 31.85512 -24.94914 1121330040 32759.68 -2153.6281 Cilla 4111 Cilla 3 587.33596 105.43981 3660 Cilla_0 41637.19 171.4182 9176.930 50 9176.927 0.1630168
2005-07-14 09:35:00 31.85694 -24.95424 1121333700 33347.02 -2048.1883 Cilla 4112 Cilla 4 -721.62814 386.34329 3600 Cilla_0 41637.19 171.4182 8582.171 50 8582.168 0.1630168
2005-07-14 10:35:00 31.85977 -24.94735 1121337300 32625.39 -1661.8450 Cilla 4113 Cilla 5 -1639.16106 -1316.48941 3540 Cilla_0 41637.19 171.4182 9196.382 50 9196.380 0.5938854
2005-07-14 11:34:00 31.84466 -24.93435 1121340840 30986.23 -2978.3344 Cilla 4114 Cilla 6 -46.20168 -50.72027 3660 Cilla_0 41637.19 171.4182 11106.934 50 11106.932 0.0185420
# Or you can use pipe operator
loc_data %>% assign_distance(buffalo) %>% assign_speed(buffalo)
knitr::kable(head(loc_data))
timestamp longitude latitude t x y identity row_name id row_no inc_x inc_y inc_t group_index median_x median_y distance_center_raw error distance_center assigned_speed
2005-07-14 05:35:00 31.88776 -24.96738 1121319300 35215.76 836.7338 Cilla 4109 Cilla 1 -3088.16863 -2466.55624 7200 Cilla_0 41637.19 171.4182 6455.809 50 6455.805 0.5489290
2005-07-14 07:35:00 31.85942 -24.94288 1121326500 32127.59 -1629.8224 Cilla 4110 Cilla 2 632.09224 -523.80571 3540 Cilla_0 41637.19 171.4182 9678.689 50 9678.686 0.2318817
2005-07-14 08:34:00 31.85512 -24.94914 1121330040 32759.68 -2153.6281 Cilla 4111 Cilla 3 587.33596 105.43981 3660 Cilla_0 41637.19 171.4182 9176.930 50 9176.927 0.1630168
2005-07-14 09:35:00 31.85694 -24.95424 1121333700 33347.02 -2048.1883 Cilla 4112 Cilla 4 -721.62814 386.34329 3600 Cilla_0 41637.19 171.4182 8582.171 50 8582.168 0.1630168
2005-07-14 10:35:00 31.85977 -24.94735 1121337300 32625.39 -1661.8450 Cilla 4113 Cilla 5 -1639.16106 -1316.48941 3540 Cilla_0 41637.19 171.4182 9196.382 50 9196.380 0.5938854
2005-07-14 11:34:00 31.84466 -24.93435 1121340840 30986.23 -2978.3344 Cilla 4114 Cilla 6 -46.20168 -50.72027 3660 Cilla_0 41637.19 171.4182 11106.934 50 11106.932 0.0185420

info summarize some basic information on data.

knitr::kable(info)
identity start end interval (hours) duration (months) points calibrated
Cilla 2005-07-14 05:35 2005-12-07 22:16 1 4.97 3527 no
Gabs 2005-04-05 05:56 2005-06-27 03:45 1 2.81 1996 no
Mvubu 2005-07-15 05:02 2005-10-29 18:49 1 3.61 2572 no
Pepper 2006-04-25 05:09 2006-12-31 14:34 2 8.48 1725 no
Queen 2005-02-17 05:05 2005-06-02 03:43 1 3.55 1756 no
Toni 2005-08-23 06:35 2006-04-22 23:09 1 8.22 5766 no

Selecting a subset

In the app we can select a subset of the full data by slecting rows in the data summary table.

To perform the same selection via package functions, we can select animal names from the identity column or numbers from the id factor column in loc_data using either data.table or data.frame syntax.

We suggest to always select a subset from full data instead of creating separate data objects for individuals, because the subset will carry the id column which still holds all animal names as levels so that a consistent color mapping can be maintained (otherwise ggplot2 will always draw the first animal in same color).

# select by identity column
loc_data_sub1 <- loc_data[identity %in% c("Gabs", "Queen")]
# select by id factor column value
loc_data[as.numeric(id) %in% c(1, 3)]
#>                 timestamp longitude  latitude          t        x          y
#>    1: 2005-07-14 05:35:00  31.88776 -24.96738 1121319300 35215.76   836.7338
#>    2: 2005-07-14 07:35:00  31.85942 -24.94288 1121326500 32127.59 -1629.8224
#>    3: 2005-07-14 08:34:00  31.85512 -24.94914 1121330040 32759.68 -2153.6281
#>    4: 2005-07-14 09:35:00  31.85694 -24.95424 1121333700 33347.02 -2048.1883
#>    5: 2005-07-14 10:35:00  31.85977 -24.94735 1121337300 32625.39 -1661.8450
#>   ---                                                                       
#> 6095: 2005-10-29 14:49:00  31.91272 -24.95795 1130597340 34515.66  3474.2789
#> 6096: 2005-10-29 15:49:00  31.91805 -24.95593 1130600940 34365.50  4037.6998
#> 6097: 2005-10-29 16:49:00  31.92247 -24.95555 1130604540 34383.81  4485.3181
#> 6098: 2005-10-29 17:49:00  31.92560 -24.95899 1130608140 34805.98  4746.6588
#> 6099: 2005-10-29 18:49:00  31.92659 -24.96000 1130611740 34930.80  4830.5502
#>       identity row_name    id row_no       inc_x       inc_y inc_t group_index
#>    1:    Cilla     4109 Cilla      1 -3088.16863 -2466.55624  7200     Cilla_0
#>    2:    Cilla     4110 Cilla      2   632.09224  -523.80571  3540     Cilla_0
#>    3:    Cilla     4111 Cilla      3   587.33596   105.43981  3660     Cilla_0
#>    4:    Cilla     4112 Cilla      4  -721.62814   386.34329  3600     Cilla_0
#>    5:    Cilla     4113 Cilla      5 -1639.16106 -1316.48941  3540     Cilla_0
#>   ---                                                                         
#> 6095:    Mvubu    11148 Mvubu   8091  -150.15422   563.42090  3600     Mvubu_0
#> 6096:    Mvubu    11150 Mvubu   8092    18.30518   447.61830  3600     Mvubu_0
#> 6097:    Mvubu    11151 Mvubu   8093   422.17554   261.34066  3600     Mvubu_0
#> 6098:    Mvubu    11152 Mvubu   8094   124.81949    83.89143  3600     Mvubu_0
#> 6099:    Mvubu    11154 Mvubu   8095          NA          NA    NA     Mvubu_0
#>       median_x median_y distance_center_raw error distance_center
#>    1: 41637.19 171.4182            6455.809    50        6455.805
#>    2: 41637.19 171.4182            9678.689    50        9678.686
#>    3: 41637.19 171.4182            9176.930    50        9176.927
#>    4: 41637.19 171.4182            8582.171    50        8582.168
#>    5: 41637.19 171.4182            9196.382    50        9196.380
#>   ---                                                            
#> 6095: 40443.01 377.6464            6687.504    50        6687.500
#> 6096: 40443.01 377.6464            7094.516    50        7094.512
#> 6097: 40443.01 377.6464            7320.312    50        7320.308
#> 6098: 40443.01 377.6464            7131.928    50        7131.925
#> 6099: 40443.01 377.6464            7086.102    50        7086.098
#>       assigned_speed
#>    1:     0.54892897
#>    2:     0.23188168
#>    3:     0.16301680
#>    4:     0.16301680
#>    5:     0.59388538
#>   ---               
#> 6095:     0.21984432
#> 6096:     0.12441133
#> 6097:     0.13789396
#> 6098:     0.04168272
#> 6099:     0.04168272

Visualization

You can reproduce most of the plots in the Visualization page of the webapp with package functions, for example:

# plot animal locations
plot_loc(loc_data)

# plot a subset only. Note the color mapping is consistent because loc_data_sub1 id column hold all animal names in levels.
plot_loc(loc_data_sub1)

# with subset and full data set both provided, subset will be drawn with full data as background. 
plot_loc(loc_data_sub1, loc_data)

# location in facet
plot_loc_facet(loc_data_sub1)

# sampling time
plot_time(loc_data_sub1)

# take the ggplot2 object to further customize it
plot_loc(loc_data_sub1, loc_data) +
  ggplot2::ggtitle("Locations of Buffalos") +
  # override the default left alignment of title and make it bigger
  ctmmweb:::CENTER_TITLE

# export plot
g <- plot_loc(loc_data_sub1, loc_data)
# use this to save as file if needed
if (interactive()) {
  ggplot2::ggsave("test.png", g)
}

Variogram

We can plot both empirical variograms alone and empirical variograms with fitted theoretical semi-variance functions superimposed on them.

For example, to plot empirical variograms from telemetry data, one would do:

vario_list <- lapply(data_sample, ctmm::variogram)
# names of vario_list are needed for figure titles
names(vario_list) <- names(data_sample)
plot_vario(vario_list)

# sometimes the default figure settings doesn't work in some systems, you can clear the plot device then use a smaller title size
# dev.off()
# plot_vario(vario_list, cex = 0.55)

Similarly, we can compare initial movement model parameter guesses to the data via variograms:

guess_list <- lapply(data_sample,
                     function(tele)
                       ctmm::ctmm.guess(tele, interactive = FALSE))
plot_vario(vario_list, guess_list)

Model summary table

We can try different models on multiple animals in parallel.

Note there could be some known errors with data.table, R 3.4 and parallel operations. You can try to restart R session if met with same error.

# try multiple models in parallel with ctmm.select. parallel mode can be turned off with parallel = FALSE
model_try_res <- par_try_models(data_sample)
#> running parallel in SOCKET cluster of 6
# `model_try_res` holds a list of items named by animal names. Each item hold the attempted models for that animal as a sub list, named by model type.
print(str(model_try_res[1:3], max.level = 2))
#> List of 3
#>  $ Cilla:List of 6
#>   ..$ OU anisotropic :Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OUF anisotropic:Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OUf anisotropic:Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OU isotropic   :Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OUF isotropic  :Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ IID anisotropic:Formal class 'ctmm' [package "ctmm"]
#> @ info
#>  $ Gabs :List of 6
#>   ..$ OU anisotropic :Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OUF anisotropic:Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OU isotropic   :Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OUF isotropic  :Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OUf anisotropic:Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ IID anisotropic:Formal class 'ctmm' [package "ctmm"]
#> @ info
#>  $ Mvubu:List of 6
#>   ..$ OU anisotropic :Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OUF anisotropic:Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OUf anisotropic:Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OU isotropic   :Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ OUF isotropic  :Formal class 'ctmm' [package "ctmm"]
#> @ info
#>   ..$ IID anisotropic:Formal class 'ctmm' [package "ctmm"]
#> @ info
#> NULL
# a data.table of models information summary
model_summary <- summary_tried_models(model_try_res)
# you can also open it with RStudio's data.frame viewer
knitr::kable(model_summary)
model_no identity model_type model_name ΔAICc DOF mean DOF area DOF speed area (km²) area CI (km²) τ[position] (days) τ[position] CI (days) τ[velocity] (hours) τ[velocity] CI (hours) speed (km/day) speed CI (km/day) τ τ[decay] τ[period]
1 Cilla OU anisotropic 1. Cilla - OU anisotropic 0.00 15.35 26.43 0.00 403.54 (264.63 - 571.29) 5.07 (2.85 - 9.04) NA NA NA NA NA NA NA
2 Cilla OUF anisotropic 2. Cilla - OUF anisotropic 2.23 15.43 26.68 0.00 402.78 (264.71 - 569.37) 5.04 (2.1 - 12.07) 0.16 (0 - 13.44) 38.88 (0 - 358.56) NA NA NA
3 Cilla OUf anisotropic 3. Cilla - OUf anisotropic 13.28 30.37 50.23 62.26 363.74 (270.17 - 471.01) NA NA NA NA 5.18 (4.32 - 6.05) 107343.94 NA NA
4 Cilla OU isotropic 4. Cilla - OU isotropic 24.57 16.89 29.29 0.00 432.94 (290.58 - 603.23) 4.58 (2.64 - 7.93) NA NA NA NA NA NA NA
5 Cilla OUF isotropic 5. Cilla - OUF isotropic 26.52 17.74 30.59 0.15 431.56 (292.39 - 597.38) 4.25 (1.77 - 10.24) 2.34 (0 - 15.21) 10.37 (0 - 31.1) NA NA NA
6 Cilla IID anisotropic 6. Cilla - IID anisotropic 150.13 100.00 98.69 0.00 385.39 (313.12 - 465.05) NA NA NA NA NA NA NA NA NA
7 Gabs OU anisotropic 7. Gabs - OU anisotropic 0.00 4.79 6.55 0.00 409.20 (158.33 - 777.09) 10.93 (2.22 - 53.84) NA NA NA NA NA NA NA
8 Gabs OUF anisotropic 8. Gabs - OUF anisotropic 2.00 5.53 11.40 0.00 347.40 (175.93 - 576.23) 9.14 (3.61 - 23.16) 0.00 (0 - 0.13) NA NA NA NA NA
9 Gabs OU isotropic 9. Gabs - OU isotropic 15.42 3.93 5.30 0.00 574.22 (194.12 - 1156.77) 14.17 (1.74 - 115.47) NA NA NA NA NA NA NA
10 Gabs OUF isotropic 10. Gabs - OUF isotropic 17.31 4.58 9.89 0.00 475.35 (226.87 - 814.17) 11.57 (4.14 - 32.38) 0.00 (0 - 0.35) NA NA NA NA NA
11 Gabs OUf anisotropic 11. Gabs - OUf anisotropic 51.09 19.79 32.65 71.32 278.53 (191.31 - 381.87) NA NA NA NA 5.18 (4.32 - 6.05) 94919.23 NA NA
12 Gabs IID anisotropic 12. Gabs - IID anisotropic 306.86 100.00 120.27 0.00 277.28 (229.94 - 328.98) NA NA NA NA NA NA NA NA NA
13 Mvubu OU anisotropic 13. Mvubu - OU anisotropic 0.00 12.62 20.08 0.00 398.19 (243.51 - 590.28) 4.56 (2.31 - 9.03) NA NA NA NA NA NA NA
14 Mvubu OUF anisotropic 14. Mvubu - OUF anisotropic 0.57 14.62 22.98 3.23 394.75 (250.17 - 571.77) 3.71 (1.57 - 8.78) 4.52 (0 - 10.94) 8.64 (4.32 - 12.96) NA NA NA
15 Mvubu OUf anisotropic 15. Mvubu - OUf anisotropic 8.20 26.93 42.48 66.93 349.79 (252.61 - 462.55) NA NA NA NA 6.91 (6.05 - 6.91) 88381.02 NA NA
16 Mvubu OU isotropic 16. Mvubu - OU isotropic 26.35 14.30 24.00 0.00 418.53 (268.15 - 601.83) 3.98 (2.16 - 7.34) NA NA NA NA NA NA NA
17 Mvubu OUF isotropic 17. Mvubu - OUF isotropic 27.77 15.70 26.19 0.91 415.21 (271.68 - 588.69) 3.48 (1.53 - 7.95) 2.97 (0 - 10.03) 10.37 (1.73 - 19.87) NA NA NA
18 Mvubu IID anisotropic 18. Mvubu - IID anisotropic 167.71 100.00 99.24 0.00 353.33 (287.25 - 426.16) NA NA NA NA NA NA NA NA NA
19 Pepper OUF anisotropic 19. Pepper - OUF anisotropic 0.00 17.81 25.69 30.66 681.98 (444.21 - 969.89) 6.62 (3.2 - 13.69) 12.05 (6.13 - 23.72) 5.18 (4.32 - 6.05) NA NA NA
20 Pepper OUf anisotropic 20. Pepper - OUf anisotropic 11.03 33.35 47.33 45.09 587.49 (432.17 - 766.29) NA NA NA NA 5.18 (4.32 - 6.05) 151019.36 NA NA
21 Pepper OU anisotropic 21. Pepper - OU anisotropic 14.03 13.13 19.42 0.00 690.12 (418.06 - 1029.25) 10.01 (4.97 - 20.18) NA NA NA NA NA NA NA
22 Pepper OUF isotropic 22. Pepper - OUF isotropic 43.40 13.43 20.88 44.82 1153.83 (713.15 - 1698.83) 9.28 (4.36 - 19.74) 12.72 (7.02 - 23.05) 5.18 (4.32 - 6.05) NA NA NA
23 Queen OU anisotropic 23. Queen - OU anisotropic 0.00 10.17 15.75 0.00 293.45 (166.89 - 455.13) 4.92 (2.27 - 10.67) NA NA NA NA NA NA NA
24 Queen OUF anisotropic 24. Queen - OUF anisotropic 0.29 11.25 17.12 3.32 298.54 (174.28 - 455.71) 4.23 (1.85 - 9.66) 2.86 (0 - 6.51) 8.64 (4.32 - 13.82) NA NA NA
25 Queen OUf anisotropic 25. Queen - OUf anisotropic 20.76 24.50 36.04 58.55 266.01 (186.35 - 359.63) NA NA NA NA 6.91 (6.05 - 7.78) 76354.22 NA NA
26 Queen OUF isotropic 26. Queen - OUF isotropic 33.98 10.60 16.38 31.13 469.21 (270.21 - 722.22) 4.43 (1.81 - 10.82) 6.12 (2.54 - 14.73) 6.91 (5.18 - 7.78) NA NA NA
27 Queen OU isotropic 27. Queen - OU isotropic 42.41 8.12 12.51 0.00 471.97 (247.76 - 767.22) 6.55 (2.68 - 16.01) NA NA NA NA NA NA NA
28 Queen IID anisotropic 28. Queen - IID anisotropic 226.94 100.00 70.40 0.00 254.93 (198.88 - 317.83) NA NA NA NA NA NA NA NA NA
29 Toni OUf anisotropic 29. Toni - OUf anisotropic 0.00 34.36 51.63 53.37 322.98 (240.95 - 416.85) NA NA NA NA 3.46 (2.59 - 3.46) 156130.99 NA NA
30 Toni OU anisotropic 30. Toni - OU anisotropic 0.05 21.80 35.35 0.00 341.21 (238.14 - 462.54) 5.75 (3.38 - 9.77) NA NA NA NA NA NA NA
31 Toni OUF anisotropic 31. Toni - OUF anisotropic 0.37 25.88 39.80 4.85 337.96 (241.22 - 450.75) 4.27 (1.56 - 11.69) 13.79 (0 - 32.43) 3.46 (2.59 - 5.18) NA NA NA
32 Toni OUF isotropic 32. Toni - OUF isotropic 0.44 23.77 38.78 2.59 358.13 (254.39 - 479.33) 4.83 (1.95 - 11.94) 11.18 (0 - 29.59) 4.32 (1.73 - 6.91) NA NA NA
33 Toni OUf isotropic 33. Toni - OUf isotropic 1.51 33.49 52.31 60.73 342.06 (255.7 - 440.78) NA NA NA NA 3.46 (2.59 - 3.46) 160366.84 NA NA
34 Toni OUΩ anisotropic 34. Toni - OUΩ anisotropic 2.24 34.36 51.63 53.37 322.98 (240.95 - 416.85) NA NA NA NA 3.46 (2.59 - 3.46) NA 156131 98338578462
35 Toni IID anisotropic 35. Toni - IID anisotropic 88.06 100.00 97.61 0.00 309.00 (250.75 - 373.24) NA NA NA NA NA NA NA NA NA

There could be multiple models attempted for each animal. In the webapp you can select a subset of models then check their variograms, home ranges and occurrences. You can also select a subset of fitted models via in a script as:

To make selecting a subset easier, we first convert this nested list structure to a flat list:

# the nested structure of model fit result
names(model_try_res)
#> [1] "Cilla"  "Gabs"   "Mvubu"  "Pepper" "Queen"  "Toni"
names(model_try_res[[1]])
#> [1] "OU anisotropic"  "OUF anisotropic" "OUf anisotropic" "OU isotropic"   
#> [5] "OUF isotropic"   "IID anisotropic"
# convert to a flat list
model_list <- flatten_models(model_try_res)
names(model_list)
#>  [1] "1. Cilla - OU anisotropic"    "2. Cilla - OUF anisotropic"  
#>  [3] "3. Cilla - OUf anisotropic"   "4. Cilla - OU isotropic"     
#>  [5] "5. Cilla - OUF isotropic"     "6. Cilla - IID anisotropic"  
#>  [7] "7. Gabs - OU anisotropic"     "8. Gabs - OUF anisotropic"   
#>  [9] "9. Gabs - OU isotropic"       "10. Gabs - OUF isotropic"    
#> [11] "11. Gabs - OUf anisotropic"   "12. Gabs - IID anisotropic"  
#> [13] "13. Mvubu - OU anisotropic"   "14. Mvubu - OUF anisotropic" 
#> [15] "15. Mvubu - OUf anisotropic"  "16. Mvubu - OU isotropic"    
#> [17] "17. Mvubu - OUF isotropic"    "18. Mvubu - IID anisotropic" 
#> [19] "19. Pepper - OUF anisotropic" "20. Pepper - OUf anisotropic"
#> [21] "21. Pepper - OU anisotropic"  "22. Pepper - OUF isotropic"  
#> [23] "23. Queen - OU anisotropic"   "24. Queen - OUF anisotropic" 
#> [25] "25. Queen - OUf anisotropic"  "26. Queen - OUF isotropic"   
#> [27] "27. Queen - OU isotropic"     "28. Queen - IID anisotropic" 
#> [29] "29. Toni - OUf anisotropic"   "30. Toni - OU anisotropic"   
#> [31] "31. Toni - OUF anisotropic"   "32. Toni - OUF isotropic"    
#> [33] "33. Toni - OUf isotropic"     "34. Toni - OU惟 anisotropic"  
#> [35] "35. Toni - IID anisotropic"

Then we can find the model names in the model summary table by model_no or model_type. The code here uses data.table syntax, but you can also use the table as data.frame if you want.

# select subset in model summary table by model_no
knitr::kable(model_summary[model_no %in% c(1, 3, 10, 11, 12, 13)])
model_no identity model_type model_name ΔAICc DOF mean DOF area DOF speed area (km²) area CI (km²) τ[position] (days) τ[position] CI (days) τ[velocity] (hours) τ[velocity] CI (hours) speed (km/day) speed CI (km/day) τ τ[decay] τ[period]
1 Cilla OU anisotropic 1. Cilla - OU anisotropic 0.00 15.35 26.43 0.00 403.54 (264.63 - 571.29) 5.07 (2.85 - 9.04) NA NA NA NA NA NA NA
3 Cilla OUf anisotropic 3. Cilla - OUf anisotropic 13.28 30.37 50.23 62.26 363.74 (270.17 - 471.01) NA NA NA NA 5.18 (4.32 - 6.05) 107343.94 NA NA
10 Gabs OUF isotropic 10. Gabs - OUF isotropic 17.31 4.58 9.89 0.00 475.35 (226.87 - 814.17) 11.57 (4.14 - 32.38) 0 (0 - 0.35) NA NA NA NA NA
11 Gabs OUf anisotropic 11. Gabs - OUf anisotropic 51.09 19.79 32.65 71.32 278.53 (191.31 - 381.87) NA NA NA NA 5.18 (4.32 - 6.05) 94919.23 NA NA
12 Gabs IID anisotropic 12. Gabs - IID anisotropic 306.86 100.00 120.27 0.00 277.28 (229.94 - 328.98) NA NA NA NA NA NA NA NA NA
13 Mvubu OU anisotropic 13. Mvubu - OU anisotropic 0.00 12.62 20.08 0.00 398.19 (243.51 - 590.28) 4.56 (2.31 - 9.03) NA NA NA NA NA NA NA
# select by model type
knitr::kable(model_summary[model_type == "OU anisotropic"])
model_no identity model_type model_name ΔAICc DOF mean DOF area DOF speed area (km²) area CI (km²) τ[position] (days) τ[position] CI (days) τ[velocity] (hours) τ[velocity] CI (hours) speed (km/day) speed CI (km/day) τ τ[decay] τ[period]
1 Cilla OU anisotropic 1. Cilla - OU anisotropic 0.00 15.35 26.43 0 403.54 (264.63 - 571.29) 5.07 (2.85 - 9.04) NA NA NA NA NA NA NA
7 Gabs OU anisotropic 7. Gabs - OU anisotropic 0.00 4.79 6.55 0 409.20 (158.33 - 777.09) 10.93 (2.22 - 53.84) NA NA NA NA NA NA NA
13 Mvubu OU anisotropic 13. Mvubu - OU anisotropic 0.00 12.62 20.08 0 398.19 (243.51 - 590.28) 4.56 (2.31 - 9.03) NA NA NA NA NA NA NA
21 Pepper OU anisotropic 21. Pepper - OU anisotropic 14.03 13.13 19.42 0 690.12 (418.06 - 1029.25) 10.01 (4.97 - 20.18) NA NA NA NA NA NA NA
23 Queen OU anisotropic 23. Queen - OU anisotropic 0.00 10.17 15.75 0 293.45 (166.89 - 455.13) 4.92 (2.27 - 10.67) NA NA NA NA NA NA NA
30 Toni OU anisotropic 30. Toni - OU anisotropic 0.05 21.80 35.35 0 341.21 (238.14 - 462.54) 5.75 (3.38 - 9.77) NA NA NA NA NA NA NA
# select first(best) model for each animal using the smallest AICc value
knitr::kable(model_summary[`ΔAICc` == 0])
model_no identity model_type model_name ΔAICc DOF mean DOF area DOF speed area (km²) area CI (km²) τ[position] (days) τ[position] CI (days) τ[velocity] (hours) τ[velocity] CI (hours) speed (km/day) speed CI (km/day) τ τ[decay] τ[period]
1 Cilla OU anisotropic 1. Cilla - OU anisotropic 0 15.35 26.43 0.00 403.54 (264.63 - 571.29) 5.07 (2.85 - 9.04) NA NA NA NA NA NA NA
7 Gabs OU anisotropic 7. Gabs - OU anisotropic 0 4.79 6.55 0.00 409.20 (158.33 - 777.09) 10.93 (2.22 - 53.84) NA NA NA NA NA NA NA
13 Mvubu OU anisotropic 13. Mvubu - OU anisotropic 0 12.62 20.08 0.00 398.19 (243.51 - 590.28) 4.56 (2.31 - 9.03) NA NA NA NA NA NA NA
19 Pepper OUF anisotropic 19. Pepper - OUF anisotropic 0 17.81 25.69 30.66 681.98 (444.21 - 969.89) 6.62 (3.2 - 13.69) 12.05 (6.13 - 23.72) 5.18 (4.32 - 6.05) NA NA NA
23 Queen OU anisotropic 23. Queen - OU anisotropic 0 10.17 15.75 0.00 293.45 (166.89 - 455.13) 4.92 (2.27 - 10.67) NA NA NA NA NA NA NA
29 Toni OUf anisotropic 29. Toni - OUf anisotropic 0 34.36 51.63 53.37 322.98 (240.95 - 416.85) NA NA NA NA 3.46 (2.59 - 3.46) 156131 NA NA
# The expression is enclosed with () to enable automatical printing of result. Both model_name and identity are selected in same filter. We need the model_name to filter the model list, and the animal name to filter the variograms
(names_sub2 <- model_summary[(model_no %in% c(1, 3, 10, 11, 12, 13)),
                             .(model_name, identity)])
#>                    model_name identity
#> 1:  1. Cilla - OU anisotropic    Cilla
#> 2: 3. Cilla - OUf anisotropic    Cilla
#> 3:   10. Gabs - OUF isotropic     Gabs
#> 4: 11. Gabs - OUf anisotropic     Gabs
#> 5: 12. Gabs - IID anisotropic     Gabs
#> 6: 13. Mvubu - OU anisotropic    Mvubu

Once you have selected the models in summary table, you can filter the actual models list with model names. Note this is a different subset from loc_data_sub1 above.

# filter model list by model names to get subset of model list.
model_list_sub2 <- model_list[names_sub2$model_name]

Now we can plot variograms with the fitted SVFs of the selected models. Note the vario_list_sub2 needs to match with model_list_sub2 in length and animal name, so they are based on same data.

# get corresponding variograms by animal names.
vario_list_sub2 <- vario_list[names_sub2$identity]
# specify a different color for model
plot_vario(vario_list_sub2, model_list_sub2, model_color = "purple")

Home range

We can estimate home range with the telemetry data and the fitted models. Note parallel mode is not used because we want to calculate all animals together to put them in same grid.

# calculate home range with ctmm::akde. 
tele_list_sub2 <- data_sample[names_sub2$identity]
hrange_list_sub2 <- akde(tele_list_sub2, CTMM = model_list_sub2)
# name by model name
names(hrange_list_sub2) <- names(model_list_sub2)
# summary of each home range. There is no summary table function here because we borrowed the model table in app to make the home range summay table. To reproduce that in functions need model table/model_try_res and the selection as parameters, which will be quite awkward. If there is a strong request from users, a summary table function can be added.
lapply(hrange_list_sub2, summary)
#> $`1. Cilla - OU anisotropic`
#> $`1. Cilla - OU anisotropic`$DOF
#>      area bandwidth 
#>  26.43320  42.92138 
#> 
#> $`1. Cilla - OU anisotropic`$CI
#>                               low      est     high
#> area (square kilometers) 238.1235 363.1199 514.0671
#> 
#> 
#> $`3. Cilla - OUf anisotropic`
#> $`3. Cilla - OUf anisotropic`$DOF
#>      area bandwidth 
#>  50.22548  67.22177 
#> 
#> $`3. Cilla - OUf anisotropic`$CI
#>                               low      est     high
#> area (square kilometers) 248.4837 334.5437 433.2031
#> 
#> 
#> $`10. Gabs - OUF isotropic`
#> $`10. Gabs - OUF isotropic`$DOF
#>      area bandwidth 
#>  9.891746  8.309818 
#> 
#> $`10. Gabs - OUF isotropic`$CI
#>                               low      est     high
#> area (square kilometers) 209.8107 439.6118 752.9657
#> 
#> 
#> $`11. Gabs - OUf anisotropic`
#> $`11. Gabs - OUf anisotropic`$DOF
#>      area bandwidth 
#>  32.65113  40.72327 
#> 
#> $`11. Gabs - OUf anisotropic`$CI
#>                              low      est     high
#> area (square kilometers) 161.681 235.3929 322.7311
#> 
#> 
#> $`12. Gabs - IID anisotropic`
#> $`12. Gabs - IID anisotropic`$DOF
#>      area bandwidth 
#> 120.27079  99.99999 
#> 
#> $`12. Gabs - IID anisotropic`$CI
#>                              low      est     high
#> area (square kilometers) 172.078 207.5022 246.1927
#> 
#> 
#> $`13. Mvubu - OU anisotropic`
#> $`13. Mvubu - OU anisotropic`$DOF
#>      area bandwidth 
#>  20.08493  32.98463 
#> 
#> $`13. Mvubu - OU anisotropic`$CI
#>                               low      est     high
#> area (square kilometers) 206.9929 338.4765 501.7671
# plot home range
plot_ud(hrange_list_sub2)

# plot home range with location overlay
plot_ud(hrange_list_sub2, tele_list = tele_list_sub2)

# plot with different level.UD values
plot_ud(hrange_list_sub2, level_vec = c(0.50, 0.95), tele_list = tele_list_sub2)

Batch script for home range

Home range is one of the most used feature in ctmm/ctmmweb. Sometimes you have lots of individuals and don’t need to fine tune each model, then you can calculate home range in batch mode. ctmmweb package provided many functions to make this process as easy as possible. Check the function documents for more details.

The script can also be used to create a reproducible example if you met any error in modeling/home range parts.

library(ctmm)
library(ctmmweb)
library(data.table)

# this read from a .rds from saved progress file. You can also import your data directly with as.telemetry
tele_list <- readRDS("/Users/xhdong/Downloads/Saved_2019-11-15_11-37-41_models/input_telemetry.rds")
# one line to try models. The result is a nested list of CTMM models
model_try_res <- par_try_models(tele_list)
save.image("test.RData")
# summary models to find the best
model_summary <- summary_tried_models(model_try_res)
# select the best model name for each individual, which is the first one with 0 AICc
best_model_names <- model_summary[, .(best_model_name = model_name[1]),
                                  by = "identity"]$best_model_name
# get a flat list of model object
model_list <- flatten_models(model_try_res)
# get the best model objects
best_models <- model_list[best_model_names]
# calculate home range in same grid, with optimal weighting on
hrange_list_same_grid <- akde(tele_list, best_models,
                              weights = TRUE)
# calculate home range separately, with optimal weighting on
hrange_list <- par_hrange_each(tele_list, best_models,
                               rep.int(list(TRUE), length(tele_list)))

Occurrence

Occurrence distribution can be calculated from the telemetry data and the fitted models in parallel.

# calculate occurrence in parallel. parallel can be disabled with fallback = TRUE
occur_list_sub2 <- par_occur(tele_list_sub2, model_list_sub2)
#> running parallel in SOCKET cluster of 6
# plot occurrence. Note tele_list is not needed here because the location overlay usually interfere with occurrence plot.
plot_ud(occur_list_sub2, level_vec = c(0.50, 0.95))

Maps

We can then build interactive online maps of the data and fitted Home Range and Occurrence distributions. For example, one can display the location data on a map via:

# loc_data_sub1 is used for visualization plot. all the model selections above are based on a different subset. We need to get corresponding loc_data subset for model selections
loc_data_sub2 <- loc_data[identity %in% names(tele_list_sub2)]
point_map(loc_data_sub2)
point_map(loc_data)