atyre2

## function to do graphical checks in one hit, given a model
## assumes that broom and gridExtra are available
## and that the model object has tidiers available
check_assumptions <- function(x, ...){
  if (inherits(x,"gam")){
    stop("check_assumptions doesn't work with gam() objects")
  }
  
  my_theme <- ggplot2::theme_classic() + 
    ggplot2::theme(text=ggplot2::element_text(size=rel(4)))

atyre2

# from https://github.com/ha0ye/personal-website/blob/master/figure%20code/causality-figures.R

library(tidyverse)
library(gganimate) # devtools::install_github("thomasp85/gganimate")

## generate model data ----

# params
t_max <- 100
h <- 0.01 # step size

atyre2

#' # Lake Ontario Cormorant matrix
cormorant_post_pop <- function(S = c(0.5, 0.87, 0.88, 0.89), 
                               B = c(0, 0.5, 0.99), 
                               F_ = c(0, 1.6, 2.4),
                               R = 0.5){
  A <- matrix(0, nrow=3, ncol=3)
  A[2,1] <- S[1]
  A[3,2:3] <- S[2:3] # ends up slightly different than pre-breeding matrix
  A[1,] <- B*F_*S[2:4]*R
  A

atyre2

dieet <- structure(list(Jaar = c(2013L, 2013L, 2013L, 2013L, 2013L, 2013L, 
2014L, 2014L, 2014L, 2014L, 2014L, 2015L, 2015L, 2015L, 2013L, 
2013L, 2013L, 2013L, 2013L, 2013L, 2013L, 2013L, 2013L, 2014L, 
2014L, 2014L, 2014L, 2014L, 2014L, 2014L, 2015L, 2015L, 2015L, 
2015L, 2015L, 2015L), Kuikenweek = c("Week 1", "Week -1", "Week 2", 
"Week -2", "Week 3", "Week -3", "Week 1", "Week 2", "Week -2", 
"Week 3", "Week -3", "Week -1", "Week -2", "Week -3", "Week 1", 
"Week -1", "Week 2", "Week -2", "Week 3", "Week -3", "Week 4", 
"Week 5", "Week 6", "Week 1", "Week -1", "Week 2", "Week 3", 
"Week -3", "Week 4", "Week 5", "Week 1", "Week -1", "Week 2", 

atyre2

library(rvest)
library(tidyverse)

webpage <- read_html("http://snr.unl.edu/aboutus/who/people/wellness/flightascendingsummary.asp")

tbls <- html_nodes(webpage, "table") %>%
  html_table() %>%
  .[[1]] %>%
  separate(Participant, into = c("Last Name", "First Name"), sep = ",") %>%
  mutate(`Flights per day` = `Total Number of Flights` / `Total Number of Days`,

atyre2

library(rvest)
library(tidyverse)

webpage <- read_html("http://snr.unl.edu/aboutus/who/people/wellness/flightascendingsummary.asp")

tbls <- html_nodes(webpage, "table") %>%
  html_table() %>%
  .[[1]] %>%
  separate(Participant, into = c("Last Name", "First Name"), sep = ",") %>%
  mutate(`Flights per day` = `Total Number of Flights` / `Total Number of Days`,

atyre2

library(tidyverse)
xx <- read_table2(" 1    family1 0.8690661709 0.751969991 0.7188319163 0.8329348587 0.887535770
     2    family1 0.2191712498 0.378020942 0.8583653572 0.4347398167 0.463576331
                     3    family1 0.6999901198 0.641613472 0.3160536611 0.6978992328 0.648399558
                     4    family3 0.7372355657 0.538900863 0.7362216445 0.7604946053 0.822828410
                     5    family1 0.7096292505 0.598671265 0.2544588205 0.6133470555 0.607725992
                     6    family1 0.6360280758 0.566566355 0.5877906233 0.6997789522 0.546953210
                     7    family1 0.5991758243 0.616002671 0.8369985343 0.8206614441 0.818448365
                     8    family1 0.8904527922 0.716208368 0.9011345475 0.5670958449 0.727816837
                     9    family2 0.1452406058 0.207274998 0.4986928476 0.1686444875 0.189754451

atyre2

n <- 200;sig <- 2
x <- 0:(n-1)/(n-1)
f <- 0.2*x^11*(10*(1-x))^6+10*(10*x)^3*(1-x)^10
e <- rnorm(n,0,sig)
for (i in 2:n) e[i] <- 0.6*e[i-1] + e[i]
y <- f + e
op <- par(mfrow=c(2,2))
## Fit model with AR1 residuals
b <- gamm(y~s(x,k=20),correlation=corAR1())
plot(b$gam);lines(x,f-mean(f),col=2)