4_1

Using sows with greater than 60 hours of data before the date of farrowing was recorded.

Loading in Libraries

library(dplyr)

Attaching package: 'dplyr'

The following objects are masked from 'package:stats':

    filter, lag

The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union

library(tidyr)
library(forcats)
library(ggplot2)
library(purrr)
library(slider)

Warning: package 'slider' was built under R version 4.4.3

library(lubridate)

Attaching package: 'lubridate'

The following objects are masked from 'package:base':

    date, intersect, setdiff, union

library(ggforce)

Warning: package 'ggforce' was built under R version 4.4.3

library(FKF)

Warning: package 'FKF' was built under R version 4.4.3

library(data.table)

Attaching package: 'data.table'

The following objects are masked from 'package:lubridate':

    hour, isoweek, mday, minute, month, quarter, second, wday, week,
    yday, year

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

    transpose

The following objects are masked from 'package:dplyr':

    between, first, last

Loading in Data

setwd("C:/Users/ipberg/OneDrive - Iowa State University/Spring26/Sensors/")
load("4_1.RData")

# setwd("/Users/IsaacBerg/Documents/Code/Sensors/")
# load("4_1.RData")

Not keeping the first 12 hours of data collected and using the next 24 as training

training_data <- focused %>% 
  group_by(Sow_ID) %>% 
  # Filter out the first 12 hours of data per sow
  # filter(sensor_ts >= min(sensor_ts) + hours(12)) %>%
  mutate(
    hours_before = as.numeric(difftime(FD_Gestal, min(sensor_ts), units="hours")),
    reaches_farrow = max(sensor_ts) >= as.POSIXct(FD_Gestal)
  ) %>% 
  filter(hours_before >= 60 & reaches_farrow == TRUE) %>% 
  filter(sensor_ts < as.POSIXct((FD_Gestal + days(2)))) %>% 
  mutate(
    time_to_farrow = -as.numeric(difftime(as.POSIXct(FD_Gestal), sensor_ts, units = "hours")),
    m_mean_24h = slide_index_dbl(m, sensor_ts, mean, .before = hours(24)),
    m_var_24h  = slide_index_dbl(m, sensor_ts, var, .before = hours(24))
  ) %>%
  filter(as.numeric(format(sensor_ts, "%M")) %% 15 == 0 & 
           as.numeric(format(sensor_ts, "%S")) == 0) %>%
  mutate(start_time = min(sensor_ts)) %>%
  # Filter for the window: between 12 and 36 hours from the start
  filter(sensor_ts >= (start_time + hours(12)) & 
         sensor_ts <  (start_time + hours(36))) %>%
  select(-hours_before, -reaches_farrow) %>% 
  nest()

Subset of data with only sows with 60 hours of data or greater before farrowing - removing the first 12 hours

sows_greater_60 <- focused %>% 
  group_by(Sow_ID) %>% 
  mutate(
    hours_before = as.numeric(difftime(FD_Gestal, min(sensor_ts), units="hours")),
    reaches_farrow = max(sensor_ts) >= as.POSIXct(FD_Gestal)
  ) %>% 
  filter(hours_before >= 60 & reaches_farrow == TRUE) %>% 
  filter(sensor_ts < as.POSIXct((FD_Gestal + days(2)))) %>% 
  mutate(
    time_to_farrow = -as.numeric(difftime(as.POSIXct(FD_Gestal), sensor_ts, units = "hours")),
    m_mean_24h = slide_index_dbl(m, sensor_ts, mean, .before = hours(24)),
    m_var_24h  = slide_index_dbl(m, sensor_ts, var,  .before = hours(24))
  ) %>%
  filter(as.numeric(format(sensor_ts, "%M")) %% 15 == 0 & 
         as.numeric(format(sensor_ts, "%S")) == 0) %>%
  mutate(start_time = min(sensor_ts)) %>%
  select(-hours_before, -reaches_farrow) %>%
  semi_join(training_data %>% unnest(data), by = "Sow_ID") %>%  # <-- only keeps sows that exist in training
  nest()

Plotting the variance of each sow on the same plot relative to their date of farrowing

# 1. Prepare data for plotting
plot_data <- sows_greater_60 %>%
  unnest(data)

# 2. Create the plot
ggplot(plot_data, aes(x = time_to_farrow, y = m_var_24h)) +
  # Add a line for each individual sow (lightly colored)
  geom_line(aes(group = Sow_ID), alpha = 0.1, color = "steelblue") +
  # Add a smoothed average line to see the general trend
  geom_smooth(color = "red", method = "gam") + 
  # Mark the farrowing event
  geom_vline(xintercept = 0, linetype = "dashed", color = "black") +
  annotate("text", x = 0.2, y = max(plot_data$m_var_24h, na.rm=T), 
           label = "Farrowing Date", angle = 90, vjust = -0.5) +
  theme_minimal() +
  labs(
    title = "Sow Activity Variance Relative to Farrowing",
    x = "Days Relative to Farrowing (0 = Start of FD_Gestal)",
    y = "24h Rolling Variance (m)"
  )

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

Warning: Removed 40 rows containing non-finite outside the scale range
(`stat_smooth()`).

Warning: Removed 40 rows containing missing values or values outside the scale range
(`geom_line()`).

using training data

# 1. Prepare data for plotting
plot_data_12 <- training_data %>%
  unnest(data) 

# 2. Create the plot
ggplot(plot_data_12, aes(x = time_to_farrow, y = m_var_24h)) +
  # Add a line for each individual sow (lightly colored)
  geom_line(aes(group = Sow_ID), alpha = 0.1, color = "steelblue") +
  # Add a smoothed average line to see the general trend
  geom_smooth(color = "red", method = "gam") + 
  # Mark the farrowing event
  geom_vline(xintercept = 0, linetype = "dashed", color = "black") +
  annotate("text", x = 0.2, y = max(plot_data_12$m_var_24h, na.rm=T), 
           label = "Farrowing Date", angle = 90, vjust = -0.5) +
  theme_minimal() +
  labs(
    title = "Sow Activity Variance Relative to Farrowing",
    x = "Days Relative to Farrowing (0 = Start of FD_Gestal)",
    y = "24h Rolling Variance (m)"
  )

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

Plotting a single sow

library(ggplot2)
library(dplyr)

# 1. Pick one sow and unnest her data
single_sow_data <- sows_greater_60 %>%
  filter(Sow_ID == "26240") %>% unnest(data) %>% 
  filter(sensor_ts >= (start_time + hours(12)))
  

# 2. Plot
ggplot(single_sow_data, aes(x = time_to_farrow, y = m_var_24h)) +
  geom_line(color = "firebrick", size = 1) +
  geom_point(size = 0.5, alpha = 0.5) + # Shows the 15-min steps clearly
  geom_vline(xintercept = 0, linetype = "dashed") +
  geom_vline(xintercept = single_sow_data$time_to_farrow[1] -12)+
  # geom_smooth(se=FALSE,span=0.2)+
  theme_minimal() +
  labs(
    title = paste("Activity Variance: Sow", unique(single_sow_data$Sow_ID)),
    subtitle = "24-hour rolling window sampled every 15 minutes",
    x = "Hours Relative to Farrowing (0 = Start of FD_Gestal)",
    y = "Variance of Magnitude (m)"
  )

Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
ℹ Please use `linewidth` instead.

plotting the training data for a single sow

library(ggplot2)
library(dplyr)

# 1. Pick one sow and unnest her data
single_sow_data_12 <- training_data %>%
  filter(Sow_ID == "26240") %>% 
  unnest(data) 

# 2. Plot
ggplot(single_sow_data_12, aes(x = time_to_farrow, y = m_var_24h)) +
  geom_line(color = "firebrick", size = 1) +
  geom_point(size = 0.5, alpha = 0.5) + # Shows the 15-min steps clearly
  geom_vline(xintercept = 0, linetype = "dashed") +
  geom_vline(xintercept = single_sow_data$time_to_farrow[1] -12)+
  # geom_smooth(se=FALSE,span=0.2)+
  theme_minimal() +
  labs(
    title = paste("Activity Variance: Sow", unique(single_sow_data$Sow_ID)),
    subtitle = "24-hour rolling window sampled every 15 minutes",
    x = "Hours Relative to Farrowing (0 = Start of FD_Gestal)",
    y = "Variance of Magnitude (m)"
  )

24 hour 15 minute step plots of sows with greater than 60hrs of data before farrowing date

Plots of all sows

Plots are of the magnitude variance

# 1. Prepare all data
all_sows_unnested <- sows_greater_60 %>%
  unnest(data)

# 2. Get unique IDs
sow_list <- unique(all_sows_unnested$Sow_ID)

# 3. Loop and Plot
for (sow in sow_list) {
  
  plot_data <- all_sows_unnested %>% filter(Sow_ID == sow) %>% 
      filter(sensor_ts >= (start_time + hours(12)))
  
  p <- ggplot(plot_data, aes(x = time_to_farrow, y = m_var_24h)) +
    geom_line(color = "firebrick", size = 0.8) +
    geom_vline(xintercept = 0, linetype = "dashed", color = "darkgrey") +
    geom_vline(xintercept = plot_data$time_to_farrow[1] -12)+
    geom_vline(xintercept = plot_data$time_to_farrow[1]+24)
    theme_minimal() +
    labs(
      title = paste("Activity Variance: Sow", sow),
      subtitle = "24-hour rolling window (15-min steps)",
      x = "Hours Relative to Farrowing",
      y = "Variance of Magnitude"
    )
  
  print(p) 
  # Optional: ggsave(paste0("Sow_", sow, ".png"), p)
}

for removing the first 12 hours

# 1. Prepare all data
all_sows_unnested_12 <- training_data %>%
  unnest(data)

# 2. Get unique IDs
sow_list <- unique(all_sows_unnested_12$Sow_ID)

# 3. Loop and Plot
for (sow in sow_list) {
  
  plot_data <- all_sows_unnested_12 %>% filter(Sow_ID == sow)%>% 
      filter(sensor_ts >= (start_time + hours(12)))
  
  p <- ggplot(plot_data, aes(x = time_to_farrow, y = m_var_24h)) +
    geom_line(color = "firebrick", size = 0.8) +
    geom_vline(xintercept = 0, linetype = "dashed", color = "darkgrey") +
    geom_vline(xintercept = plot_data$time_to_farrow[1] -12)+
    theme_minimal() +
    labs(
      title = paste("Activity Variance: Sow", sow),
      subtitle = "24-hour rolling window (15-min steps)",
      x = "Hours Relative to Farrowing",
      y = "Variance of Magnitude"
    )
  
  print(p) 
  # Optional: ggsave(paste0("Sow_", sow, ".png"), p)
}

Applying the Kalmon filtering

Using the first 48hrs of data as “baseline” for model to compare with

y <- training_data %>% 
  unnest(cols=everything()) %>% 
  filter(Sow_ID==26240) %>% 
  select(m_var_24h)

Adding missing grouping variables: `Sow_ID`

y <- y[[2]]

a0 <- y[1]
P0 <- matrix(1)
dt <- ct <- matrix(0)
Zt <- Tt <- matrix(1)

fit.fkf <- optim(c(HHt = var(y, na.rm = TRUE) * .5,
                   GGt = var(y, na.rm = TRUE) * .5),
                 fn = function(par, ...)
                 -fkf(HHt = matrix(par[1]), GGt = matrix(par[2]), ...)$logLik,
                 yt = rbind(y), a0 = a0, P0 = P0, dt = dt, ct = ct,
                 Zt = Zt, Tt = Tt)

HHt <- as.numeric(fit.fkf$par[1])
GGt <- as.numeric(fit.fkf$par[2])
HHt; GGt

[1] 3.998639e-09

[1] -7.965795e-10

y_fkf <- fkf(a0, P0, dt, ct, Tt, Zt,
             HHt = matrix(HHt), GGt = matrix(GGt),
             yt = rbind(y))


data <- data.table(x=single_sow_data$time_to_farrow,
                  y=as.numeric(single_sow_data$m_var_24h),
                   y_kalman = as.numeric(y_fkf$att))

Warning in as.data.table.list(x, keep.rownames = keep.rownames, check.names =
check.names, : Item 3 has 89 rows but longest item has 344; recycled with
remainder.

ggplot(data, aes(x = x, y = y)) +
  geom_line() +
  geom_line(data = data, aes(x = x, y = y_kalman), col = "blue") +
  expand_limits(y=0)+
  geom_vline(xintercept = 0)+
  geom_vline(xintercept = data$x[1]-12)+
  geom_vline(xintercept = data$x[1]+24)+
  xlab("time relative to farrowing")  +
  ggtitle("24h window w/ 15 min steps with first 12hrs removed and Kalman Filtering") +
  theme_bw()

y_full <- as.numeric(single_sow_data$m_var_24h)


# Manual Test for Smoothing
y_fkf_smooth <- fkf(
  a0 = a0, P0 = P0, dt = dt, ct = ct, Tt = Tt, Zt = Zt,
  HHt = matrix(0.000001), # System Noise (Lower = Smoother)
  GGt = matrix(0.001),   # Measurement Noise (Higher = Smoother)
  yt = rbind(y_full)
)

# Plot this version
plot_smooth <- data.table(
  x = single_sow_data$time_to_farrow,
  y = y_full,
  y_kalman = as.numeric(y_fkf_smooth$att[1, ])
)

ggplot(plot_smooth, aes(x = x)) +
  geom_line(aes(y = y), color = "black", alpha = 0.3) + # Raw data in background
  geom_line(aes(y = y_kalman), color = "blue", size = 1) +
  expand_limits(y=0)+
  geom_vline(xintercept = 0)+
  geom_vline(xintercept = data$x[1]-12)+
  geom_vline(xintercept = data$x[1]+24)+
  xlab("time relative to farrowing")  +
  ggtitle("24h window w/ 15 min steps with first 12hrs removed and Kalman Filtering Sow 26240") +
  theme_bw()

Additional Sow

single_sow_data_2 <- sows_greater_60 %>%
  filter(Sow_ID == "28738") %>% unnest(data) %>% 
  filter(sensor_ts >= (start_time + hours(12)))


y <- training_data %>% 
  unnest(cols=everything()) %>% 
  filter(Sow_ID==28738) %>% 
  select(m_var_24h)

Adding missing grouping variables: `Sow_ID`

y <- y[[2]]

a0 <- y[1]
P0 <- matrix(1)
dt <- ct <- matrix(0)
Zt <- Tt <- matrix(1)

fit.fkf <- optim(c(HHt = var(y, na.rm = TRUE) * .5,
                   GGt = var(y, na.rm = TRUE) * .5),
                 fn = function(par, ...)
                 -fkf(HHt = matrix(par[1]), GGt = matrix(par[2]), ...)$logLik,
                 yt = rbind(y), a0 = a0, P0 = P0, dt = dt, ct = ct,
                 Zt = Zt, Tt = Tt)

HHt <- as.numeric(fit.fkf$par[1])
GGt <- as.numeric(fit.fkf$par[2])
HHt; GGt

[1] 4.278644e-09

[1] -5.018446e-10

y_fkf <- fkf(a0, P0, dt, ct, Tt, Zt,
             HHt = matrix(HHt), GGt = matrix(GGt),
             yt = rbind(y))


data <- data.table(x=single_sow_data_2$time_to_farrow,
                  y=as.numeric(single_sow_data_2$m_var_24h),
                   y_kalman = as.numeric(y_fkf$att))

Warning in as.data.table.list(x, keep.rownames = keep.rownames, check.names =
check.names, : Item 3 has 89 rows but longest item has 344; recycled with
remainder.

ggplot(data, aes(x = x, y = y)) +
  geom_line() +
  geom_line(data = data, aes(x = x, y = y_kalman), col = "blue") +
  expand_limits(y=0)+
  geom_vline(xintercept = 0)+
  geom_vline(xintercept = data$x[1]-12)+
  geom_vline(xintercept = data$x[1]+24)+
  xlab("time relative to farrowing")  +
  ggtitle("24h window w/ 15 min steps with first 12hrs removed and Kalman Filtering Sow 28738") +
  theme_bw()

y_full <- as.numeric(single_sow_data_2$m_var_24h)


# Manual Test for Smoothing
y_fkf_smooth <- fkf(
  a0 = a0, P0 = P0, dt = dt, ct = ct, Tt = Tt, Zt = Zt,
  HHt = matrix(0.000001), # System Noise (Lower = Smoother)
  GGt = matrix(0.001),   # Measurement Noise (Higher = Smoother)
  yt = rbind(y_full)
)

# Plot this version
plot_smooth <- data.table(
  x = single_sow_data_2$time_to_farrow,
  y = y_full,
  y_kalman = as.numeric(y_fkf_smooth$att[1, ])
)

ggplot(plot_smooth, aes(x = x)) +
  geom_line(aes(y = y), color = "black", alpha = 0.3) + # Raw data in background
  geom_line(aes(y = y_kalman), color = "blue", size = 1) +
  expand_limits(y=0)+
  geom_vline(xintercept = 0)+
  geom_vline(xintercept = data$x[1]-12)+
  geom_vline(xintercept = data$x[1]+24)+
  xlab("time relative to farrowing")  +
  ggtitle("24h window w/ 15 min steps with first 12hrs removed and Kalman Filtering Sow 28738") +
  theme_bw()

single_sow_data_3 <- sows_greater_60 %>%
  filter(Sow_ID == "34083") %>% unnest(data) %>% 
  filter(sensor_ts >= (start_time + hours(12)))


y <- training_data %>% 
  unnest(cols=everything()) %>% 
  filter(Sow_ID==34083) %>% 
  select(m_var_24h)

Adding missing grouping variables: `Sow_ID`

y <- y[[2]]

a0 <- y[1]
P0 <- matrix(1)
dt <- ct <- matrix(0)
Zt <- Tt <- matrix(1)

fit.fkf <- optim(c(HHt = var(y, na.rm = TRUE) * .5,
                   GGt = var(y, na.rm = TRUE) * .5),
                 fn = function(par, ...)
                 -fkf(HHt = matrix(par[1]), GGt = matrix(par[2]), ...)$logLik,
                 yt = rbind(y), a0 = a0, P0 = P0, dt = dt, ct = ct,
                 Zt = Zt, Tt = Tt)

HHt <- as.numeric(fit.fkf$par[1])
GGt <- as.numeric(fit.fkf$par[2])
HHt; GGt

[1] 6.994863e-09

[1] 2.367943e-11

y_fkf <- fkf(a0, P0, dt, ct, Tt, Zt,
             HHt = matrix(HHt), GGt = matrix(GGt),
             yt = rbind(y))


data <- data.table(x=single_sow_data_3$time_to_farrow,
                  y=as.numeric(single_sow_data_3$m_var_24h),
                   y_kalman = as.numeric(y_fkf$att))

Warning in as.data.table.list(x, keep.rownames = keep.rownames, check.names =
check.names, : Item 3 has 88 rows but longest item has 491; recycled with
remainder.

ggplot(data, aes(x = x, y = y)) +
  geom_line() +
  geom_line(data = data, aes(x = x, y = y_kalman), col = "blue") +
  expand_limits(y=0)+
  geom_vline(xintercept = 0)+
  geom_vline(xintercept = data$x[1]-12)+
  geom_vline(xintercept = data$x[1]+24)+
  xlab("time relative to farrowing")  +
  ggtitle("24h window w/ 15 min steps with first 12hrs removed and Kalman Filtering Sow 34083") +
  theme_bw()

y_full <- as.numeric(single_sow_data_3$m_var_24h)


# Manual Test for Smoothing
y_fkf_smooth <- fkf(
  a0 = a0, P0 = P0, dt = dt, ct = ct, Tt = Tt, Zt = Zt,
  HHt = matrix(0.000001), # System Noise (Lower = Smoother)
  GGt = matrix(0.001),   # Measurement Noise (Higher = Smoother)
  yt = rbind(y_full)
)

# Plot this version
plot_smooth <- data.table(
  x = single_sow_data_3$time_to_farrow,
  y = y_full,
  y_kalman = as.numeric(y_fkf_smooth$att[1, ])
)

ggplot(plot_smooth, aes(x = x)) +
  geom_line(aes(y = y), color = "black", alpha = 0.3) + # Raw data in background
  geom_line(aes(y = y_kalman), color = "blue", size = 1) +
  expand_limits(y=0)+
  geom_vline(xintercept = 0)+
  geom_vline(xintercept = data$x[1]-12)+
  geom_vline(xintercept = data$x[1]+24)+
  xlab("time relative to farrowing")  +
  ggtitle("24h window w/ 15 min steps with first 12hrs removed and Kalman Filtering Sow 34083") +
  theme_bw()

# --- Create the training baseline lookup ---
train_lookup <- training_data %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  # We store the variance values as a list so they can be passed to the filter function
  summarise(train_y = list(m_var_24h), .groups = "drop")

Adding CI

library(tidyverse)

Warning: package 'tidyverse' was built under R version 4.4.3

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ readr   2.1.5     ✔ tibble  3.2.1
✔ stringr 1.5.1     
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ data.table::between()   masks dplyr::between()
✖ dplyr::filter()         masks stats::filter()
✖ data.table::first()     masks dplyr::first()
✖ data.table::hour()      masks lubridate::hour()
✖ data.table::isoweek()   masks lubridate::isoweek()
✖ dplyr::lag()            masks stats::lag()
✖ data.table::last()      masks dplyr::last()
✖ data.table::mday()      masks lubridate::mday()
✖ data.table::minute()    masks lubridate::minute()
✖ data.table::month()     masks lubridate::month()
✖ data.table::quarter()   masks lubridate::quarter()
✖ data.table::second()    masks lubridate::second()
✖ data.table::transpose() masks purrr::transpose()
✖ data.table::wday()      masks lubridate::wday()
✖ data.table::week()      masks lubridate::week()
✖ data.table::yday()      masks lubridate::yday()
✖ data.table::year()      masks lubridate::year()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(zoo)

Attaching package: 'zoo'

The following objects are masked from 'package:data.table':

    yearmon, yearqtr

The following objects are masked from 'package:base':

    as.Date, as.Date.numeric

library(ggplot2)

  # adjust to taste

plot_smooth <- plot_smooth %>%
  mutate(
    resid    = y - y_kalman,
    se       = rollapply(resid, width = 96, FUN = sd, fill = NA, align = "right"),
    ci_lower = y_kalman - qt(0.975, df = 95) * se,
    ci_upper = y_kalman + qt(0.975, df = 95) * se
  ) %>% 
  select(x,y,y_kalman,ci_lower,ci_upper)

ggplot(plot_smooth, aes(x = x)) +
  geom_line(aes(y = y_kalman), color = "red") +
  geom_line(aes(y = y),        color = "blue") +
  geom_line(aes(y = ci_upper), color = "orange") +
  geom_line(aes(y = ci_lower), color = "brown")

Warning: Removed 95 rows containing missing values or values outside the scale range
(`geom_line()`).

Warning: Removed 95 rows containing missing values or values outside the scale range
(`geom_line()`).

Function for kalman filtering

# --- Kalman filter function ---
fit_kalman <- function(train_y, full_y) {
  
  a0 <- train_y[1]
  P0 <- matrix(1)
  dt <- ct <- matrix(0)
  Zt <- Tt <- matrix(1)
  
  fit <- optim(
    c(HHt = var(train_y, na.rm = TRUE) * .5,
      GGt = var(train_y, na.rm = TRUE) * .5),
    fn = function(par, ...)
      -fkf(HHt = matrix(0.000001), GGt = matrix(0.001), ...)$logLik,
    yt = rbind(train_y), a0 = a0, P0 = P0,
    dt = dt, ct = ct, Zt = Zt, Tt = Tt
  )
  
  fkf(
    a0 = a0, P0 = P0, dt = dt, ct = ct, Tt = Tt, Zt = Zt,
    HHt = matrix(0.000001),
    GGt = matrix(0.001),
    yt  = rbind(full_y)
  )$att[1, ] |> as.numeric()
}

# --- apply to all sows ---
results95 <- sows_greater_60 %>%
  unnest(data) %>%
  filter(sensor_ts >= (start_time + hours(12))) %>%
  group_by(Sow_ID) %>%
  nest() %>%
  left_join(train_lookup, by = "Sow_ID") %>%
  mutate(
    data = map2(train_y, data, function(train, df) {
      df %>% 
        mutate(
          # 1. Standard Kalman and CI calculations
          y_kalman = fit_kalman(train, m_var_24h),
          resid    = m_var_24h - y_kalman,
          se       = zoo::rollapply(resid, width = 96, FUN = sd, fill = NA, align = "right"),
          ci_lower = y_kalman - qt(0.975, df = 95) * se,
          ci_upper = y_kalman + qt(0.975, df = 95) * se,
          
          # 2. Define the breach flags
          above_upper = m_var_24h > ci_upper,
          below_lower = m_var_24h < ci_lower,
          
          # 3. Create the index for the FIRST spike
          # We use if_else to handle sows that never spike (returning a very large number)
          first_spike_idx = if(any(above_upper, na.rm = TRUE)) min(which(above_upper == TRUE)) else Inf,
          
          # 4. Now we can safely reference first_spike_idx
          is_post_spike = row_number() > first_spike_idx,
          crossed_lower_after_upper = is_post_spike & below_lower,
          
          # 5. Status Column
          violation_type = case_when(
            m_var_24h > ci_upper ~ 1,
            m_var_24h < ci_lower ~ -1,
            TRUE ~ 0
          )
        )
    })
  ) %>%
  select(-train_y)

results95 %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  group_walk(~ {
    baseline_end <- min(.x$time_to_farrow, na.rm = TRUE) - 12
    p <- ggplot(.x, aes(x = time_to_farrow)) +
      geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), fill = "steelblue", alpha = 0.2) +
      geom_line(aes(y = m_var_24h), color = "black", alpha = 0.3) +
      geom_line(aes(y = y_kalman),  color = "blue", linewidth = 0.8) +
      geom_vline(xintercept = 0, linetype = "dashed") +
      geom_vline(xintercept = baseline_end)+
      expand_limits(y = 0) +
      theme_bw() +
      labs(
        title = paste("Sow", .y$Sow_ID),
        x     = "Hours relative to farrowing",
        y     = "24h variance of magnitude"
      )
    print(p)
  })

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

For a total of 94 sows, the kalman filter was applied with a manualy fixed HHt and GGt matrices for smoothing

Quantitative results

alarm_performance_all95 <- results95 %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  summarise(
    # First Alarm: First time crossing ABOVE Upper CI
    first_alarm = min(time_to_farrow[above_upper == TRUE], na.rm = TRUE),
    
    # Second Alarm: First time crossing BELOW Lower CI (after the first spike)
    second_alarm = min(time_to_farrow[crossed_lower_after_upper == TRUE], na.rm = TRUE),
    
    .groups = "drop"
  )

Warning: There were 59 warnings in `summarise()`.
The first warning was:
ℹ In argument: `first_alarm = min(time_to_farrow[above_upper == TRUE], na.rm =
  TRUE)`.
ℹ In group 13: `Sow_ID = 28442`.
Caused by warning in `min()`:
! no non-missing arguments to min; returning Inf
ℹ Run `dplyr::last_dplyr_warnings()` to see the 58 remaining warnings.

# View the full list
print(alarm_performance_all95)

# A tibble: 94 × 3
   Sow_ID first_alarm second_alarm
   <fct>        <dbl>        <dbl>
 1 28918       -12.5         Inf  
 2 26441       -27.8         Inf  
 3 31719         7.75        Inf  
 4 34585       -45.5         -21.5
 5 29938       -15.8          22.8
 6 29934       -41.8         Inf  
 7 34083       -32.2         -31.2
 8 34042        32           Inf  
 9 26240       -14           Inf  
10 28738         9.25         35.2
# ℹ 84 more rows

# 1. Ensure the per-sow summary is fresh
dt_perf_all95 <- results95 %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  summarise(
    first_alarm  = min(time_to_farrow[above_upper == TRUE], na.rm = TRUE),
    second_alarm = min(time_to_farrow[crossed_lower_after_upper == TRUE], na.rm = TRUE),
    .groups = "drop"
  )

Warning: There were 59 warnings in `summarise()`.
The first warning was:
ℹ In argument: `first_alarm = min(time_to_farrow[above_upper == TRUE], na.rm =
  TRUE)`.
ℹ In group 13: `Sow_ID = 28442`.
Caused by warning in `min()`:
! no non-missing arguments to min; returning Inf
ℹ Run `dplyr::last_dplyr_warnings()` to see the 58 remaining warnings.

# 2. Build the final results table with lead-time tiers
quant_res_final95 <- data.table(
  # --- OVERALL COUNTS ---
  n_total_sows               = nrow(dt_perf_all95),
  n_first_alarm              = sum(is.finite(dt_perf_all95$first_alarm)),
  n_no_first_alarm           = sum(!is.finite(dt_perf_all95$first_alarm)),
  n_second_alarm             = sum(is.finite(dt_perf_all95$second_alarm)),
  n_no_second_alarm          = sum(!is.finite(dt_perf_all95$second_alarm)),
  
  # --- MEANS ---
  mean_first_alarm           = mean(dt_perf_all95$first_alarm[is.finite(dt_perf_all95$first_alarm)], na.rm = TRUE),
  mean_second_alarm          = mean(dt_perf_all95$second_alarm[is.finite(dt_perf_all95$second_alarm)], na.rm = TRUE),

  # --- FIRST ALARM DETAILED BREAKDOWN (Sums to 78) ---
  n_first_alarm_gt_48h_before = sum(dt_perf_all95$first_alarm < -48, na.rm = TRUE),
  n_first_alarm_24_to_48h_pre = sum(dt_perf_all95$first_alarm < -24 & dt_perf_all95$first_alarm >= -48, na.rm = TRUE),
  n_first_alarm_0_to_24h_pre  = sum(dt_perf_all95$first_alarm < 0 & dt_perf_all95$first_alarm >= -24, na.rm = TRUE),
  n_first_alarm_on_FD         = sum(dt_perf_all95$first_alarm >= 0 & dt_perf_all95$first_alarm <= 24, na.rm = TRUE),
  n_first_alarm_after_FD      = sum(dt_perf_all95$first_alarm > 24 & is.finite(dt_perf_all95$first_alarm), na.rm = TRUE),
  
  # --- SECOND ALARM BREAKDOWN (Sums to 51) ---
  n_second_alarm_before_FD    = sum(dt_perf_all95$second_alarm < 0, na.rm = TRUE),
  n_second_alarm_on_FD        = sum(dt_perf_all95$second_alarm >= 0 & dt_perf_all95$second_alarm <= 24, na.rm = TRUE),
  n_second_alarm_after_FD     = sum(dt_perf_all95$second_alarm > 24 & is.finite(dt_perf_all95$second_alarm), na.rm = TRUE)
)

# Transpose for final presentation
print(t(quant_res_final95))

                                 [,1]
n_total_sows                 94.00000
n_first_alarm                78.00000
n_no_first_alarm             16.00000
n_second_alarm               51.00000
n_no_second_alarm            43.00000
mean_first_alarm            -21.72756
mean_second_alarm             2.27451
n_first_alarm_gt_48h_before  11.00000
n_first_alarm_24_to_48h_pre  27.00000
n_first_alarm_0_to_24h_pre   27.00000
n_first_alarm_on_FD           7.00000
n_first_alarm_after_FD        6.00000
n_second_alarm_before_FD     22.00000
n_second_alarm_on_FD         17.00000
n_second_alarm_after_FD      12.00000

Results:

For a 95% confidence interval with a total of 94 sows with greater than 60 hours of data before the date of farrowing. There was 78 sows that raised a first alarm and of those 78 with an average of 21.7 hours before the date of farrowing, 51 then raised a second alarm with an average of 2.2 hours into the farrowing date. There were 65 sows that raised the first alarm before the date of farrowing and 7 sows during the farrowing date and 6 sows after the date of farrowing. There were 22 sows that raised a second alarm before the date of farrowing, 17 that raised a second alarm on the date of farrowing and 12 with a second alarm after the date of farrowing. There is 27 sows that raised the first alarm less than 24hrs before the date of farrowing, 27 sows 24-48hrs before the farrowing date, 11 sows that raised the first alarm greater than 48hrs before the date of farrowing, 7 sows that raised the first alarm on the farrowing date, and 6 sows that raised a first alarm after the date of farrowing

Now to look at a 68% CI interval

# --- apply to all sows ---
results68 <- sows_greater_60 %>%
  unnest(data) %>%
  filter(sensor_ts >= (start_time + hours(12))) %>%
  group_by(Sow_ID) %>%
  nest() %>%
  left_join(train_lookup, by = "Sow_ID") %>%
  mutate(
    data = map2(train_y, data, function(train, df) {
      df %>% 
        mutate(
          # 1. Standard Kalman and CI calculations
          y_kalman = fit_kalman(train, m_var_24h),
          resid    = m_var_24h - y_kalman,
          se       = zoo::rollapply(resid, width = 96, FUN = sd, fill = NA, align = "right"),
          ci_lower = y_kalman - qt(0.84, df = 95) * se,
          ci_upper = y_kalman + qt(0.84, df = 95) * se,
          
          # 2. Define the breach flags
          above_upper = m_var_24h > ci_upper,
          below_lower = m_var_24h < ci_lower,
          
          # 3. Create the index for the FIRST spike
          # We use if_else to handle sows that never spike (returning a very large number)
          first_spike_idx = if(any(above_upper, na.rm = TRUE)) min(which(above_upper == TRUE)) else Inf,
          
          # 4. Now we can safely reference first_spike_idx
          is_post_spike = row_number() > first_spike_idx,
          crossed_lower_after_upper = is_post_spike & below_lower,
          
          # 5. Status Column
          violation_type = case_when(
            m_var_24h > ci_upper ~ 1,
            m_var_24h < ci_lower ~ -1,
            TRUE ~ 0
          )
        )
    })
  ) %>%
  select(-train_y)

results68 %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  group_walk(~ {
    baseline_end <- min(.x$time_to_farrow, na.rm = TRUE) - 12
    p <- ggplot(.x, aes(x = time_to_farrow)) +
      geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), fill = "steelblue", alpha = 0.2) +
      geom_line(aes(y = m_var_24h), color = "black", alpha = 0.3) +
      geom_line(aes(y = y_kalman),  color = "blue", linewidth = 0.8) +
      geom_vline(xintercept = 0, linetype = "dashed") +
      geom_vline(xintercept = baseline_end)+
      expand_limits(y = 0) +
      theme_bw() +
      labs(
        title = paste("Sow", .y$Sow_ID),
        x     = "Hours relative to farrowing",
        y     = "24h variance of magnitude"
      )
    print(p)
  })

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

For a total of 94 sows, the kalman filter was applied with a manualy fixed HHt and GGt matrices for smoothing

Quantitative results

alarm_performance_all68 <- results68 %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  summarise(
    # First Alarm: First time crossing ABOVE Upper CI
    first_alarm = min(time_to_farrow[above_upper == TRUE], na.rm = TRUE),
    
    # Second Alarm: First time crossing BELOW Lower CI (after the first spike)
    second_alarm = min(time_to_farrow[crossed_lower_after_upper == TRUE], na.rm = TRUE),
    
    .groups = "drop"
  )

Warning: There were 22 warnings in `summarise()`.
The first warning was:
ℹ In argument: `first_alarm = min(time_to_farrow[above_upper == TRUE], na.rm =
  TRUE)`.
ℹ In group 21: `Sow_ID = 28225`.
Caused by warning in `min()`:
! no non-missing arguments to min; returning Inf
ℹ Run `dplyr::last_dplyr_warnings()` to see the 21 remaining warnings.

# View the full list
print(alarm_performance_all68)

# A tibble: 94 × 3
   Sow_ID first_alarm second_alarm
   <fct>        <dbl>        <dbl>
 1 28918       -12.5          11.5
 2 26441       -32.2          31.8
 3 31719       -46.8         -24.2
 4 34585       -45.5         -21.5
 5 29938       -17.5          17.8
 6 29934       -48.5          11  
 7 34083       -36.2         -31.2
 8 34042       -26           -17.2
 9 26240       -14.2          18.5
10 28738         3.25         34.2
# ℹ 84 more rows

library(data.table)

# 1. Ensure the per-sow summary is fresh
dt_perf_all68 <- results68 %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  summarise(
    first_alarm  = min(time_to_farrow[above_upper == TRUE], na.rm = TRUE),
    second_alarm = min(time_to_farrow[crossed_lower_after_upper == TRUE], na.rm = TRUE),
    .groups = "drop"
  )

Warning: There were 22 warnings in `summarise()`.
The first warning was:
ℹ In argument: `first_alarm = min(time_to_farrow[above_upper == TRUE], na.rm =
  TRUE)`.
ℹ In group 21: `Sow_ID = 28225`.
Caused by warning in `min()`:
! no non-missing arguments to min; returning Inf
ℹ Run `dplyr::last_dplyr_warnings()` to see the 21 remaining warnings.

# 2. Build the final results table with lead-time tiers
quant_res_final68 <- data.table(
  # --- OVERALL COUNTS ---
  n_total_sows               = nrow(dt_perf_all68),
  n_first_alarm              = sum(is.finite(dt_perf_all68$first_alarm)),
  n_no_first_alarm           = sum(!is.finite(dt_perf_all68$first_alarm)),
  n_second_alarm             = sum(is.finite(dt_perf_all68$second_alarm)),
  n_no_second_alarm          = sum(!is.finite(dt_perf_all68$second_alarm)),
  
  # --- MEANS ---
  mean_first_alarm           = mean(dt_perf_all68$first_alarm[is.finite(dt_perf_all68$first_alarm)], na.rm = TRUE),
  mean_second_alarm          = mean(dt_perf_all68$second_alarm[is.finite(dt_perf_all68$second_alarm)], na.rm = TRUE),

  # --- FIRST ALARM DETAILED BREAKDOWN (Sums to 78) ---
  n_first_alarm_gt_48h_before = sum(dt_perf_all68$first_alarm < -48, na.rm = TRUE),
  n_first_alarm_24_to_48h_pre = sum(dt_perf_all68$first_alarm < -24 & dt_perf_all68$first_alarm >= -48, na.rm = TRUE),
  n_first_alarm_0_to_24h_pre  = sum(dt_perf_all68$first_alarm < 0 & dt_perf_all68$first_alarm >= -24, na.rm = TRUE),
  n_first_alarm_on_FD         = sum(dt_perf_all68$first_alarm >= 0 & dt_perf_all68$first_alarm <= 24, na.rm = TRUE),
  n_first_alarm_after_FD      = sum(dt_perf_all68$first_alarm > 24 & is.finite(dt_perf_all68$first_alarm), na.rm = TRUE),
  
  # --- SECOND ALARM BREAKDOWN (Sums to 51) ---
  n_second_alarm_before_FD    = sum(dt_perf_all68$second_alarm < 0, na.rm = TRUE),
  n_second_alarm_on_FD        = sum(dt_perf_all68$second_alarm >= 0 & dt_perf_all68$second_alarm <= 24, na.rm = TRUE),
  n_second_alarm_after_FD     = sum(dt_perf_all68$second_alarm > 24 & is.finite(dt_perf_all68$second_alarm), na.rm = TRUE)
)

# Transpose for final presentation
print(t(quant_res_final68))

                                  [,1]
n_total_sows                 94.000000
n_first_alarm                87.000000
n_no_first_alarm              7.000000
n_second_alarm               79.000000
n_no_second_alarm            15.000000
mean_first_alarm            -29.597701
mean_second_alarm            -7.344937
n_first_alarm_gt_48h_before  21.000000
n_first_alarm_24_to_48h_pre  34.000000
n_first_alarm_0_to_24h_pre   24.000000
n_first_alarm_on_FD           3.000000
n_first_alarm_after_FD        5.000000
n_second_alarm_before_FD     46.000000
n_second_alarm_on_FD         26.000000
n_second_alarm_after_FD       7.000000

Now trying 90% CI

# --- apply to all sows ---
results90 <- sows_greater_60 %>%
  unnest(data) %>%
  filter(sensor_ts >= (start_time + hours(12))) %>%
  group_by(Sow_ID) %>%
  nest() %>%
  left_join(train_lookup, by = "Sow_ID") %>%
  mutate(
    data = map2(train_y, data, function(train, df) {
      df %>% 
        mutate(
          # 1. Standard Kalman and CI calculations
          y_kalman = fit_kalman(train, m_var_24h),
          resid    = m_var_24h - y_kalman,
          se       = zoo::rollapply(resid, width = 96, FUN = sd, fill = NA, align = "right"),
          ci_lower = y_kalman - qt(0.95, df = 95) * se,
          ci_upper = y_kalman + qt(0.95, df = 95) * se,
          
          # 2. Define the breach flags
          above_upper = m_var_24h > ci_upper,
          below_lower = m_var_24h < ci_lower,
          
          # 3. Create the index for the FIRST spike
          # We use if_else to handle sows that never spike (returning a very large number)
          first_spike_idx = if(any(above_upper, na.rm = TRUE)) min(which(above_upper == TRUE)) else Inf,
          
          # 4. Now we can safely reference first_spike_idx
          is_post_spike = row_number() > first_spike_idx,
          crossed_lower_after_upper = is_post_spike & below_lower,
          
          # 5. Status Column
          violation_type = case_when(
            m_var_24h > ci_upper ~ 1,
            m_var_24h < ci_lower ~ -1,
            TRUE ~ 0
          )
        )
    })
  ) %>%
  select(-train_y)

results90 %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  group_walk(~ {
    baseline_end <- min(.x$time_to_farrow, na.rm = TRUE) - 12
    p <- ggplot(.x, aes(x = time_to_farrow)) +
      geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), fill = "steelblue", alpha = 0.2) +
      geom_line(aes(y = m_var_24h), color = "black", alpha = 0.3) +
      geom_line(aes(y = y_kalman),  color = "blue", linewidth = 0.8) +
      geom_vline(xintercept = 0, linetype = "dashed") +
      geom_vline(xintercept = baseline_end)+
      expand_limits(y = 0) +
      theme_bw() +
      labs(
        title = paste("Sow", .y$Sow_ID),
        x     = "Hours relative to farrowing",
        y     = "24h variance of magnitude"
      )
    print(p)
  })

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning
-Inf

For a total of 94 sows, the kalman filter was applied with a manualy fixed HHt and GGt matrices for smoothing

Quantitative results

alarm_performance_all90 <- results90 %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  summarise(
    # First Alarm: First time crossing ABOVE Upper CI
    first_alarm = min(time_to_farrow[above_upper == TRUE], na.rm = TRUE),
    
    # Second Alarm: First time crossing BELOW Lower CI (after the first spike)
    second_alarm = min(time_to_farrow[crossed_lower_after_upper == TRUE], na.rm = TRUE),
    
    .groups = "drop"
  )

Warning: There were 41 warnings in `summarise()`.
The first warning was:
ℹ In argument: `first_alarm = min(time_to_farrow[above_upper == TRUE], na.rm =
  TRUE)`.
ℹ In group 21: `Sow_ID = 28225`.
Caused by warning in `min()`:
! no non-missing arguments to min; returning Inf
ℹ Run `dplyr::last_dplyr_warnings()` to see the 40 remaining warnings.

# View the full list
print(alarm_performance_all90)

# A tibble: 94 × 3
   Sow_ID first_alarm second_alarm
   <fct>        <dbl>        <dbl>
 1 28918        -12.5        Inf  
 2 26441        -30.5        Inf  
 3 31719        -26.8        -23.8
 4 34585        -45.5        -21.5
 5 29938        -16.5         22  
 6 29934        -41.8         11  
 7 34083        -32.8        -31.2
 8 34042         32          Inf  
 9 26240        -14           19.8
10 28738          6.5         34.5
# ℹ 84 more rows

# 1. Ensure the per-sow summary is fresh
dt_perf_all90 <- results90 %>%
  unnest(data) %>%
  group_by(Sow_ID) %>%
  summarise(
    first_alarm  = min(time_to_farrow[above_upper == TRUE], na.rm = TRUE),
    second_alarm = min(time_to_farrow[crossed_lower_after_upper == TRUE], na.rm = TRUE),
    .groups = "drop"
  )

Warning: There were 41 warnings in `summarise()`.
The first warning was:
ℹ In argument: `first_alarm = min(time_to_farrow[above_upper == TRUE], na.rm =
  TRUE)`.
ℹ In group 21: `Sow_ID = 28225`.
Caused by warning in `min()`:
! no non-missing arguments to min; returning Inf
ℹ Run `dplyr::last_dplyr_warnings()` to see the 40 remaining warnings.

# 2. Build the final results table with lead-time tiers
quant_res_final90 <- data.table(
  # --- OVERALL COUNTS ---
  n_total_sows               = nrow(dt_perf_all90),
  n_first_alarm              = sum(is.finite(dt_perf_all90$first_alarm)),
  n_no_first_alarm           = sum(!is.finite(dt_perf_all90$first_alarm)),
  n_second_alarm             = sum(is.finite(dt_perf_all90$second_alarm)),
  n_no_second_alarm          = sum(!is.finite(dt_perf_all90$second_alarm)),
  
  # --- MEANS ---
  mean_first_alarm           = mean(dt_perf_all90$first_alarm[is.finite(dt_perf_all90$first_alarm)], na.rm = TRUE),
  mean_second_alarm          = mean(dt_perf_all90$second_alarm[is.finite(dt_perf_all90$second_alarm)], na.rm = TRUE),

  # --- FIRST ALARM DETAILED BREAKDOWN (Sums to 78) ---
  n_first_alarm_gt_48h_before = sum(dt_perf_all90$first_alarm < -48, na.rm = TRUE),
  n_first_alarm_24_to_48h_pre = sum(dt_perf_all90$first_alarm < -24 & dt_perf_all90$first_alarm >= -48, na.rm = TRUE),
  n_first_alarm_0_to_24h_pre  = sum(dt_perf_all90$first_alarm < 0 & dt_perf_all90$first_alarm >= -24, na.rm = TRUE),
  n_first_alarm_on_FD         = sum(dt_perf_all90$first_alarm >= 0 & dt_perf_all90$first_alarm <= 24, na.rm = TRUE),
  n_first_alarm_after_FD      = sum(dt_perf_all90$first_alarm > 24 & is.finite(dt_perf_all90$first_alarm), na.rm = TRUE),
  
  # --- SECOND ALARM BREAKDOWN (Sums to 51) ---
  n_second_alarm_before_FD    = sum(dt_perf_all90$second_alarm < 0, na.rm = TRUE),
  n_second_alarm_on_FD        = sum(dt_perf_all90$second_alarm >= 0 & dt_perf_all90$second_alarm <= 24, na.rm = TRUE),
  n_second_alarm_after_FD     = sum(dt_perf_all90$second_alarm > 24 & is.finite(dt_perf_all90$second_alarm), na.rm = TRUE)
)

# Transpose for final presentation
print(t(quant_res_final90))

                                   [,1]
n_total_sows                 94.0000000
n_first_alarm                83.0000000
n_no_first_alarm             11.0000000
n_second_alarm               64.0000000
n_no_second_alarm            30.0000000
mean_first_alarm            -23.4698795
mean_second_alarm             0.3398438
n_first_alarm_gt_48h_before  16.0000000
n_first_alarm_24_to_48h_pre  28.0000000
n_first_alarm_0_to_24h_pre   27.0000000
n_first_alarm_on_FD           4.0000000
n_first_alarm_after_FD        8.0000000
n_second_alarm_before_FD     28.0000000
n_second_alarm_on_FD         21.0000000
n_second_alarm_after_FD      15.0000000

Saving data and moving to a new qmd to compare results

save(results90,results68,results95,quant_res_final68,quant_res_final90,quant_res_final95, file='C:/Users/ipberg/OneDrive - Iowa State University/Spring26/Sensors/4-13.RData')

Add 2 things confidence interval bands alarms from there

Apply this to all animals is it done in a friendly that can be applied to all animals

can a data set to quantify results

CI sow specific data - check with graphics from data for add alarm loop through all 50

Read paper BLUP is a good thing