gistfile1.txt

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))
plot(m1_gam_comp2$gam)
m1_gam_comp2_aug<-augment(m1_gam_comp2$gam, confint=TRUE)
m1_gam_comp2_aug%>%filter(lcn_id%in%(901:909))%>%
  ggplot(., aes(x=comp2_provided, y=.fitted, colour=.resid))+
  geom_point()+
  theme_bw()+
  facet_wrap(~lcn_id+newField)+
  geom_abline(intercept=0, slope=1, lty=2, col="black")+
  scale_color_viridis_c(option="turbo")
testing$comp2_provided<-ifelse(is.na(testing$comp2_provided),0, 
                               testing$comp2_provided)
testing$pred<-predict(m1_gam_comp2$gam, 
                      newdata=testing)
testing%>%filter(lcn_id%in%(901:920))%>%
  ggplot(., aes(x=curve_day, group=newField))+
  geom_point(aes(y=comp2_provided), alpha=0.5, color="black")+
  geom_point(aes(y=pred, colour=newField), lty=2, lwd=1)+
  theme_bw()+
  facet_wrap(~lcn_id+newField)
g1<-testing%>%
  group_by(ogn_id, curve_day)%>%
  summarize(prov = sum(comp2_provided),
            pred = sum(pred))%>%
  ggplot(., aes(x=curve_day))+
  geom_area(aes(y=prov, 
                fill=as.factor(ogn_id)), 
            alpha=0.5)+
  geom_point(aes(y=pred), alpha=0.5)+
  geom_line(aes(y=pred))+
  theme_bw()+
  labs(x="Day in Curve", 
       y="comp2", 
       fill="Farm", 
       title="Predicted vs Observed on Component 2 - TEST data")+
  facet_wrap(~ogn_id)
g2<-m1_gam_comp2_aug%>%
  group_by(ogn_id,curve_day)%>%
  summarize(prov = sum(comp2_provided, na.rm = TRUE),
            pred = sum(.fitted, na.rm = TRUE))%>%
  ggplot(., aes(x=curve_day))+
  geom_area(aes(y=prov, 
                fill=as.factor(ogn_id)), 
            alpha=0.5)+
  geom_point(aes(y=pred), alpha=0.5)+
  geom_line(aes(y=pred))+
  theme_bw()+
  labs(x="Day in Curve", 
       y="comp2", 
       fill="Farm", 
       title="Predicted vs Observed on Component 2 - TRAIN data")+
  facet_wrap(~ogn_id)
grid.arrange(g2,g1,ncol=1)