Getting Started with soReta

What soReta does

soReta takes your camera-trap dataset and returns datasets ready for statistical analysis: GLMM, GAMM, occupancy, kernel density and circular activity analysis, temporal interactions between species, capture-mark-recapture -- without you having to write the aggregation code yourself.

To do this, it needs a recordTable featuring a station (Station), a date/time (DateTimeOriginal), and a species (Species) column formatted in camtrapR's standard style (Niedballa et al. 2016), along with a matching camOp effort matrix. Additionally, some functions require an individual-count column; unlike the core variables, this column does not need to follow camtrapR formatting.

This vignette walks through the main families of functions using soReta's own bundled example data: camOp_soReta (the effort matrix, already built via camtrapR::cameraOperation() from camtraps_soReta, a station table -- see ?camtraps_soReta for its own format), recordTable_soReta (species detections), and recordTableIndividuals_soReta (individually-identifiable detections, for the capture-mark-recapture section) -- all entirely synthetic; see ?camOp_soReta, ?recordTable_soReta, and ?recordTableIndividuals_soReta for details on how they were built.

library(soReta)

head(camOp_soReta[, 1:6])   # only the first 6 of 120 days, or the table gets too wide
#>      2026-01-01 2026-01-02 2026-01-03 2026-01-04 2026-01-05 2026-01-06
#> S_01        0.5        1.0        1.0        1.0        1.0          1
#> S_02         NA        0.5        1.0        1.0        1.0          1
#> S_03         NA         NA        0.5        1.0        1.0          1
#> S_04         NA         NA         NA        0.5        1.0          1
#> S_05         NA         NA         NA         NA        0.5          1
head(recordTable_soReta)
#>   Station    DateTimeOriginal  Species N_individuals
#> 1    S_01 2026-01-10 16:59:57  red fox             1
#> 2    S_01 2026-01-17 19:18:45     wolf             1
#> 3    S_01 2026-01-18 20:03:14  red fox             1
#> 4    S_01 2026-01-20 22:42:13 red deer             1
#> 5    S_01 2026-01-20 22:56:47 red deer             1
#> 6    S_01 2026-01-21 00:13:52 red deer             2
head(recordTableIndividuals_soReta)
#>   Station  Species Individual    DateTimeOriginal
#> 1    S_05 red deer      RD_01 2026-01-08 07:32:36
#> 2    S_02 red deer      RD_01 2026-01-20 23:01:07
#> 3    S_01 red deer      RD_01 2026-01-21 00:13:52
#> 4    S_01 red deer      RD_01 2026-01-21 00:20:58
#> 5    S_03 red deer      RD_01 2026-01-24 04:18:29
#> 6    S_02 red deer      RD_01 2026-01-27 07:46:41

1. GLMM/GAMM-ready datasets: the build_site_* and build_day_total family

These functions aggregate independent events into counts (N events) and Relative Abundance Index (RAI, O'Brien et al. 2003), at whatever temporal grain you need. Thay all use a configurable threshold, in minutes (threshold_min), to determine which events count as independent. All of them recompute event independence from the raw timestamps (never relying on a pre-set threshold baked into the recordTable).

# 1 row = 1 site x 1 day (no RAI column here: with n_days_active always
# equal to 1, RAI would just be N x 100 -- no extra information over N)
ds_day <- build_site_day(recordTable_soReta, camOp_soReta, threshold_min = 30)
head(ds_day[ds_day$N_sp > 0, ])
#>    Station       Date N_sp red_deer_N red_fox_N wild_boar_N wolf_N
#> 10    S_01 2026-01-10    1          0         1           0      0
#> 17    S_01 2026-01-17    1          0         0           0      1
#> 18    S_01 2026-01-18    1          0         1           0      0
#> 20    S_01 2026-01-20    1          1         0           0      0
#> 21    S_01 2026-01-21    1          1         0           0      0
#> 22    S_01 2026-01-22    1          1         0           0      0

Choosing how events are grouped: independence_method and require_uninterrupted

Before going further, it's worth understanding exactly how threshold_min turns raw timestamps into "independent events" -- because two further arguments (independence_method, require_uninterrupted), present on build_site_day() and on every other function in this package that takes a threshold_min, let you control that precisely. "chain"/FALSE are set as the defaults. This choice is usually unsettled in the published literature, though, I encourage you to report which independence_method/require_uninterrupted combination you used if you publish results built with soReta, the same way you'd report threshold_min itself.

independence_method ("chain", default, or "window"): "chain" compares each photo only to the one immediately before it -- a close-together run of photos can, in principle, span far longer than threshold_min, as long as no single gap between consecutive photos exceeds it. "window" instead compares each photo to the start of the current bout: a bout can never last longer than threshold_min. The two methods only disagree on bouts long enough to contain more than one sub-threshold gap in a row -- illustrated below on a real sequence from recordTable_soReta itself (a wolf pack at station S_02, three photos 16 and 25 minutes apart -- "chain" counts this as one event, "window" as two, because the first 16-minute gap already used up most of the 30-minute budget before the second photo arrived):

wolf_burst <- recordTable_soReta[
  recordTable_soReta$Station == "S_02" &
    recordTable_soReta$Species == "wolf" &
    format(recordTable_soReta$DateTimeOriginal, "%Y-%m-%d") == "2026-02-14",
]
wolf_burst
#>     Station    DateTimeOriginal Species N_individuals
#> 102    S_02 2026-02-14 21:10:00    wolf             2
#> 103    S_02 2026-02-14 21:26:00    wolf             1
#> 104    S_02 2026-02-14 21:51:00    wolf             6
wolf_burst_camOp <- camOp_soReta["S_02", "2026-02-14", drop = FALSE]

build_site_day(wolf_burst, wolf_burst_camOp, threshold_min = 30, independence_method = "chain")$wolf_N
#> [1] 1
build_site_day(wolf_burst, wolf_burst_camOp, threshold_min = 30, independence_method = "window")$wolf_N
#> [1] 2

require_uninterrupted (FALSE, default, or TRUE): implements the third independence criterion of O'Brien et al. (2003) -- the paper most commonly cited for the whole convention. Two photos of the same species are only subject to the time threshold at all if they are literally consecutive in the station's raw, all-species timeline: if a photo of any other species falls chronologically between them, the two are independent regardless of elapsed time. Another real sequence from recordTable_soReta shows this directly (station S_04: a wolf, then a red fox 10 minutes later, then a wolf again 10 minutes after that):

interrupted_burst <- recordTable_soReta[
  recordTable_soReta$Station == "S_04" &
    format(recordTable_soReta$DateTimeOriginal, "%Y-%m-%d") == "2026-03-10",
]
interrupted_burst
#>     Station    DateTimeOriginal Species N_individuals
#> 285    S_04 2026-03-10 22:15:00    wolf             3
#> 286    S_04 2026-03-10 22:25:00 red fox             1
#> 287    S_04 2026-03-10 22:35:00    wolf             5

interrupted_burst_camOp <- camOp_soReta["S_04", "2026-03-10", drop = FALSE]

build_site_day(interrupted_burst, interrupted_burst_camOp, threshold_min = 30, require_uninterrupted = FALSE)$wolf_N
#> [1] 1
build_site_day(interrupted_burst, interrupted_burst_camOp, threshold_min = 30, require_uninterrupted = TRUE)$wolf_N
#> [1] 2

With require_uninterrupted = FALSE the two wolf photos are 20 minutes apart (below the 30-minute threshold) and collapse into one event. With require_uninterrupted = TRUE the red fox in between breaks the run, so both wolf photos count as independent regardless of the 20-minute gap.

independence_method and require_uninterrupted combine freely, giving four possible conventions in total. Every function below that takes a threshold_min supports both -- from here on, this vignette just points back to this section rather than re-explaining them each time.

O'Brien et al. (2003) actually define three criteria for independence: (1) consecutive photos of different individuals of the same or different species, (2) consecutive photos of the same species more than 0.5 hours apart, and (3) non-consecutive photos of the same species. independence_method implements criterion (2) (as "chain" or "window"), and require_uninterrupted implements criterion (3). Criterion (1) is deliberately not implemented: individual recognition is only feasible for some species (e.g. those with strong coat-pattern variation or sexual dimorphism) and not others, and applying it selectively -- crediting extra independent events only where individuals happen to be distinguishable -- would bias N and RAI inconsistently across species within the same study. Since most studies now monitor several species at once rather than a single one, this package does not implement criterion (1) for any species.

The rest of the build_site_* family

The other functions in the build_site_* family perform the same task as build_site_day(), but aggregate events at a different temporal grain:

All three also take min_days (default NULL, no filtering): set it to drop any site x interval combination where the station was active for fewer than that many days -- useful for excluding a short, unreliable block or month from the dataset rather than keeping it in with an unreliably low sample size.

# 1 row = 1 site x 1 fixed N-day block
ds_week <- build_site_block(recordTable_soReta, camOp_soReta, block_days = 7, threshold_min = 30)
head(ds_week)
#>   Station block_start  block_end    mid_day n_days_active N_sp red_deer_N
#> 1    S_01  2026-01-01 2026-01-07 2026-01-04             7    0          0
#> 2    S_01  2026-01-08 2026-01-14 2026-01-11             7    1          0
#> 3    S_01  2026-01-15 2026-01-21 2026-01-18             7    3          2
#> 4    S_01  2026-01-22 2026-01-28 2026-01-25             7    2          4
#> 5    S_01  2026-01-29 2026-02-04 2026-02-01             7    1          2
#> 6    S_01  2026-02-05 2026-02-11 2026-02-08             7    2          2
#>   red_deer_RAI red_fox_N red_fox_RAI wild_boar_N wild_boar_RAI wolf_N wolf_RAI
#> 1         0.00         0        0.00           0             0      0     0.00
#> 2         0.00         1       14.29           0             0      0     0.00
#> 3        28.57         1       14.29           0             0      1    14.29
#> 4        57.14         0        0.00           0             0      1    14.29
#> 5        28.57         0        0.00           0             0      0     0.00
#> 6        28.57         0        0.00           0             0      1    14.29

# 1 row = 1 site x 1 calendar month
ds_month <- build_site_month(recordTable_soReta, camOp_soReta, threshold_min = 30)
head(ds_month)
#> # A tibble: 6 × 17
#>   Station  year month month_start month_end  mid_day    n_days_active
#>   <chr>   <dbl> <dbl> <date>      <date>     <date>             <int>
#> 1 S_01     2026     1 2026-01-01  2026-01-31 2026-01-16            31
#> 2 S_01     2026     2 2026-02-01  2026-02-28 2026-02-14            28
#> 3 S_01     2026     3 2026-03-01  2026-03-31 2026-03-16            16
#> 4 S_01     2026     4 2026-04-01  2026-04-30 2026-04-15            30
#> 5 S_02     2026     1 2026-01-01  2026-01-31 2026-01-16            30
#> 6 S_02     2026     2 2026-02-01  2026-02-28 2026-02-14            28
#> # ℹ 10 more variables: n_days_in_month <int>, N_sp <dbl>, red_deer_N <int>,
#> #   red_deer_RAI <dbl>, red_fox_N <int>, red_fox_RAI <dbl>, wild_boar_N <int>,
#> #   wild_boar_RAI <dbl>, wolf_N <int>, wolf_RAI <dbl>

# 1 row = 1 site x 1 named period, recurring every year -- our example
# data only spans January to April 2026, so a two-period split fits
# better here than a full four-season year
ds_period <- build_site_period(
  recordTable_soReta, camOp_soReta,
  period_names  = c("early", "late"),
  period_starts = c("01/01/2026", "01/03/2026"),
  threshold_min = 30
)
head(ds_period)
#>   Station period_name year period_start period_end    mid_day n_days_active
#> 1    S_01       early 2026        01/01      28/02 2026-01-30            59
#> 2    S_01        late 2026        01/03      31/12 2026-07-31            46
#> 3    S_02       early 2026        01/01      28/02 2026-01-30            58
#> 4    S_02        late 2026        01/03      31/12 2026-07-31            61
#> 5    S_03       early 2026        01/01      28/02 2026-01-30            37
#> 6    S_03        late 2026        01/03      31/12 2026-07-31            21
#>   period_length_days N_sp red_deer_N red_deer_RAI red_fox_N red_fox_RAI
#> 1                 59    4         13        22.03         6       10.17
#> 2                306    4          3         6.52        11       23.91
#> 3                 59    4         16        27.59         4        6.90
#> 4                306    4         10        16.39        16       26.23
#> 5                 59    4          5        13.51         5       13.51
#> 6                306    4          5        23.81         2        9.52
#>   wild_boar_N wild_boar_RAI wolf_N wolf_RAI
#> 1           2          3.39      7    11.86
#> 2           7         15.22      4     8.70
#> 3           5          8.62     10    17.24
#> 4           5          8.20      8    13.11
#> 5           4         10.81      7    18.92
#> 6           5         23.81      8    38.10

build_site_block(), build_site_month(), and build_site_period() all include also a mid_day column: the day roughly halfway between the start and end of each block/month/period, rounded down when the span is even so it always lands on a real day. It's there to make joining an external daily covariate straightforward -- lunar illumination, average temperature, whatever your analysis needs -- without you having to compute a representative date yourself:

# example: joining a per-day covariate you already have, e.g. lunar fraction
ds_month |> dplyr::left_join(my_lunar_fraction, by = c("mid_day" = "Date"))

Station-level covariates (elevation, habitat type, distance to the nearest road...) join just as directly -- every function in this package returns a Station column, so a single left_join() is all it takes:

ds_month |> dplyr::left_join(my_station_covariates, by = "Station")

To exclude species you're not interested in (domestic animals, unclear identification, etc.) from every downstream calculation, filter recordTable before calling any build_*() function -- RAI and richness are recomputed correctly from whatever species remain:

recordTable_soReta |> dplyr::filter(!Species %in% c("...")) |> build_site_month(camOp_soReta, threshold_min = 30)

Finally, two related functions collapse a different dimension: build_site_total() gives one row per site over the whole survey -- handy for data exploration, community analyses (e.g. vegan::vegdist()), and modelling -- while build_day_total() gives one row per calendar day, summed across all stations. build_day_total() also returns the number of active stations that day and RAI values, calculated as the sum of events for each species divided by the number of active stations that day, times 100, analogous to the classic RAI proposed by O'Brien et al. (2003). This shape is specifically suited to GLMM/GAMM analyses of behaviour in response to day-to-day environmental variation, such as the lunar cycle, and can be readily further aggregated for tests such as the chi-square test. Both functions also take independence_method/require_uninterrupted, as explained above.

# 1 row = 1 site
ds_site_tot <- build_site_total(recordTable_soReta, camOp_soReta, threshold_min = 30)
head(ds_site_tot)
#> # A tibble: 5 × 11
#>   Station n_days_active  N_sp red_deer_N red_deer_RAI red_fox_N red_fox_RAI
#>   <chr>           <int> <dbl>      <int>        <dbl>     <int>       <dbl>
#> 1 S_01              105     4         16         15.2        17       16.2 
#> 2 S_02              119     4         26         21.8        20       16.8 
#> 3 S_03               58     4         10         17.2         7       12.1 
#> 4 S_04               94     4         16         17.0         4        4.26
#> 5 S_05              106     4         11         10.4        20       18.9 
#> # ℹ 4 more variables: wild_boar_N <int>, wild_boar_RAI <dbl>, wolf_N <int>,
#> #   wolf_RAI <dbl>

# 1 row = 1 day
ds_day_tot <- build_day_total(recordTable_soReta, camOp_soReta, threshold_min = 30)
head(ds_day_tot)
#> # A tibble: 6 × 11
#>   Date       n_stations_active  N_sp red_deer_N red_deer_RAI red_fox_N
#>   <date>                 <int> <dbl>      <int>        <dbl>     <int>
#> 1 2026-01-01                 1     0          0            0         0
#> 2 2026-01-02                 2     1          1           50         0
#> 3 2026-01-03                 3     0          0            0         0
#> 4 2026-01-04                 4     1          1           25         0
#> 5 2026-01-05                 5     2          0            0         0
#> 6 2026-01-06                 5     0          0            0         0
#> # ℹ 5 more variables: red_fox_RAI <dbl>, wild_boar_N <int>,
#> #   wild_boar_RAI <dbl>, wolf_N <int>, wolf_RAI <dbl>

2. Group size

For species where the number of animals per photo matters -- herd size, pack size, flock size -- the build_group_size_* family mirrors the build_site_* family, but tracks group size instead of just counting events. All seven functions take a countCol argument: the name of the column in recordTable holding the number of animals per photo. In recordTable_soReta that column is called "N_individuals". They also support independence_method and require_uninterrupted, exactly as explained in Section 1: the boundary of a "bout" is defined identically here and in build_site_*(), so the event counts and the group sizes always agree on where one passage ends and the next begins. countCol is only needed for this family -- every other function in the package works without it.

build_group_size_events() is the base of the family: one row per independent event, with no temporal aggregation at all, not even daily. It has no counterpart in build_site_*() -- there, an event is always worth exactly 1 towards the count, so a raw, unaggregated listing of events wouldn't add any information beyond build_site_day() itself. Here, the raw group size per event is itself the finest-grained information there is, before any averaging or summing collapses it.

Within a single bout, group_size is always the maximum count across the photos in that bout, never the sum: summing would treat repeated photos of the same passing group as if they were separate individuals.

gr_size_ev <- build_group_size_events(recordTable_soReta, camOp_soReta, countCol = "N_individuals",
                         independence_method = "window", threshold_min = 30)
head(gr_size_ev)
#> # A tibble: 6 × 6
#>   Station Species  DateTime            Date       Time     group_size
#>   <chr>   <chr>    <dttm>              <date>     <chr>         <dbl>
#> 1 S_01    red deer 2026-01-20 22:42:13 2026-01-20 22:42:13          1
#> 2 S_01    red deer 2026-01-21 00:13:52 2026-01-21 00:13:52          2
#> 3 S_01    red deer 2026-01-22 19:13:45 2026-01-22 19:13:45          2
#> 4 S_01    red deer 2026-01-28 07:27:29 2026-01-28 07:27:29          2
#> 5 S_01    red deer 2026-01-28 08:39:37 2026-01-28 08:39:37          1
#> 6 S_01    red deer 2026-01-28 09:14:28 2026-01-28 09:14:28          1

build_group_size_day(), build_group_size_block(), build_group_size_month(), build_group_size_period(), build_group_size_total(), and build_group_day_total() aggregate build_group_size_events() at the same temporal grains as the build_site_* family, each producing three columns per species: <species>_mean_group_size, <species>_max_group_size, and <species>_sum_group_size -- plus n_days_active/n_stations_active, N_sp, and <species>_RAI_individuals (= sum_group_size / n_days_active * 100, the individual-based counterpart of the event-based RAI in build_site_*()) on every level except the daily one, where they would be redundant. Unlike build_site_*(), a species with no events in a given day/period gets NA in the group-size columns, not 0: a group size of zero never happens, and using 0 would silently distort any later averaging across days or periods.

# most site x day x species combinations are genuinely empty
# at this grain, so filter for the informative rows
ds_group_day <- build_group_size_day(recordTable_soReta, camOp_soReta,
                                       countCol = "N_individuals",
                                       threshold_min = 30)
ds_group_day[!is.na(ds_group_day$red_deer_mean_group_size), ] |>
  head() |>
  print(width = 90)
#> # A tibble: 6 × 15
#>   Station Date        N_sp red_deer_mean_group_size red_deer_max_group_size
#>   <chr>   <date>     <dbl>                    <dbl>                   <dbl>
#> 1 S_01    2026-01-20     1                     1                          1
#> 2 S_01    2026-01-21     1                     2                          2
#> 3 S_01    2026-01-22     1                     2                          2
#> 4 S_01    2026-01-28     1                     1.33                       2
#> 5 S_01    2026-02-01     1                     2.5                        4
#> 6 S_01    2026-02-11     1                     1.5                        2
#>   red_deer_sum_group_size red_fox_mean_group_size red_fox_max_group_size
#>                     <dbl>                   <dbl>                  <dbl>
#> 1                       1                      NA                     NA
#> 2                       2                      NA                     NA
#> 3                       2                      NA                     NA
#> 4                       4                      NA                     NA
#> 5                       5                      NA                     NA
#> 6                       3                      NA                     NA
#> # ℹ 7 more variables: red_fox_sum_group_size <dbl>, wild_boar_mean_group_size <dbl>,
#> #   wild_boar_max_group_size <dbl>, wild_boar_sum_group_size <dbl>,
#> #   wolf_mean_group_size <dbl>, wolf_max_group_size <dbl>, wolf_sum_group_size <dbl>
ds_group_total <- build_group_size_total(recordTable_soReta, camOp_soReta,
                                           countCol = "N_individuals",
                                           threshold_min = 30)
print(ds_group_total, width = 100)
#> # A tibble: 5 × 19
#>   Station n_days_active  N_sp red_deer_mean_group_size red_deer_max_group_size
#>   <chr>           <int> <dbl>                    <dbl>                   <dbl>
#> 1 S_01              105     4                     1.56                       4
#> 2 S_02              119     4                     1.92                       4
#> 3 S_03               58     4                     1.6                        3
#> 4 S_04               94     4                     1.81                       4
#> 5 S_05              106     4                     1.45                       3
#>   red_deer_sum_group_size red_deer_RAI_individuals red_fox_mean_group_size
#>                     <dbl>                    <dbl>                   <dbl>
#> 1                      25                     23.8                    1.24
#> 2                      50                     42.0                    1.1 
#> 3                      16                     27.6                    1.14
#> 4                      29                     30.8                    1.25
#> 5                      16                     15.1                    1.15
#> # ℹ 11 more variables: red_fox_max_group_size <dbl>, red_fox_sum_group_size <dbl>,
#> #   red_fox_RAI_individuals <dbl>, wild_boar_mean_group_size <dbl>, wild_boar_max_group_size <dbl>,
#> #   wild_boar_sum_group_size <dbl>, wild_boar_RAI_individuals <dbl>, wolf_mean_group_size <dbl>,
#> #   wolf_max_group_size <dbl>, wolf_sum_group_size <dbl>, wolf_RAI_individuals <dbl>

3. Occupancy detection histories

build_occupancy_day() and build_occupancy_block() return, for each species, a site x occasion matrix of 0/1/NA -- ready for unmarked, ubms, or spOccupancy. NA marks an occasion with no active camera effort at all; 0 means the station was active but the species was not detected. build_occupancy_day() creates a daily matrix; build_occupancy_block() aggregates days over a specified interval (block_days), as several functions in the build_site_* and build_group_size_* families also do. Only build_occupancy_block() takes min_days -- and, unlike elsewhere in the package, it defaults to 1 here rather than NA: a block needs at least one active day to be trusted as a real 0 (absence) rather than left as NA (effort too low to draw any conclusion). These two functions have no independence threshold at all (presence/absence doesn't need one), so independence_method/require_uninterrupted don't apply here.

occ_day <- build_occupancy_day(recordTable_soReta, camOp_soReta)
head(occ_day[["wolf"]])
#> # A tibble: 5 × 121
#>   Station `2026-01-01` `2026-01-02` `2026-01-03` `2026-01-04` `2026-01-05`
#>   <chr>          <int>        <int>        <int>        <int>        <int>
#> 1 S_01               0            0            0            0            0
#> 2 S_02              NA            0            0            0            0
#> 3 S_03              NA           NA            0            0            0
#> 4 S_04              NA           NA           NA            0            1
#> 5 S_05              NA           NA           NA           NA            0
#> # ℹ 115 more variables: `2026-01-06` <int>, `2026-01-07` <int>,
#> #   `2026-01-08` <int>, `2026-01-09` <int>, `2026-01-10` <int>,
#> #   `2026-01-11` <int>, `2026-01-12` <int>, `2026-01-13` <int>,
#> #   `2026-01-14` <int>, `2026-01-15` <int>, `2026-01-16` <int>,
#> #   `2026-01-17` <int>, `2026-01-18` <int>, `2026-01-19` <int>,
#> #   `2026-01-20` <int>, `2026-01-21` <int>, `2026-01-22` <int>,
#> #   `2026-01-23` <int>, `2026-01-24` <int>, `2026-01-25` <int>, …

occ_week <- build_occupancy_block(recordTable_soReta, camOp_soReta, block_days = 7, min_days = 4)
head(occ_week[["wolf"]])
#> # A tibble: 5 × 19
#>   Station `2026-01-01--2026-01-07` 2026-01-08--2026-01-…¹ 2026-01-15--2026-01-…²
#>   <chr>                      <int>                  <int>                  <int>
#> 1 S_01                           0                      0                      1
#> 2 S_02                           0                      1                      0
#> 3 S_03                           0                      0                      0
#> 4 S_04                           1                     NA                     NA
#> 5 S_05                          NA                      1                      1
#> # ℹ abbreviated names: ¹​`2026-01-08--2026-01-14`, ²​`2026-01-15--2026-01-21`
#> # ℹ 15 more variables: `2026-01-22--2026-01-28` <int>,
#> #   `2026-01-29--2026-02-04` <int>, `2026-02-05--2026-02-11` <int>,
#> #   `2026-02-12--2026-02-18` <int>, `2026-02-19--2026-02-25` <int>,
#> #   `2026-02-26--2026-03-04` <int>, `2026-03-05--2026-03-11` <int>,
#> #   `2026-03-12--2026-03-18` <int>, `2026-03-19--2026-03-25` <int>,
#> #   `2026-03-26--2026-04-01` <int>, `2026-04-02--2026-04-08` <int>, …

4. Kernel density and circular activity analysis

extract_radians() converts detection times into radians (clock time, 0 to 2*pi) -- the format overlap::densityPlot(), overlap::overlapPlot(), and activity::fitact() expect, one vector per species. It computes independence over the entire detection history before any subsetting, so a threshold-based cut never splits one continuous visit into two. Like the functions above, it also supports independence_method/require_uninterrupted (Section 1).

rad <- extract_radians(recordTable_soReta, threshold_min = 30)
names(rad)
#> [1] "red deer"  "red fox"   "wild boar" "wolf"
# overlap::densityPlot(rad[["wolf"]], xcenter = "midnight")

With no further arguments, each species' result is one flat vector, covering the whole dataset. The group_col argument splits it further -- one nested list level per column you name, in that order -- useful if you want a separate activity-pattern estimate per station, per season, or both, without calling the function once per subgroup yourself:

recordTable_soReta$bimonth <- paste0("bim", ceiling(lubridate::month(recordTable_soReta$DateTimeOriginal) / 2))

# one level: split by station
rad_by_station <- extract_radians(recordTable_soReta, threshold_min = 30, group_col = "Station")
rad_by_station[["wolf"]][["S_01"]]
#>  [1] 5.056000677 5.199263120 1.632319183 3.655519396 3.941573272 5.895940379
#>  [7] 6.085239133 0.580103810 0.747882346 5.995642313 0.007503861

# two levels: station, then two-month period within station
rad_by_station_bimonth <- extract_radians(recordTable_soReta, threshold_min = 30,
                                            group_col = c("Station", "bimonth"))
names(rad_by_station_bimonth[["wolf"]][["S_01"]])   # check which periods actually
#> [1] "bim1" "bim2"
                                                       # exist for this species/station
                                                       # before indexing further
rad_by_station_bimonth[["wolf"]][["S_01"]][["bim1"]]
#> [1] 5.056001 5.199263 1.632319 3.655519 3.941573 5.895940 6.085239

Don't confuse group_col with the separate group_cols argument. group_cols decides which events get compared against each other when checking independence -- by default, only events at the same station and of the same species are compared, which is the right choice for almost every study. group_col, instead, only decides how the already-computed result gets organised into nested lists for you to browse -- it never changes which photos count as one event.

One situation where you might genuinely want to change group_cols: two camera traps placed a few meters apart, effectively watching the same spot (a narrow trail, a den entrance). An animal walking past could trigger both cameras within seconds of each other -- two separate "events" at two different Station values, even though it was really one single passage. Merging that pair into one combined station name before calling extract_radians(), and passing that combined name via group_cols, treats the two cameras as one:

recordTable_soReta$Cluster <- ifelse(
  recordTable_soReta$Station %in% c("S_01", "S_02"),
  "S_01_S_02_cluster",
  recordTable_soReta$Station
)

# group_cols changes what counts as independent: the two real stations
# are now merged for this purpose
rad_clustered <- extract_radians(recordTable_soReta, threshold_min = 30,
                                   group_cols = c("Cluster", "Species"))

# group_col, unchanged in meaning: still just splits the finished
# result -- here, by the same clustered column, just to display it
rad_clustered_split <- extract_radians(recordTable_soReta, threshold_min = 30,
                                         group_cols = c("Cluster", "Species"),
                                         group_col  = "Cluster")
rad_clustered_split[["wolf"]][["S_01_S_02_cluster"]]
#>  [1] 0.002108940 0.784525499 5.056000677 2.385501477 5.199263120 0.116718894
#>  [7] 1.632319183 3.655519396 3.941573272 5.541420375 1.249583022 1.778451538
#> [13] 5.895940379 6.085239133 0.086830130 0.320010432 0.633355099 4.615086875
#> [19] 0.580103810 0.747882346 0.320486084 4.940954390 5.237437829 3.908446694
#> [25] 0.042106068 5.995642313 0.007503861 5.245950677 1.431969929

One more detail worth knowing before plotting: the radians this function returns are in clock time, not solar time. For comparisons across seasons or across sites at different latitudes, convert the result with activity::solartime() or overlap::sunTime() before handing it to densityPlot()/fitact() -- extract_radians() itself does not do this conversion.

If you specifically need one separate object per species/group combination in your environment -- to reuse an older script that was already written that way -- radians_to_env() takes the (possibly nested) list extract_radians() returns and creates one object per vector, named by prefix + species + group levels:

Rad_obj <- radians_to_env(rad_by_station, prefix = "Rad_", sep = "-")
#> 20 objects created: Rad_red deer-S_01, Rad_red deer-S_02, Rad_red deer-S_03, Rad_red deer-S_04, Rad_red deer-S_05, ...
Rad_obj
#>  [1] "Rad_red deer-S_01"  "Rad_red deer-S_02"  "Rad_red deer-S_03" 
#>  [4] "Rad_red deer-S_04"  "Rad_red deer-S_05"  "Rad_red fox-S_01"  
#>  [7] "Rad_red fox-S_02"   "Rad_red fox-S_03"   "Rad_red fox-S_04"  
#> [10] "Rad_red fox-S_05"   "Rad_wild boar-S_01" "Rad_wild boar-S_02"
#> [13] "Rad_wild boar-S_03" "Rad_wild boar-S_04" "Rad_wild boar-S_05"
#> [16] "Rad_wolf-S_01"      "Rad_wolf-S_02"      "Rad_wolf-S_03"     
#> [19] "Rad_wolf-S_04"      "Rad_wolf-S_05"
str(get(Rad_obj[1]))
#>  num [1:16] 5.9438 0.0605 5.0342 1.9525 2.2673 ...

5. Temporal interactions between species

build_species_pair_intervals() computes the AB/BA time intervals between two species (Parsons et al. 2016, Niedballa et al. 2019): the time from a detection of species A to the next detection of species B at the same station, and vice versa. Intervals with no follow-up detection before the end of a station's monitoring period are not discarded -- they are flagged censored = TRUE, so you can hand them to a proper survival model (survival::survreg()) instead of losing that information. It also supports independence_method/require_uninterrupted (Section 1).

pair_data <- build_species_pair_intervals(
  recordTable_soReta, camOp_soReta,
  speciesA = "wolf", speciesB = "wild boar",
  threshold_min = 30
)
head(pair_data)
#>   Station direction           time_from delta_hours censored
#> 1    S_01        AB 2026-01-17 19:18:45   841.29917    FALSE
#> 2    S_01        AB 2026-01-27 19:51:35   600.75194    FALSE
#> 3    S_01        AB 2026-02-05 06:14:06   398.37667    FALSE
#> 4    S_01        AB 2026-02-12 13:57:47   222.64861    FALSE
#> 5    S_01        AB 2026-02-12 15:03:20   221.55597    FALSE
#> 6    S_01        BA 2026-02-21 20:36:42    49.90917    FALSE

Its output is exactly what all three of Niedballa et al. (2019)'s AB/BA-level methods need, unchanged:

# 1. linear model (log-transformed, as Niedballa et al. 2019 did to meet
#    linear model assumptions)
mod <- lm(log(delta_hours) ~ direction, data = pair_data[!pair_data$censored, ])
summary(mod)
#> 
#> Call:
#> lm(formula = log(delta_hours) ~ direction, data = pair_data[!pair_data$censored, 
#>     ])
#> 
#> Residuals:
#>     Min      1Q  Median      3Q     Max 
#> -3.4909 -0.7489  0.1749  0.7839  2.5096 
#> 
#> Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)    
#> (Intercept)  4.71018    0.15481  30.426   <2e-16 ***
#> directionBA -0.08118    0.23747  -0.342    0.733    
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 1.286 on 118 degrees of freedom
#> Multiple R-squared:  0.0009894,  Adjusted R-squared:  -0.007477 
#> F-statistic: 0.1169 on 1 and 118 DF,  p-value: 0.7331

# 2. Mann-Whitney U-test 
wilcox.test(delta_hours ~ direction, data = pair_data[!pair_data$censored, ])
#> 
#>  Wilcoxon rank sum test with continuity correction
#> 
#> data:  delta_hours by direction
#> W = 1828, p-value = 0.7181
#> alternative hypothesis: true location shift is not equal to 0

# 3. permutation test (shuffle species labels, keeping real timestamps
#    and each species' total count fixed, recompute the AB/BA ratio
#    each time)
set.seed(1)
n_perm <- 999
rt_pair <- recordTable_soReta[recordTable_soReta$Species %in% c("wolf", "wild boar"), ]
obs_ratio <- median(pair_data$delta_hours[pair_data$direction == "AB" & !pair_data$censored]) /
  median(pair_data$delta_hours[pair_data$direction == "BA" & !pair_data$censored])
null_ratio <- replicate(n_perm, {
  rt_perm <- rt_pair
  rt_perm$Species <- sample(rt_perm$Species)
  out_perm <- build_species_pair_intervals(rt_perm, camOp_soReta, speciesA = "wolf", speciesB = "wild boar",
                                            threshold_min = 30)
  median(out_perm$delta_hours[out_perm$direction == "AB" & !out_perm$censored]) /
    median(out_perm$delta_hours[out_perm$direction == "BA" & !out_perm$censored])
})
mean(null_ratio >= obs_ratio, na.rm = TRUE)  # empirical p-value
#> [1] 0.6246246

A methodological caution before interpreting any of these results: with two species detected at very different rates (a common situation), a strong AB/BA asymmetry can appear even with no real behavioural avoidance -- it's largely an artefact of how often each species is detected at all (Dymit & Levi 2025). Before drawing conclusions, compare the two species' overall detection rates (build_site_total() is the quickest way); the permutation test above already accounts for this by holding each species' total count fixed under the null.

build_species_pair_interruptions() covers the complementary AA/BB/ABA/BAB intervals from the same framework (also supports independence_method/require_uninterrupted). Unlike AB/BA, Niedballa et al. (2019) only applied the permutation test to these four interval types -- not the linear model or Mann-Whitney U-test used above, which in their paper were reserved for AB/BA specifically:

interruptions <- build_species_pair_interruptions(
  recordTable_soReta,
  speciesA = "wolf", speciesB = "wild boar",
  threshold_min = 30
)
table(interruptions$type)
#> 
#>  AA ABA BAB  BB 
#>  46  20  21  43

# permutation test on AA vs BB (same logic as above, applied to the
# interruption-type intervals instead of AB/BA)
obs_ratio_aabb <- median(interruptions$delta_hours[interruptions$type == "AA"]) /
  median(interruptions$delta_hours[interruptions$type == "BB"])
null_ratio_aabb <- replicate(n_perm, {
  rt_perm <- rt_pair
  rt_perm$Species <- sample(rt_perm$Species)
  out_perm <- build_species_pair_interruptions(rt_perm, speciesA = "wolf", speciesB = "wild boar",
                                                threshold_min = 30)
  median(out_perm$delta_hours[out_perm$type == "AA"]) /
    median(out_perm$delta_hours[out_perm$type == "BB"])
})
mean(null_ratio_aabb >= obs_ratio_aabb, na.rm = TRUE)
#> [1] 0.4604605

6. Hierarchical diel activity models

build_diel_binomial_block(), build_diel_binomial_month(), and build_diel_binomial_period() turn detection times into a binomial success/failure dataset -- how many active days a species was detected in each time-bin (bin_hours sets the width of each bin, e.g. 1 for hourly bins), at whatever temporal grain you choose -- ready for the trigonometric GLMMs and cyclic cubic spline HGAMs described in Iannarilli et al. (2024)'s tutorial for hierarchical diel activity models. There is no independence threshold here by design (see the caution in Section 4): detected is already a per-day, per-bin binary indicator, so more photos in the same bin don't change the result.

diel_month <- build_diel_binomial_month(recordTable_soReta, camOp_soReta, bin_hours = 1, min_days = 10)
diel_wolf  <- diel_month[diel_month$Species == "wolf", ]
diel_wolf$Station <- factor(diel_wolf$Station)

# trigonometric GLMM (Iannarilli et al. 2024, section 3.3)
trig_model <- GLMMadaptive::mixed_model(
  fixed  = cbind(success, failure) ~ cos(2 * pi * Time / 24) + sin(2 * pi * Time / 24) +
                                      cos(2 * pi * Time / 12) + sin(2 * pi * Time / 12),
  random = ~ 1 | Station,
  data   = diel_wolf,
  family = binomial()
)
summary(trig_model)
#> 
#> Call:
#> GLMMadaptive::mixed_model(fixed = cbind(success, failure) ~ cos(2 * 
#>     pi * Time/24) + sin(2 * pi * Time/24) + cos(2 * pi * Time/12) + 
#>     sin(2 * pi * Time/12), random = ~1 | Station, data = diel_wolf, 
#>     family = binomial())
#> 
#> Data Descriptives:
#> Number of Observations: 408
#> Number of Groups: 5 
#> 
#> Model:
#>  family: binomial
#>  link: logit 
#> 
#> Fit statistics:
#>    log.Lik      AIC      BIC
#>  -183.0795 378.1591 375.8157
#> 
#> Random effects covariance matrix:
#>                StdDev
#> (Intercept) 0.2293207
#> 
#> Fixed effects:
#>                       Estimate Std.Err  z-value p-value
#> (Intercept)            -5.3926  0.2056 -26.2310 < 1e-04
#> cos(2 * pi * Time/24)   1.3531  0.2726   4.9640 < 1e-04
#> sin(2 * pi * Time/24)  -0.2029  0.2071  -0.9797 0.32722
#> cos(2 * pi * Time/12)  -0.2080  0.2029  -1.0249 0.30539
#> sin(2 * pi * Time/12)   0.1433  0.1982   0.7231 0.46965
#> 
#> Integration:
#> method: adaptive Gauss-Hermite quadrature rule
#> quadrature points: 11
#> 
#> Optimization:
#> method: hybrid EM and quasi-Newton
#> converged: TRUE

7. Classic capture-mark-recapture

For individually identifiable species, build_cmr_day() and build_cmr_block() build an individual x occasion capture history from a recordTable with an individual-ID column (the same shape as camtrapR's own recordTableIndividual() output). Station plays no role in the result. as_capture_strings() converts that into the capture-history string format used by RMark, MARK, and marked. Neither of these two functions has an independence threshold: distinct individual sightings on the same day are already unambiguous.

recordTableIndividuals_soReta provides this directly: 8 individually-recognizable red deer ("RD_01" to "RD_08"), built by tagging the "red deer" detections already present in recordTable_soReta -- same stations, same camOp_soReta used throughout this vignette, no separate dataset or dependency needed.

ch <- build_cmr_block(recordTableIndividuals_soReta, camOp_soReta, block_days = 7, min_days = 4)
as_capture_strings(ch)
#> $`red deer`
#>                RD_01                RD_02                RD_03 
#> "01111010110010110." "11111101100011110." "10101001101101110." 
#>                RD_04                RD_05                RD_06 
#> "10110011111000010." "00010110000010110." "00010000111100110." 
#>                RD_07                RD_08 
#> "01000010001000010." "01110000000010010."

Where to go from here

Every function's help page (?build_site_month, ?build_species_pair_intervals, ...) documents its parameters and return value in full. This vignette only shows the default behaviour -- most functions have a threshold_min, min_days, or similar argument worth tuning to your own study design.

References