ggircal: Raw Accelerometer Data Analysis

Documented in g.report.part4

g.report.part4 = function(datadir=c(),metadatadir=c(),loglocation = c(),f0=c(),f1=c(),storefolderstructure=TRUE) {
  # description: function to merge generate report from milestone data generated by g.part4 (if store.ms was set to TRUE)
  ms4.out = "/meta/ms4.out"
  if (file.exists(paste(metadatadir,ms4.out,sep=""))) {
    if (length(dir(paste(metadatadir,ms4.out,sep=""))) == 0) {
      try.generate.report = FALSE
    } else {
      try.generate.report = TRUE
    }
  } else {
    try.generate.report = FALSE
  }
  if (try.generate.report == TRUE) {
    resultfolder = metadatadir #resultfolder = "Q:/studies/sleep/output_pi_sleep_wrist"
    meta.sleep.folder = paste(metadatadir,"/meta/ms3.out",sep="")
    #   only.use.sleeplog = TRUE #do you want to fall back on generic assumptions about sleep if sleep log is not available?
    if (length(loglocation) > 0) {
      only.use.sleeplog = TRUE
    } else {
      only.use.sleeplog = FALSE
    }
    #============================================================================
    # Go through accelerometer datafiles and merge with sleep log data
    fnames = dir(meta.sleep.folder)
    if (f1 > length(fnames)) {
      print(paste("f1 changed from, ",f1," to ",length(fnames),sep=""))
      f1 = length(fnames)
    }
    if (length(f1) == 0 | f1 > length(fnames))  f1 = length(fnames)
    if (storefolderstructure == TRUE) {
      filelist = FALSE
      if (length(datadir) == 1) { #could be a directory or one file
        if (length(unlist(strsplit(datadir,"[.]bi")))>1) filelist = TRUE
        if (length(unlist(strsplit(datadir,"[.]cs")))>1) filelist = TRUE
      } else { #multiple files
        filelist = TRUE    
      }
      if (filelist == FALSE) {
        fnamesfull = c(dir(datadir,recursive=TRUE,pattern="[.]csv"),dir(datadir,recursive=TRUE,pattern="[.]bin"))
      } else {
        fnamesfull = datadir
      }
      f16 = function(X) {
        out = unlist(strsplit(X,"/"))
        f16 = out[length(out)]
      }
      f17 = function(X) {
        out = unlist(strsplit(X,"/"))
        f17 = out[(length(out)-1)]
      }
      ffd = ffp = rep("",length(fnamesfull))
      if (length(fnamesfull) > 0) {
        fnamesshort = apply(X=as.matrix(fnamesfull),MARGIN=1,FUN=f16)
        foldername = apply(X=as.matrix(fnamesfull),MARGIN=1,FUN=f17)
        for (i in 1:length(fnames)) { #
          ff = as.character(unlist(strsplit(fnames[i],".RDa"))[1])
          if (length(which(fnamesshort == ff)) > 0) {
            ffd[i] = fnamesfull[which(fnamesshort == ff)]
            ffp[i] = foldername[which(fnamesshort == ff)]
          }
        }
      }
    }
    #-----------------------------------------------------
    nightsummary2 = as.data.frame(matrix(0,0,29)) 
    colnames(nightsummary2) = c("id", "night","acc_onset", "acc_wake", "acc_timeinbed", "acc_def", 
                                "sleeplog_onset", "sleeplog_wake", "sleeplog_timeinbed",
                                "error_onset", "error_wake", "error_dur",
                                "fraction_night_invalid",
                                "acc_dur_noc","acc_dur_sibd","acc_n_noc","acc_n_sibd",
                                "acc_onset_ts","acc_wake_ts","sleeplog_onset_ts", "sleeplog_wake_ts",
                                "page","daysleeper","weekday","calendardate","filename",
                                "cleaningcode","sleeplog_used","acc_available")
    sumi = 1
    sleeplog_used = rep(" ",((f1-f0)+1))
    fnames.ms4 = list.files(paste0(metadatadir,ms4.out),full.names=TRUE)
    if(length(fnames.ms4) < f1) f1 = length(fnames.ms4)
    cat(" loading all the milestone data from part 4 this can take a few minutes\n")
    myfun = function(x) { 
      load(file=x)
      cut = which(nightsummary[,1] == "")
      if (length(cut) > 0 & length(cut) < nrow(nightsummary)) {
        nightsummary = nightsummary[-cut,]
      }
      out = as.matrix(nightsummary)
    }
    
    nightsummary2 = as.data.frame(do.call(rbind,lapply(fnames.ms4,myfun)),stringsAsFactors=FALSE)
    #----------------
    nightsummary = nightsummary2
    # next 4 lines moved here on 25/11/2015
    pko = which(nightsummary$acc_onset == 0 & nightsummary$acc_wake == 0 & nightsummary$acc_timeinbed == 0)
    if (length(pko) > 0) {
      nightsummary = nightsummary[-pko,]
    }
    #####################################################
    # COLLAPSING nightsummary TO A ONELINE summary PER PARTICIPANT
    if (nrow(nightsummary) == 0) {
      print("report not stored, because no results available")
    } else {
      write.csv(nightsummary,file=paste(resultfolder,"/results/QC/part4_nightsummary_sleep_full.csv",sep=""),row.names=FALSE)
      nightsummary_bu = nightsummary
    }
    ####
    summarynames_backup = c()
    for (dotwice in 1:2) { #store data twice, once full and once cleaned
      if (dotwice == 2) {
        # ignore nights that were derived without sleep log?
        if (only.use.sleeplog == TRUE) {
          del = which(nightsummary$cleaningcode > 0 | nightsummary$sleeplog_used == "FALSE")
        } else {
          #only delete nights with no or no valid accelerometer data data, but consider nigths with missing sleep log data
          del = which(nightsummary$cleaningcode > 1) 
        }
        if (length(del) > 0) {
          nightsummary = nightsummary_bu[-del,]
        }
      }
      NIDS = length(unique(nightsummary$id))
      NDEF = length(unique(nightsummary$acc_def))
      uuu = unique(nightsummary$acc_def)
      
      rem = which(uuu == 0 | uuu == "0" | is.na(uuu) == TRUE)
      if (length(rem) > 0) {
        uuu = uuu[-rem]
        NDEF = length(uuu)
      }
      if (storefolderstructure==TRUE) {
        summary = matrix(0,NIDS,((NDEF*3*19)+13+(6*3)))
      } else {
        summary = matrix(0,NIDS,((NDEF*3*19)+11+(6*3)))
      }
      uid = unique(nightsummary$id)
      
      if (nrow(nightsummary) > 0) {
        for (i in 1:length(uid)) {
          summarynames = c() #moved here on 3/12/2014
          # fully cleaned from nights that need to be deleted
          # nightsummary.tmp = nightsummary[which(nightsummary$id == uid[i] & nightsummary$cleaningcode == 0),]
          nightsummary.tmp = nightsummary[which(nightsummary$id == uid[i]),] #back up
          udef = as.character(unique(nightsummary.tmp$acc_def))
          if(length(which(as.character(udef)=="0") > 0)) udef = udef[-c(which(as.character(udef)=="0"))]
          udefn = udef
          #-------------------------------------------
          # general info about file
          summary[i,1] = uid[i] #id
          summarynames = c(summarynames,"id")
          summary[i,2] = as.character(nightsummary$filename[which(nightsummary$id == uid[i])][1]) #filename
          if (length(unlist(strsplit(as.character(summary[i,2]),".RDa"))) > 1) summary[i,2] = unlist(strsplit(summary[i,2],".RDa"))[1]
          summarynames = c(summarynames,"filename")
          cntt = 2
          summary[i,cntt+1] = as.character(nightsummary$calendardate[which(nightsummary$id == uid[i])][1]) #date
          summarynames = c(summarynames,"calendardate")
          summary[i,cntt+2] = nightsummary$weekday[which(nightsummary$id == uid[i])][1] #date
          summarynames = c(summarynames,"weekday")
          # sleep log used
          summary[i,cntt+3] = as.character(nightsummary.tmp$sleeplog_used[1])
          summarynames = c(summarynames,paste("sleeplog_used",sep=""))
          # total number of nights with acceleration and accelerometer worn
          summary[i,cntt+4] = length(which(
            nightsummary.tmp$acc_available[which(nightsummary.tmp$acc_def == udef[1])] == "TRUE" &
              nightsummary.tmp$cleaningcode[which(nightsummary.tmp$acc_def == udef[1])] != 2))
          summarynames = c(summarynames,paste("n_nights_acc",sep=""))
          # total number of nights with sleep log
          TOTALNIGHTS = length(nightsummary.tmp$night[which(nightsummary.tmp$acc_def == udef[1] & nightsummary.tmp$cleaningcode != 1 )]) # number of total nights
          summary[i,cntt+5] = TOTALNIGHTS
          summarynames = c(summarynames,paste("n_nights_sleeplog",sep=""))
          # total number of complete weekend and week nights
          th3 = nightsummary.tmp$weekday[which(nightsummary.tmp$acc_def == udef[1])]
          summary[i,cntt+6] = length(which(nightsummary.tmp$cleaningcode[which(nightsummary.tmp$acc_def == udef[1])] == 0 &
                                             (th3 == "Friday" | th3 == "Saturday")))
          summary[i,cntt+7] = length(which(nightsummary.tmp$cleaningcode[which(nightsummary.tmp$acc_def == udef[1])] == 0 &
                                             (th3 == "Monday" | th3 == "Tuesday" | th3 == "Wednesday" | th3 == "Thursday" | th3 == "Sunday")))
          summarynames = c(summarynames,paste("n_WE_nights_complete",sep=""),
                           paste("n_WD_nights_complete",sep=""))
          # number of days with sleep during the day
          summary[i,cntt+8] = length(which(nightsummary.tmp$daysleep[which(nightsummary.tmp$acc_def == udef[1])] == 1 &
                                             (th3 == "Friday" | th3 == "Saturday")))
          summary[i,cntt+9] = length(which(nightsummary.tmp$daysleep[which(nightsummary.tmp$acc_def == udef[1])] == 1 &
                                             (th3 == "Monday" | th3 == "Tuesday" | th3 == "Wednesday" | th3 == "Thursday" | th3 == "Sunday")))
          summarynames = c(summarynames,paste("n_WEnights_daysleeper",sep=""),
                           paste("n_WDnights_daysleeper",sep=""))
          cnt = cntt+9
          #-------------------------------------------
          # sleep log summary
          nightsummary.tmp$sleeplog_timeinbed = as.numeric(nightsummary.tmp$sleeplog_timeinbed)
          nightsummary.tmp$sleeplog_onset = as.numeric(nightsummary.tmp$sleeplog_onset)
          nightsummary.tmp$sleeplog_wake = as.numeric(nightsummary.tmp$sleeplog_wake)
          nightsummary.tmp$acc_timeinbed = as.numeric(nightsummary.tmp$acc_timeinbed)
          nightsummary.tmp$acc_onset = as.numeric(nightsummary.tmp$acc_onset)
          nightsummary.tmp$acc_wake = as.numeric(nightsummary.tmp$acc_wake)
                                                  
          nightsummary.tmp$acc_dur_noc = as.numeric(nightsummary.tmp$acc_dur_noc)
          nightsummary.tmp$acc_n_noc = as.numeric(nightsummary.tmp$acc_n_noc)
          nightsummary.tmp$acc_dur_sibd = as.numeric(nightsummary.tmp$acc_dur_sibd)
          nightsummary.tmp$acc_n_sibd = as.numeric(nightsummary.tmp$acc_n_sibd)
          weekday = nightsummary.tmp$weekday[which(nightsummary.tmp$acc_def == udef[1])]
          for (k in 1:3) {
            if (k == 1) {
              TW = "AD"
              Seli = 1:length(weekday)
            } else if (k == 2) {
              TW = "WD"
              Seli = which(weekday == "Monday" | weekday == "Tuesday" | weekday == "Wednesday" | weekday == "Thursday" | weekday == "Sunday")
            } else if (k == 3) {
              TW = "WE"
              Seli = which(weekday == "Friday" | weekday == "Saturday")              
            }
            summary[i,cnt+1] = mean(nightsummary.tmp$sleeplog_timeinbed[which(nightsummary.tmp$acc_def == udef[1])[Seli]],na.rm=TRUE)
            summary[i,cnt+2] = sd(nightsummary.tmp$sleeplog_timeinbed[which(nightsummary.tmp$acc_def == udef[1])[Seli]],na.rm=TRUE)
            summarynames = c(summarynames,paste("sleeplog_timeinbed_",TW,"_mn",sep=""),
                             paste("sleeplog_timeinbed_",TW,"_sd",sep=""))
            summary[i,cnt+3] = mean(nightsummary.tmp$sleeplog_onset[which(nightsummary.tmp$acc_def == udef[1])[Seli]],na.rm=TRUE)
            summary[i,cnt+4] = sd(nightsummary.tmp$sleeplog_onset[which(nightsummary.tmp$acc_def == udef[1])[Seli]],na.rm=TRUE)
            summarynames = c(summarynames,paste("sleeplog_onset_",TW,"_mn",sep=""),
                             paste("sleeplog_onset_",TW,"_sd",sep=""))
            summary[i,cnt+5] = mean(nightsummary.tmp$sleeplog_wake[which(nightsummary.tmp$acc_def == udef[1])[Seli]],na.rm=TRUE)
            summary[i,cnt+6] = sd(nightsummary.tmp$sleeplog_wake[which(nightsummary.tmp$acc_def == udef[1])[Seli]],na.rm=TRUE)
            summarynames = c(summarynames,paste("sleeplog_wake_",TW,"_mn",sep=""),
                             paste("sleeplog_wake_",TW,"_sd",sep=""))
            cnt = cnt + 6
          }
          nightsummary$cleaningcode = as.numeric(nightsummary$cleaningcode)
          nightsummary$id = as.character(nightsummary$id)
          uid = as.character(uid)
          #-------------------------------------------
          # accelerometer summary
          #----------------------------------------------
          if (only.use.sleeplog == FALSE) { #when sleep log is not available
            if (dotwice == 2) {
              CRIT = which(nightsummary$id == uid[i] & (nightsummary$cleaningcode == 0 | nightsummary$cleaningcode == 1))
            } else {
              CRIT = which(nightsummary$id == uid[i])
            }
          } else {
            CRIT = which(nightsummary$id == uid[i] & nightsummary$cleaningcode == 0) #when sleep log is available
          }
          # summarynames_backup = c()
          if (length(CRIT) > 0) { #summarise data if there is data
            #-----------------------------------------------
            #         if (length(which(nightsummary$id == uid[i] & nightsummary$cleaningcode == 0)) > 0) { #summarise data if there is data
            for (j in 1:length(udef)) {
              weekday = nightsummary.tmp$weekday[which(nightsummary.tmp$acc_def == udef[j])]
              for (k in 1:3) {
                if (ncol(summary) < (cnt + 22)) { # expand summary matrix if there is a change that is not big enough
                  expansion = as.matrix(summary[,(cnt+1):ncol(summary)])
                  if (nrow(expansion) != nrow(summary)) expansion = t(expansion)
                  summary = cbind(summary,expansion)
                }
                
                if (k == 1) {
                  TW = "AD"
                  Seli = 1:length(weekday)
                } else if (k == 2) {
                  TW = "WD"
                  Seli = which(weekday == "Monday" | weekday == "Tuesday" | weekday == "Wednesday" | weekday == "Thursday" | weekday == "Sunday")
                } else if (k == 3) {
                  TW = "WE"
                  Seli = which(weekday == "Friday" | weekday == "Saturday")              
                }
                indexUdef = which(nightsummary.tmp$acc_def == udef[j])[Seli]
                summary[i,(cnt+1)] = mean(nightsummary.tmp$acc_timeinbed[indexUdef],na.rm=TRUE)
                summary[i,(cnt+2)] = sd(nightsummary.tmp$acc_timeinbed[indexUdef],na.rm=TRUE)
                summarynames = c(summarynames,paste("acc_timeinbed_",TW,"_",udefn[j],"_mn",sep=""),
                                 paste("acc_timeinbed_",TW,"_",udefn[j],"_sd",sep=""))
                summary[i,(cnt+3)] = mean(nightsummary.tmp$acc_dur_noc[indexUdef],na.rm=TRUE)
                summary[i,(cnt+4)] = sd(nightsummary.tmp$acc_dur_noc[indexUdef],na.rm=TRUE)
                summarynames = c(summarynames,paste("acc_dur_noc_",TW,"_",udefn[j],"_mn",sep=""),
                                 paste("acc_dur_noc_",TW,"_",udefn[j],"_sd",sep=""))
                summary[i,(cnt+5)] = mean(nightsummary.tmp$acc_dur_noc[indexUdef] /
                                            nightsummary.tmp$acc_timeinbed[indexUdef],na.rm=TRUE)
                summary[i,(cnt+6)] = sd(nightsummary.tmp$acc_dur_noc[indexUdef] /
                                          nightsummary.tmp$acc_timeinbed[indexUdef],na.rm=TRUE)
                summarynames = c(summarynames,paste("acc_eff_",TW,"_",udefn[j],"_mn",sep=""),
                                 paste("acc_eff_",TW,"_",udefn[j],"_sd",sep=""))
                summary[i,(cnt+7)] = mean(nightsummary.tmp$acc_dur_sibd[indexUdef],na.rm=TRUE)
                summary[i,(cnt+8)] = sd(nightsummary.tmp$acc_dur_sibd[indexUdef],na.rm=TRUE)
                summarynames = c(summarynames,paste("acc_dur_sibd_",TW,"_",udefn[j],"_mn",sep=""),
                                 paste("acc_dur_sibd_",TW,"_",udefn[j],"_sd",sep=""))
                summary[i,(cnt+9)] = mean(nightsummary.tmp$acc_n_noc[indexUdef],na.rm=TRUE)
                summary[i,(cnt+10)] = sd(nightsummary.tmp$acc_n_noc[indexUdef],na.rm=TRUE)
                summarynames = c(summarynames,paste("acc_n_noc_",TW,"_",udefn[j],"_mn",sep=""),
                                 paste("acc_n_noc_",TW,"_",udefn[j],"_sd",sep=""))
                summary[i,(cnt+11)] = mean(nightsummary.tmp$acc_n_sibd[indexUdef],na.rm=TRUE)
                summary[i,(cnt+12)] = sd(nightsummary.tmp$acc_n_sibd[indexUdef],na.rm=TRUE)
                summarynames = c(summarynames,paste("acc_n_sibd_",TW,"_",udefn[j],"_mn",sep=""),
                                 paste("acc_n_sibd_",TW,"_",udefn[j],"_sd",sep=""))
                # average sibd per night
                AVEsibdDUR = c(nightsummary.tmp$acc_dur_sibd[indexUdef] /
                                 nightsummary.tmp$acc_n_sibd[indexUdef])
                if (length(which(nightsummary.tmp$acc_n_sibd[indexUdef] == 0))) {
                  AVEsibdDUR[which(nightsummary.tmp$acc_n_sibd[indexUdef] == 0)] = 0
                }
                summary[i,(cnt+13)] = mean(AVEsibdDUR,na.rm=TRUE)
                summary[i,(cnt+14)] = sd(AVEsibdDUR,na.rm=TRUE)
                summarynames = c(summarynames,paste("acc_dur_msibd_",TW,"_",udefn[j],"_mn",sep=""),
                                 paste("acc_dur_msibd_",TW,"_",udefn[j],"_sd",sep=""))
                NDAYsibd = length(which(nightsummary.tmp$acc_n_sibd[indexUdef] > 0))
                if (length(NDAYsibd) == 0) NDAYsibd = 0
                summary[i,(cnt+15)] = NDAYsibd
                summarynames = c(summarynames,paste("acc_n_days_w_sibds_",TW,"_",udefn[j],sep=""))
                summary[i,(cnt+16)] = mean(nightsummary.tmp$acc_onset[indexUdef],na.rm=TRUE)
                summary[i,(cnt+17)] = sd(nightsummary.tmp$acc_onset[indexUdef],na.rm=TRUE)
                summarynames = c(summarynames,paste("acc_onset_",TW,"_",udefn[j],"_mn",sep=""),
                                 paste("acc_onset_",TW,"_",udefn[j],"_sd",sep=""))
                summary[i,(cnt+18)] = mean(nightsummary.tmp$acc_wake[indexUdef],na.rm=TRUE)
                summary[i,(cnt+19)] = sd(nightsummary.tmp$acc_wake[indexUdef],na.rm=TRUE)
                summarynames = c(summarynames,paste("acc_wake_",TW,"_",udefn[j],"_mn",sep=""),
                                 paste("acc_wake_",TW,"_",udefn[j],"_sd",sep=""))
                cnt = cnt + 19
              }
            }
            summarynames_backup = summarynames #if (length(summarynames) >= 29)
          }
        }
        summary = as.data.frame(summary)
        if (length(summarynames) != ncol(summary)) {
          if (length(summarynames_backup) > 0) {
            names(summary) = summarynames_backup
          } else {
            if (length(summarynames) >ncol(summary)) {
              names(summary)[1:length(summarynames)] = summarynames
            } else {
               names(summary) = summarynames[1:ncol(summary)]
            }
            
          }
        } else {
          names(summary) =summarynames
        }
        
        if (storefolderstructure==TRUE) {
          colnames(summary)[length(colnames(summary))-1] = "filename_dir"
          colnames(summary)[length(colnames(summary))] = "foldername"
          summary$filename_dir = as.character(summary$filename_dir)
          summary$foldername = as.character(summary$foldername)
          for (iii in 1:length(fnames)) { 
            tyi = as.character(unlist(strsplit(as.character(ffd[iii]),"/")))
            tyi = tyi[length(tyi)]
            if (length(tyi) >0) {
              indexx = which(as.character(summary$filename) == tyi)
              if (length(indexx) > 0) {
                summary$filename_dir[indexx] = ffd[iii]
                summary$foldername[indexx] = ffp[iii]
                
              }
            }
          }
        }
      }
      #######################################################
      if (nrow(nightsummary) == 0) {
        print("report not stored, because no results available")
      } else {
        if (dotwice == 1) {
          write.csv(nightsummary,file=paste(resultfolder,"/results/QC/part4_nightsummary_sleep_full.csv",sep=""),row.names=FALSE)
          write.csv(summary,file=paste(resultfolder,"/results/QC/part4_summary_sleep_full.csv",sep=""),row.names=FALSE)
        } else {
          write.csv(nightsummary,file=paste(resultfolder,"/results/part4_nightsummary_sleep_cleaned.csv",sep=""),row.names=FALSE)
          write.csv(summary,file=paste(resultfolder,"/results/part4_summary_sleep_cleaned.csv",sep=""),row.names=FALSE)
        }
      }
    }  
    SI = sessionInfo() 
    save(SI,file=paste(resultfolder,"/results/QC/sessioninfo_part4.RData",sep=""))
  }
}