RELIVAL/CODE/functions.R

# Functions

get_msd=function(x) paste0("M = ", round(mean(x),2),", SD = ", round(sd(x),2))



define_contiguous<-function(mydat){
  #this function creates a dataset with one line per paired recordings, with session_id & next_session being the paired recs
  paired <- mydat %>% 
    #the next line sorts the table by child id then age
    arrange(child_id, age) %>% 
    #the next line keeps only one line per combination of experiment and child_id
    group_by(child_id) %>% 
    #the next line, mutate, defines 3 new variables in the dataset, n_rec, age_dist_next_rec, and next_session
    mutate( #n_rec = n(), #this var isn't used later
      age_dist_next_rec = lead(age) - age, 
      #this gets the diff corresponding cell in the preceding row and a given row
      next_session = lead(session_id)) %>% 
    #the next line keeps only recs that are less than 2 months away
    filter(age_dist_next_rec<2) 
  
  return(paired)
}

get_type<-function(mytab){
  temp_tab<-mytab
  temp_tab$Type="Output"
  temp_tab$Type[grep("fem", temp_tab$met)] <- "Female"
  temp_tab$Type[grep("mal", temp_tab$met)] <- "Male"
  temp_tab$Type[grep("och", temp_tab$met)] <- "Other children"
  temp_tab$Type[grep("adu", temp_tab$met)] <- "Adults"
  temp_tab$Type[grep("CTC", temp_tab$met)] <- "Adults"
  temp_tab$Type[grep("AWC", temp_tab$met)] <- "Adults"
  temp_tab$Type=factor(temp_tab$Type)
  temp_tab$Type
}

check_small_var<-function(x,y,i) round(x[i],3)==round(y[i],3) & round(y[i],3) == 0

fit_child_model<-function(dataframe, metric){
  # Fit formula where experiment is removed
  formula <- as.formula(paste0(metric, "~ age_s + (1|child_id)"))
  model <- lmer(formula, data=dataframe)
  
  return (model)
}

extract_chi_variables<-function(model){
  icc.result.mixed <- c(icc(model)$ICC_adjusted,icc(model)$ICC_conditional)
  icc.result.split <- c(as.data.frame(icc(model, by_group=TRUE))$ICC, NA)
  
  r2.nakagawa <- r2_nakagawa(model)
  
  ranefs_vars <- t(as.data.frame(VarCorr(model))["vcov"])
  ranefs_stdv <- t(as.data.frame(VarCorr(model))["sdcor"])
  # chi, NA, residual
  ranefs_vars <-c(ranefs_vars[1],NA,ranefs_vars[2])
  ranefs_stdv <-c(ranefs_stdv[1],NA,ranefs_stdv[2])
  # chi, NA
  ns<-c(unlist(summary(model)$n),NA)
  return (c(coefficients(summary(model))["age_s",],
            icc.result.mixed,
            icc.result.split,
            ranefs_vars,
            ranefs_stdv,
            r2.nakagawa,
            nobs(model),ns))
}

extract_full_variables<-function(model){
  icc.result.mixed <- c(icc(model)$ICC_adjusted,icc(model)$ICC_conditional)
  icc.result.split <- t(as.data.frame(icc(model, by_group=TRUE))$ICC)
  r2.nakagawa <- r2_nakagawa(model)
  
  ranefs_vars <- t(as.data.frame(VarCorr(model))["vcov"])
  ranefs_stdv <- t(as.data.frame(VarCorr(model))["sdcor"])
  
  ns<-t(data.frame(summary(model)$n))
  
  return (c(  coefficients(summary(model))["age_s",],
              icc.result.mixed,
              icc.result.split,
              ranefs_vars,
              ranefs_stdv,
              r2.nakagawa,
              nobs(model),ns))
}

new_fit_models<-function(dataframe, data_set, metric, age, fit_full = TRUE){
  #dataframe=mydat ; age = NA ; fit_full = TRUE
  #metric="wc_adu_ph"
  #to check age: dataframe<-thiscordata ; age <- thisage ; fit_full = TRUE
  iqr = quantile(dataframe[,metric],.75,na.rm=T)-quantile(dataframe[,metric],.25,na.rm=T)
  
  if(fit_full){
    # Fit full model
    formula <- as.formula(paste0(metric, "~ age_s + (1|experiment/child_id)"))
    model <- lmer(formula, data=dataframe)
  }
  
  if(fit_full && !isSingular(model)) # Fitted full model short circuit and
  {
    form="full"
    
    sw=shapiro.test(resid(model))$p
    mod_variables <- extract_full_variables(model)
    
    
  } else {
    model <- fit_child_model(dataframe, metric)
    if(!isSingular(model)){
      
      form = "no_exp"
      sw=shapiro.test(resid(model))$p
      mod_variables <- extract_chi_variables(model)
      
    } else {
      form='no_chi_effect'
      sw = NA
      mod_variables = c(
        NA,NA,NA, # c(coefficients(summary(model))["age_s",]
        NA,NA,    # icc.result.mixed
        NA,NA,    # icc.result.split
        NA,NA,NA, # ranefs_vars
        NA,NA,NA, # raners_std
        NA,NA,    # r2
        NA,NA,NA) # nobs (child, corpus), nobs
    }
  }
  
  icc.row = c(data_set, age, metric, iqr, mod_variables, form, sw)
  return (icc.row)
}