Exercises week 5

Load packages

library(readr)
library(tidyverse)

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.4     ✔ purrr     1.1.0
✔ forcats   1.0.1     ✔ stringr   1.5.2
✔ ggplot2   4.0.0     ✔ tibble    3.3.0
✔ lubridate 1.9.4     ✔ tidyr     1.3.1
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(lubridate)
library(sf)

Linking to GEOS 3.13.0, GDAL 3.8.5, PROJ 9.5.1; sf_use_s2() is TRUE

library(dplyr)
library(ggplot2)

Demo

testfun <- function(sometext) {print(sometext)} 

testfun("this function does slighlty more but still not much ")

[1] "this function does slighlty more but still not much "

my_age <- function(birthday, output_unit ="auto") {
  today <- Sys.Date() 
  difftime(today,birthday, units = output_unit)
}

my_age("1996-02-11", "weeks")

Time difference of 1578.577 weeks

Task 1 : Write your own functions

# BMI function
bmi_func <- function(weight_kg, height_m) {
  weight_kg / height_m^2
}

bmi_func(weight_kg = 55, height_m = 1.63)

[1] 20.70082

# Celsius to Fahrenheit function
celsius_to_fahrenheit <- function(celsius) {
  celsius * 9/5 + 32
}

celsius_to_fahrenheit(celsius = 20)

[1] 68

# Euclidean distance function
euclidean_distance <- function(x1, y1, x2, y2) {
  sqrt((x2 - x1)^2 + (y2 - y1)^2)
}

euclidean_distance(x1 = 0, y1 = 0, x2 = 3, y2 = 4)

[1] 5

Task 2 : Prepare Analysis

#import data
wildschwein <- read_delim("wildschwein_BE_2056.csv")

Rows: 51246 Columns: 6
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr  (2): TierID, TierName
dbl  (3): CollarID, E, N
dttm (1): DatetimeUTC

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

#filter 
wildschwein_filter <- wildschwein |>
  st_as_sf(coords = c("E", "N"), crs = 2056, remove = FALSE) |> 
  mutate(DatetimeUTC = ymd_hms(DatetimeUTC)) |>
  filter(
    TierName %in% c("Rosa", "Sabi"),
    DatetimeUTC >= ymd_hms("2015-04-01 00:00:00"),
    DatetimeUTC <= ymd_hms("2015-04-15 23:59:59")
  )

#check
wildschwein_filter |> 
  count(TierName)

Simple feature collection with 2 features and 2 fields
Geometry type: MULTIPOINT
Dimension:     XY
Bounding box:  xmin: 2569715 ymin: 1202620 xmax: 2574355 ymax: 1205669
Projected CRS: CH1903+ / LV95
# A tibble: 2 × 3
  TierName     n                                                        geometry
* <chr>    <int>                                                <MULTIPOINT [m]>
1 Rosa      1420 ((2570177 1204497), (2570185 1204495), (2570186 1204493), (257…
2 Sabi      1440 ((2569715 1205295), (2569732 1205288), (2569734 1205262), (256…

range(wildschwein_filter$DatetimeUTC)

[1] "2015-04-01 00:00:10 UTC" "2015-04-15 23:47:56 UTC"

Task 3: Create Join Key

# round the tiemstamps to the nearest multiple of 15 
wildschwein_filter <- wildschwein_filter %>%
  mutate(DatetimeRound = round_date(DatetimeUTC, unit = "15 minutes"))

wildschwein_filter %>%
  select(TierID, TierName, CollarID, DatetimeUTC, E, N, DatetimeRound) %>%
  head()

Simple feature collection with 6 features and 7 fields
Geometry type: POINT
Dimension:     XY
Bounding box:  xmin: 2570309 ymin: 1205237 xmax: 2570372 ymax: 1205313
Projected CRS: CH1903+ / LV95
# A tibble: 6 × 8
  TierID TierName CollarID DatetimeUTC              E      N DatetimeRound      
  <chr>  <chr>       <dbl> <dttm>               <dbl>  <dbl> <dttm>             
1 002A   Sabi        12275 2015-04-01 00:00:11 2.57e6 1.21e6 2015-04-01 00:00:00
2 002A   Sabi        12275 2015-04-01 00:15:22 2.57e6 1.21e6 2015-04-01 00:15:00
3 002A   Sabi        12275 2015-04-01 00:30:11 2.57e6 1.21e6 2015-04-01 00:30:00
4 002A   Sabi        12275 2015-04-01 00:45:16 2.57e6 1.21e6 2015-04-01 00:45:00
5 002A   Sabi        12275 2015-04-01 01:00:44 2.57e6 1.21e6 2015-04-01 01:00:00
6 002A   Sabi        12275 2015-04-01 01:15:17 2.57e6 1.21e6 2015-04-01 01:15:00
# ℹ 1 more variable: geometry <POINT [m]>

Task 4: Measuring distance at concurrent locations

#split inro one data.frame per animal 
sabi <- wildschwein_filter |> 
  filter(TierName == "Sabi")

rosa <- wildschwein_filter |> 
  filter(TierName == "Rosa")


#join by key created
boar_joined <- inner_join(
  st_drop_geometry(sabi),
  st_drop_geometry(rosa),
  by = "DatetimeRound",
  suffix = c("_sabi", "_rosa")
)


# calculate the distance between the 2 objects 

boar_joined <- boar_joined |> 
  mutate(
    distance_m = euclidean_distance(
      x1 = E_sabi,
      y1 = N_sabi,
      x2 = E_rosa,
      y2 = N_rosa
    )
  )

#Define meetings
#A meeting happens when distance is less than or equal to 100 meters.
boar_joined <- boar_joined |> 
  mutate(
    meet = distance_m <= 100
  )

# how many meetings there were ? 
boar_joined |> 
  count(meet)

# A tibble: 2 × 2
  meet      n
  <lgl> <int>
1 FALSE  1409
2 TRUE     11

#view meeting rows

meets <- boar_joined |> 
  filter(meet == TRUE)

meets

# A tibble: 11 × 15
   TierID_sabi TierName_sabi CollarID_sabi DatetimeUTC_sabi      E_sabi   N_sabi
   <chr>       <chr>                 <dbl> <dttm>                 <dbl>    <dbl>
 1 002A        Sabi                  12275 2015-04-02 04:45:16 2570422. 1205393.
 2 002A        Sabi                  12275 2015-04-02 05:00:11 2570387. 1205381.
 3 002A        Sabi                  12275 2015-04-02 18:00:15 2570368. 1205370.
 4 002A        Sabi                  12275 2015-04-02 18:30:23 2570409. 1205318.
 5 002A        Sabi                  12275 2015-04-03 04:00:13 2570298. 1205391.
 6 002A        Sabi                  12275 2015-04-03 04:15:09 2570302. 1205386.
 7 002A        Sabi                  12275 2015-04-03 04:45:10 2570325. 1205377.
 8 002A        Sabi                  12275 2015-04-03 05:45:12 2570302. 1205377.
 9 002A        Sabi                  12275 2015-04-03 06:00:43 2570294. 1205386.
10 002A        Sabi                  12275 2015-04-03 06:15:18 2570298. 1205388.
11 002A        Sabi                  12275 2015-04-03 06:30:14 2570304. 1205388.
# ℹ 9 more variables: DatetimeRound <dttm>, TierID_rosa <chr>,
#   TierName_rosa <chr>, CollarID_rosa <dbl>, DatetimeUTC_rosa <dttm>,
#   E_rosa <dbl>, N_rosa <dbl>, distance_m <dbl>, meet <lgl>

### Visualize Data

ggplot() +
  # Regular locations: Rosa
  geom_point(
    data = rosa,
    aes(x = E, y = N, color = "rosa"),
    alpha = 0.25,
    size = 1.5
  ) +
  # Regular locations: Sabi
  geom_point(
    data = sabi,
    aes(x = E, y = N, color = "sabi"),
    alpha = 0.25,
    size = 1.5
  ) +
  # Meets: Rosa
  geom_point(
    data = meets,
    aes(x = E_rosa, y = N_rosa, shape = "rosa"),
    color = "red",
    size = 2
  ) +
  # Meets: Sabi
  geom_point(
    data = meets,
    aes(x = E_sabi, y = N_sabi, shape = "sabi"),
    color = "darkcyan",
    size = 2
  ) +
  labs(
    title = "Meeting locations of Rosa and Sabi",
    x = "E",
    y = "N",
    color = "Regular Locations",
    shape = "Meets"
  ) +
  coord_equal() +
  theme_minimal()

### Task 6 (optional): Visualize data as timecube with plotly

library(plotly)

Attaching package: 'plotly'

The following object is masked from 'package:ggplot2':

    last_plot

The following object is masked from 'package:stats':

    filter

The following object is masked from 'package:graphics':

    layout

plot_ly() |> 
  # Sabi trajectory
  add_trace(
    data = sabi,
    x = ~E,
    y = ~N,
    z = ~DatetimeUTC,
    type = "scatter3d",
    mode = "lines",
    name = "Sabi"
  ) |> 
  # Rosa trajectory
  add_trace(
    data = rosa,
    x = ~E,
    y = ~N,
    z = ~DatetimeUTC,
    type = "scatter3d",
    mode = "lines",
    name = "Rosa"
  ) |> 
  # meet locations based on Sabi coordinates
  add_trace(
    data = meets,
    x = ~E_sabi,
    y = ~N_sabi,
    z = ~DatetimeRound,
    type = "scatter3d",
    mode = "markers",
    name = "Meets",
    marker = list(size = 5)
  ) |> 
  layout(
    title = "Space-time cube of meeting patterns between Sabi and Rosa",
    scene = list(
      xaxis = list(title = "E"),
      yaxis = list(title = "N"),
      zaxis = list(title = "Time")
    )
  )