4 Light map
In the advanced tutorial, we construct the same three products with a Great Reed Warbler (18LX), this time including light and wind data.
This first chapter focuses on constructing a likelihood map from light data. This optional step is helpful in most cases, particularly for short stopovers bracketed by long flights. It can also be used to cross-check the pressure map and reduce the computational cost of building the graph.
Here, we use a threshold-based approach calibrated in terms of solar zenith angle rather than sunrise and sunset times. This approach is presented in more detail in section 2.4 of Nussbaumer et al. (2023).
A more thorough introduction to light-level geolocation is available in the geolocation manual (Lisovski et al. 2020). Other methods can also produce likelihood maps, including those presented in Basson et al. (2016) and Bindoff et al. (2018).
4.1 Basic tutorial catch up
Before working with light data, we need to create and label the tag object and define its map. We repeat the steps from Tag object and Pressure map in a few lines using the native pipe operator |>. This assumes that labelling has already been completed.
tag <- tag_create(
"18LX",
crop_start = "2017-06-20",
crop_end = "2018-05-02",
quiet = TRUE
) |>
tag_label(quiet = TRUE) |>
tag_set_map(
extent = c(-16, 23, 0, 50),
scale = 4,
known = data.frame(
stap_id = 1,
known_lat = 48.9,
known_lon = 17.05
)
) |>
geopressure_map(quiet = TRUE)4.2 Estimate twilights
We identify daily sunrise and sunset events with twilight_create(). It performs the same task as TwGeos::FindTwilight(), but uses a matrix representation that is faster, though potentially less general.
By default, the light threshold twl_thr is defined from the first and last non-zero light measurements each day (i.e., tag$light$value > 0). The twl_offset parameter centres the night and day in the matrix representation. Good centring is necessary to identify the correct twilights.
tag <- twilight_create(tag, twl_thr = NULL, twl_offset = NULL)We can visualise the twilights and check that the day is centred correctly.
plot(tag, type = "twilight")Re-run twilight_create() with a different twl_offset until the night and day are properly centred.
Because we have already labelled stationary periods, the colour of each point indicates its period. This helps distinguish outliers from a change in the bird’s position.
4.3 Manual labelling of twilight
Use the GeoLightViz Shiny app included in GeoPressureR to discard twilight outliers. See the GeoLightViz chapter for an overview of its controls and workflows.
geolightviz(tag)Twilight outliers can often be identified visually because they do not follow a smooth line.
The point colours help distinguish outliers from a change in the bird’s position. Twilight times within a stationary period should change smoothly, whereas changes between positions can be abrupt.
Changing a twilight’s stationary-period label has no effect later in the workflow, because twilight_label_read() uses only the "discard" label.
Avoid over-editing calibration periods. Twilight variability is needed to build a calibration that represents the uncertainty of individual twilights. Outliers are easier to identify during long stationary periods than during short ones, which can leave more variability in short periods and bias their estimated positions.
When you have finished labelling, read and check the twilight annotations.
tag <- twilight_label_read(tag)
plot(tag, type = "twilight")4.4 Compute likelihood map
The light likelihood maps are built in three steps:
- Calibrate the solar zenith angles compatible with the recorded twilights with
geolight_map_calibrate(). - Evaluate that calibration at every map location for each twilight with
geolight_map_likelihood(). - Aggregate twilight maps by stationary period with
geolight_map_aggregate().
These three steps are generally performed with geolight_map().
4.4.1 Calibrate twilight error
A recorded sunrise or sunset does not always occur at exactly the same solar zenith angle. Local shading, weather, habitat, and tag behaviour can all affect when the light threshold is crossed. We therefore calibrate the twilight observations in terms of solar zenith angle rather than a timing error. This makes the calibration more transferable across locations and seasons.
Here, we use hybrid calibration. The calibration density is built from twilights recorded at the known breeding location, supplied when we created the tag object, and from long stationary periods whose locations are fitted from the twilight data. Combining these sources helps the calibration represent light conditions across the annual cycle rather than relying on a single site.
tag <- geolight_map_calibrate(
tag,
twl_calib_adjust = 1.4,
fitted_location_duration = 30,
refine_fitted_location_max_iter = 2
)The function converts retained twilights at these calibration locations into zenith angles and fits a smooth density to their distribution using stats::density(). twl_calib_adjust controls the smoothness of this density. Since the calibration is applied throughout the annual cycle, a slightly broader fit is generally safer (twl_calib_adjust > 1).
For each eligible stationary period, GeoPressureR assumes that the bird remained at one location. It then finds the longitude, latitude, and solar zenith angle that best reconcile the retained twilight observations from that period. A weak prior on the zenith angle helps stabilise fits when the twilight geometry alone is not sufficiently informative.
The fitted location is then refined by constructing a temporary light likelihood map from the same twilights and moving the location towards its local mode. Here, refine_fitted_location_max_iter = 2 performs two refinement iterations; known locations remain fixed.
fitted_location_duration = 30 makes periods of at least 30 days eligible for fitting. These fitted locations are calibration anchors, not final, uncertainty-aware location estimates. They can be used on their own when no known calibration location is available, or together with known locations as in this example.
You should always check calibration before computing the map:
plot_twl_calib(tag)You can also compare the calibration against a provisional path:
plot_twl_calib(tag, path = tag2path(tag, interp = 2))Warning: ! First and/or last modelled stationary periods (1 and 29) are shorter than 2
days but cannot be interpolated.
→ They will not be interpolated.For each stationary period, this diagnostic compares the distribution of observed twilight zenith angles along the provisional path with the calibrated distribution. Periods highlighted as out of range may indicate a problem with the twilight labels, calibration settings, or provisional path. Use this plot iteratively to refine the labels and calibration before modelling the final trajectory.
4.4.2 Compute likelihood maps for individual twilights
We can now compute a likelihood map for each twilight event. For each map cell, GeoPressureR calculates the solar zenith angle expected at the observed twilight time and evaluates its compatibility with the calibrated zenith-angle distribution.
tag <- geolight_map_likelihood(
tag,
compute_known = FALSE
)These maps are stored in tag$map_light_twl and are the inputs for aggregation by stationary period. With compute_known = FALSE, GeoPressureR does not calculate twilight maps for periods with known coordinates, reducing computation. Those periods are retained later as fixed locations.
4.4.3 Aggregate twilight maps by stationary period
The trajectory model works at the stationary-period level, rather than treating every twilight as an independent location estimate. We therefore combine all retained twilight maps within each stationary period.
tag <- geolight_map_aggregate(
tag,
twl_llp = \(n) log(n) / n,
compute_known = FALSE
)The twl_llp argument sets the log-linear pooling weight. The default, \(n) log(n) / n, down-weights each twilight as their number (n) increases. This helps limit overconfidence when consecutive twilights share errors caused by persistent shading, weather, habitat, or behaviour.
The output of this step is tag$map_light, which is the map used in the rest of the trajectory workflow.
Finally, we can visualise the aggregated likelihood map for each stationary period:
plot(tag, type = "map_light")Before building the graph, compare the light and pressure likelihood maps. They should show meaningful spatial agreement. Strong disagreement is a reason to revisit the tag or twilight labels and the calibration settings.
This task is best performed with GeoPressureViz, presented in its dedicated chapter. For now, we can simply visualise the combined pressure and light likelihood maps.
plot(tag, type = "map")4.5 Check light label
As with pressure labels, light labels should also be checked. We first estimate a trajectory, calculate the twilights expected along it, and compare them with the observed twilights. Here, tag2path() finds the most likely position for each stationary period, regardless of the flight-duration or movement model. We interpolate periods shorter than two days to avoid unrealistic position estimates.
path <- tag2path(tag, interp = 2)This path can be visualised with plot_path(path).
We can calculate the theoretical twilights expected for a bird following this path with path2twilight().
twilight_line <- path2twilight(path)We can compare these theoretical twilights with the observed twilights using plot_tag_twilight().
plot_tag_twilight(tag, twilight_line = twilight_line, plot_plotly = T)