MJacobs1985

df_nest<-day%>%
  dplyr::select(date_interval,location_identification, egg_count)%>%
  group_by(location_identification)%>%
  filter(!location_identification%in%c(4))%>%
  nest()
df_ts<-df_nest%>%mutate(data.ts = map(.x       = data, 
                              .f       = tk_ts, 
                              select   = -date_interval, 
                              start    = 2021,
                              freq     = 12))

MJacobs1985

ggplot(day, aes(x=date_interval))+
  geom_point(aes(y=egg_count/animals_actual, colour=as.factor(location_identification)))+
  theme_bw()+  theme(legend.position="none")

ggplot(day, aes(x=date_interval))+
  geom_point(aes(y=(egg_count*average_weight_eggs)/1000000, 
                 colour=as.factor(location_identification)))+
  theme_bw()+  theme(legend.position="none")

ggplot(day, aes(x=date_interval))+

MJacobs1985

DataExplorer::plot_missing(day)
DataExplorer::plot_missing(hour)

g1<-ggplot(Norm, aes(x=`Age in Weeks`))+
  geom_point(aes(y=`Cumm egg prod (HH)`))+
  geom_point(aes(y=`Cumm egg produced (HD)`), col="red")+
  theme_bw()
g2<-ggplot(Norm, aes(x=`Age in Weeks`))+
  geom_point(aes(y=`Egg weight (g)`))+
  geom_point(aes(y=`Egg mass per day (g)`), col="red")+

MJacobs1985

rm(list = ls())
options("scipen"=100, "digits"=10)

#### LIBRARIES ####
library(readr)
library(readxl)
library(DataExplorer)
library(tidyverse)
library(ggExtra)
library(gridExtra)

MJacobs1985

dfnew%>%filter(!is.na(comp3_provided))%>%
  dplyr::select(ogn_id,curve_day,
                feed_provided,
                comp1_provided,
                comp2_provided,comp3_provided)

### Forecast summed feed curves over time
testdata2<-dfnew%>%
  group_by(date_interval, lcn_id, newField) %>%
  mutate(comp2_provided   = ifelse(is.na(comp2_provided),0,comp2_provided),

MJacobs1985

m1_gam_comp1 <- gamm4(comp1_provided~
                        s(curve_day)+
                        s(animals_actual)+
                        s(feed_provided_cum)+
                        s(comp2_provided)+
                        ogn_id+
                        newField,
                      data=training,
                      random=~(1+curve_day|lcn_id))
m1_gam_comp2 <- gamm4(comp2_provided~

MJacobs1985

m1_gam_comp2 <- gamm4(comp2_provided~
                        s(curve_day)+
                        s(animals_actual)+
                        s(feed_provided_cum)+
                        s(comp1_provided)+
                        ogn_id+
                        newField,
                      data=training,
                      random=~(1+curve_day|lcn_id))
par(mfrow = c(2, 2))

MJacobs1985

set.seed(1123)
inTraining <- createDataPartition(dfnew$newField, 
                                  p = .50, 
                                  list = FALSE)
training   <- dfnew[ inTraining,]
testing    <- dfnew[-inTraining,]

m1_gam_comp1 <- gamm4(comp1_provided~
                        s(curve_day)+
                        s(animals_actual)+

MJacobs1985

dfnew<-as.data.frame(dfnew)
inTraining <- createDataPartition(dfnew$newField, 
                                  p = .75, 
                                  list = FALSE)

training   <- dfnew[ inTraining,]
testing    <- dfnew[-inTraining,]

training%>%
  filter(lcn_id%in%(901:906))%>%

MJacobs1985

dfnew2%>%
  filter(lcn_id%in%(901:909))%>%
  xyplot(feed_provided_cum ~ curve_day|lcn_id*newField,., 
       strip = FALSE, aspect = "xy", pch = 16, cex=0.5, grid = TRUE,
       panel = function(x, y, ..., fit, subscripts) {
         panel.xyplot(x, y, ...)
         ypred1 <- fitted(fit.1)[subscripts]
         panel.xyplot(x, ypred1, lty=1, col = "orange")
         ypred2 <- fitted(fit.2)[subscripts]
         panel.xyplot(x, ypred2, lty=1, col = "purple")