functions.R 4.3 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132
  1. # Functions
  2. get_msd=function(x) paste0("M = ", round(mean(x),2),", SD = ", round(sd(x),2))
  3. define_contiguous<-function(mydat){
  4. #this function creates a dataset with one line per paired recordings, with session_id & next_session being the paired recs
  5. paired <- mydat %>%
  6. #the next line sorts the table by child id then age
  7. arrange(child_id, age) %>%
  8. #the next line keeps only one line per combination of experiment and child_id
  9. group_by(child_id) %>%
  10. #the next line, mutate, defines 3 new variables in the dataset, n_rec, age_dist_next_rec, and next_session
  11. mutate( #n_rec = n(), #this var isn't used later
  12. age_dist_next_rec = lead(age) - age,
  13. #this gets the diff corresponding cell in the preceding row and a given row
  14. next_session = lead(session_id)) %>%
  15. #the next line keeps only recs that are less than 2 months away
  16. filter(age_dist_next_rec<2)
  17. return(paired)
  18. }
  19. get_type<-function(mytab){
  20. temp_tab<-mytab
  21. temp_tab$Type="Output"
  22. temp_tab$Type[grep("fem", temp_tab$met)] <- "Female"
  23. temp_tab$Type[grep("mal", temp_tab$met)] <- "Male"
  24. temp_tab$Type[grep("och", temp_tab$met)] <- "Other children"
  25. temp_tab$Type[grep("adu", temp_tab$met)] <- "Adults"
  26. temp_tab$Type[grep("CTC", temp_tab$met)] <- "Adults"
  27. temp_tab$Type[grep("AWC", temp_tab$met)] <- "Adults"
  28. temp_tab$Type=factor(temp_tab$Type)
  29. temp_tab$Type
  30. }
  31. check_small_var<-function(x,y,i) round(x[i],3)==round(y[i],3) & round(y[i],3) == 0
  32. fit_child_model<-function(dataframe, metric){
  33. # Fit formula where experiment is removed
  34. formula <- as.formula(paste0(metric, "~ age_s + (1|child_id)"))
  35. model <- lmer(formula, data=dataframe)
  36. return (model)
  37. }
  38. extract_chi_variables<-function(model){
  39. icc.result.mixed <- c(icc(model)$ICC_adjusted,icc(model)$ICC_conditional)
  40. icc.result.split <- c(as.data.frame(icc(model, by_group=TRUE))$ICC, NA)
  41. r2.nakagawa <- r2_nakagawa(model)
  42. ranefs_vars <- t(as.data.frame(VarCorr(model))["vcov"])
  43. ranefs_stdv <- t(as.data.frame(VarCorr(model))["sdcor"])
  44. # chi, NA, residual
  45. ranefs_vars <-c(ranefs_vars[1],NA,ranefs_vars[2])
  46. ranefs_stdv <-c(ranefs_stdv[1],NA,ranefs_stdv[2])
  47. # chi, NA
  48. ns<-c(unlist(summary(model)$n),NA)
  49. return (c(coefficients(summary(model))["age_s",],
  50. icc.result.mixed,
  51. icc.result.split,
  52. ranefs_vars,
  53. ranefs_stdv,
  54. r2.nakagawa,
  55. nobs(model),ns))
  56. }
  57. extract_full_variables<-function(model){
  58. icc.result.mixed <- c(icc(model)$ICC_adjusted,icc(model)$ICC_conditional)
  59. icc.result.split <- t(as.data.frame(icc(model, by_group=TRUE))$ICC)
  60. r2.nakagawa <- r2_nakagawa(model)
  61. ranefs_vars <- t(as.data.frame(VarCorr(model))["vcov"])
  62. ranefs_stdv <- t(as.data.frame(VarCorr(model))["sdcor"])
  63. ns<-t(data.frame(summary(model)$n))
  64. return (c( coefficients(summary(model))["age_s",],
  65. icc.result.mixed,
  66. icc.result.split,
  67. ranefs_vars,
  68. ranefs_stdv,
  69. r2.nakagawa,
  70. nobs(model),ns))
  71. }
  72. new_fit_models<-function(dataframe, data_set, metric, age, fit_full = TRUE){
  73. #dataframe=mydat ; age = NA ; fit_full = TRUE
  74. #to check age: dataframe<-thiscordata ; age <- thisage ; fit_full = TRUE
  75. iqr = quantile(dataframe[,metric],.75,na.rm=T)-quantile(dataframe[,metric],.25,na.rm=T)
  76. if(fit_full){
  77. # Fit full model
  78. formula <- as.formula(paste0(metric, "~ age_s + (1|experiment/child_id)"))
  79. model <- lmer(formula, data=dataframe)
  80. }
  81. if(fit_full && !isSingular(model)) # Fitted full model short circuit and
  82. {
  83. form="full"
  84. sw=shapiro.test(resid(model))$p
  85. mod_variables <- extract_full_variables(model)
  86. } else {
  87. model <- fit_child_model(dataframe, metric)
  88. if(!isSingular(model)){
  89. form = "no_exp"
  90. sw=shapiro.test(resid(model))$p
  91. mod_variables <- extract_chi_variables(model)
  92. } else {
  93. form='no_chi_effect'
  94. sw = NA
  95. mod_variables = c(
  96. NA,NA,NA, # c(coefficients(summary(model))["age_s",]
  97. NA,NA, # icc.result.mixed
  98. NA,NA, # icc.result.split
  99. NA,NA,NA, # ranefs_vars
  100. NA,NA,NA, # raners_std
  101. NA,NA, # r2
  102. NA,NA,NA) # nobs (child, corpus), nobs
  103. }
  104. }
  105. icc.row = c(data_set, age, metric, iqr, mod_variables, form, sw)
  106. return (icc.row)
  107. }