R/reusable.R

In faridcher/futils: Farid utility functions

#' ggplot mosaic
#'
#' @param var1
#' @param var2
#'
#' @return
#' @export
#'
#' @examples
#' with(ch1_rsf, ggMMplot(factor(human_foot), factor(used)))
ggMMplot <- function(var1, var2){
  require(ggplot2)
  levVar1 <- length(levels(var1))
  levVar2 <- length(levels(var2))

  jointTable <- prop.table(table(var1, var2))
  plotData <- as.data.frame(jointTable)
  plotData$marginVar1 <- prop.table(table(var1))
  plotData$var2Height <- plotData$Freq / plotData$marginVar1
  plotData$var1Center <- c(0, cumsum(plotData$marginVar1)[1:levVar1 -1]) +
    plotData$marginVar1 / 2

  ggplot(plotData, aes(var1Center, var2Height)) +
    geom_bar(stat = "identity", aes(width = marginVar1, fill = var2), col = "Black") +
    geom_text(aes(label = as.character(var1), x = var1Center, y = 1.05)) +
    scale_fill_discrete(name="Used",labels=c("No","Yes"))+
    xlab("Human foot print index")+
    ylab("")+
    scale_x_continuous(breaks=NULL)+
    scale_y_continuous(breaks=NULL)+
    theme_bw()
}

#' Lag of a vector
#'
#' @param x
#' @param shift
#'
#' @return
#' @export
#'
#' @examples
#' Lag(1:5,-1)
#' Lag(1:5,1)
#' Hmsic:Lag(1:5) #same effect but with NA. we want to avoid the NAs
Lag <- function(x,shift=1)
{
  if(shift==0)
    x
  else
    c(tail(x,-shift),head(x,shift))
}
#' keep only numeric columns of dataframe
#'
#' This is for usage in pairs function
#'
#' @param df
#'
#' @return
#' @export
#'
#' @examples
numeric_col <- function(df)
{
  df %>%
    sapply(mode) %in%
    c("numeric","logical","integer","double") %>%
    df[,.]
}

#' reset graphical parameters
#'
#' @return
#' @export
#'
#' @examples
par_reset <- function(...){
  op <- structure(list(xlog = FALSE, ylog = FALSE, adj = 0.5, ann = TRUE,
                       ask = FALSE, bg = "transparent", bty = "o", cex = 1, cex.axis = 1,
                       cex.lab = 1, cex.main = 1.2, cex.sub = 1, col = "black",
                       col.axis = "black", col.lab = "black", col.main = "black",
                       col.sub = "black", crt = 0, err = 0L, family = "", fg = "black",
                       fig = c(0, 1, 0, 1), fin = c(6.99999895833333, 6.99999895833333
                       ), font = 1L, font.axis = 1L, font.lab = 1L, font.main = 2L,
                       font.sub = 1L, lab = c(5L, 5L, 7L), las = 0L, lend = "round",
                       lheight = 1, ljoin = "round", lmitre = 10, lty = "solid",
                       lwd = 1, mai = c(1.02, 0.82, 0.82, 0.42), mar = c(5.1, 4.1,
                                                                         4.1, 2.1), mex = 1, mfcol = c(1L, 1L), mfg = c(1L, 1L, 1L,
                                                                                                                        1L), mfrow = c(1L, 1L), mgp = c(3, 1, 0), mkh = 0.001, new = FALSE,
                       oma = c(0, 0, 0, 0), omd = c(0, 1, 0, 1), omi = c(0, 0, 0,
                                                                         0), pch = 1L, pin = c(5.75999895833333, 5.15999895833333),
                       plt = c(0.117142874574832, 0.939999991071427, 0.145714307397962,
                               0.882857125425167), ps = 12L, pty = "m", smo = 1, srt = 0,
                       tck = NA_real_, tcl = -0.5, usr = c(0.568, 1.432, 0.568,
                                                           1.432), xaxp = c(0.6, 1.4, 4), xaxs = "r", xaxt = "s", xpd = FALSE,
                       yaxp = c(0.6, 1.4, 4), yaxs = "r", yaxt = "s", ylbias = 0.2), .Names = c("xlog",
                                                                                                "ylog", "adj", "ann", "ask", "bg", "bty", "cex", "cex.axis",
                                                                                                "cex.lab", "cex.main", "cex.sub", "col", "col.axis", "col.lab",
                                                                                                "col.main", "col.sub", "crt", "err", "family", "fg", "fig", "fin",
                                                                                                "font", "font.axis", "font.lab", "font.main", "font.sub", "lab",
                                                                                                "las", "lend", "lheight", "ljoin", "lmitre", "lty", "lwd", "mai",
                                                                                                "mar", "mex", "mfcol", "mfg", "mfrow", "mgp", "mkh", "new", "oma",
                                                                                                "omd", "omi", "pch", "pin", "plt", "ps", "pty", "smo", "srt",
                                                                                                "tck", "tcl", "usr", "xaxp", "xaxs", "xaxt", "xpd", "yaxp", "yaxs",
                                                                                                "yaxt", "ylbias"))
  par(op)
  par(...)
}

#' @export
par_fav <- function(...){
  par(
    mar = c(3, 3, 1, 1),
    mgp = c(1.7, .75, 0),
    tck = -.01,
    cex.axis = 1.25,
    cex.lab = 1.25,
    bty='l'
  )
  par(...)
}

#' @export
png_inch <- function(...,
                     width = 7,
                     height = 7,
                     units = 'in',
                     res= 300,
                     bg = 'transparent',
                     antialias = 'cleartype'){
  li <- list(...,
             width = width,
             height = height,
             units = units,
             res = res,
             bg = bg,
             antialias = antialias
  )
  do.call(png, li)
}

#' @export
alpha_col <- function(col,alpha=1)
  (col2rgb(col %>% unclass)/256) %>%
  t %>%
  data.frame() %$%
  rgb(red,green,blue,alpha=alpha)

# improved list of objects
.ls.objects <- function (pos = 1, pattern, order.by,
                         decreasing=FALSE, head=FALSE, n=5) {
  napply <- function(names, fn) sapply(names, function(x)
    fn(get(x, pos = pos)))
  names <- ls(pos = pos, pattern = pattern)
  obj.class <- napply(names, function(x) as.character(class(x))[1])
  obj.mode <- napply(names, mode)
  obj.type <- ifelse(is.na(obj.class), obj.mode, obj.class)
  obj.prettysize <- napply(names, function(x) {
    capture.output(format(utils::object.size(x), units = "auto")) })
  obj.size <- napply(names, object.size)
  obj.dim <- t(napply(names, function(x)
    as.numeric(dim(x))[1:2]))
  vec <- is.na(obj.dim)[, 1] & (obj.type != "function")
  obj.dim[vec, 1] <- napply(names, length)[vec]
  out <- data.frame(obj.type, obj.size, obj.prettysize, obj.dim)
  names(out) <- c("Type", "Size", "PrettySize", "Rows", "Columns")
  if (!missing(order.by))
    out <- out[order(out[[order.by]], decreasing=decreasing), ]
  if (head)
    out <- head(out, n)
  out
}

#' @export
lsos <- function(..., n=10) {
  .ls.objects(..., order.by="Size", decreasing=TRUE, head=TRUE, n=n)
}

#' @export
panel.smooth_abc <- function (x, y, col = par("col"), bg = NA, pch = par("pch"),
                              cex = 1, col.smooth = "red", span = 2/3, iter = 3, ...)
{
  points(x, y, pch = pch, col = col, bg = bg, cex = cex)
  #ok <- is.finite(x) & is.finite(y)
  #if (any(ok))
  #  lines(stats::lowess(x[ok], y[ok], f = span, iter = iter),
  #        col = col.smooth, ...)

  mtext(letters[assign('i',
                       get('i',
                           env = .GlobalEnv) + 1,
                       env = .GlobalEnv)] %>%
          paste(")"), adj = 0)
  grid()
}

#' @export
rm_global_functions <- function()
  rm(list=names(Filter(is.function,
                       sapply(ls(.GlobalEnv),
                              get,envir=.GlobalEnv))),
     envir = .GlobalEnv)

#' @export
#'
rename <- function(df, value, column=names(df),add=F){
  x <- names(df)

  if (is.numeric(column) | is.logical(column))
    column <- x[column]  #Take numeric input by indexing

  if (add)
    names(df)[x %in% column] <- paste0(column,value)
  else
    names(df)[x %in% column] <- value              #What you're interested in
  return(df)
}

#' caterpillar
#'
#' @param Z a complex vector of locations ordered by time
#' @param do_plot whether to plot the result
#'
#' @return a matrix with a dimension of (n-2)*(n-2) of pseudo absences
#' (n-2)*n.pseudos
#' @export
#'
#' @examples
#' #simulate
#' X <- cumsum(arima.sim(n=3000, model=list(ar=.7)))
#' Y <- cumsum(arima.sim(n=3000, model=list(ar=.7)))
#' Z <- X + 1i*Y
#'
#' op <- par(mfrow=c(1,2))
#' znull<-caterpillar_elie(Z[1:10],T)
#' znull<-caterpillar(Z[1:10],T)
#' par(op)
#'
#' print(system.time(caterpillar_elie(Z)))
#' print(system.time(caterpillar(Z)))
caterpillar <- function(used_sfc , time = NA_real_, n.pseudos=10) {
  #if move data aren't regularly sampled must divide by difftime
  #caterpillar only works with projected for now
  #if(!inherits(used_sf,c("sfc")))
  #  stop("used_sfc must be of type sf")

  geom <- sf::st_geometry(used_sfc)
  xy <- do.call(rbind,geom)

  Z <- complex(re = xy[,1], imag=xy[,2])
  n <- length(Z)
  S <- Mod(diff(Z))
  Phi <- Arg(diff(Z))
  Theta <- diff(Phi)
  RelSteps <- complex(mod = S[-1], arg=Theta)

  Z_inner <- Z[-c(1,n)] #2:9
  Azimuths <- complex(mod = 1, arg=Phi[-(n-1)] ) # [-9]

  #Add the rotated steps to the last step
  z_null <-apply(outer(Azimuths,RelSteps),2,'+',Z_inner)
  #the first and n-1 sample of pseudo-absences
  z_null <- subset(z_null,select=sample(2:ncol(z_null), n.pseudos))

  # +1 because, stratum starts from 1 in used_df
  stratum <- rep(1:nrow(z_null), ncol(z_null)) + 1
  if (all(!is.na(time)))
    time <- rep(time[-c(1,length(time))],ncol(z_null))

  #as.vector is by col. use t to make it by row
  z_null <- as.vector(z_null)
  cater <- data.frame(stratum=stratum, x=Re(z_null), y=Im(z_null),time=time )

  return(sf::st_as_sf(cater,coords=c('x','y'), crs=st_crs(geom) ))
}


#' plot caterpillar
#'
#' calling this function on the result of the ssf_ready, plot the sampling of caterpillar i.e. the ssf sampling.
#' @param ssf_ready
#'
#' @return
#' @export
#'
#' @examples
plot_caterpillar <- function(geometry,stratum,used){
  geometry %>%
    subset(used==1) %>%
    plot(
      col=alpha_col(1,.5),
      pch=19,cex=1,type='o',lwd=2)

  xy <- do.call(rbind,args = geometry)
  data.frame(stratum, used, x=xy[,1],y=xy[,2]) %>%
    merge(.,.,by='stratum',suffixes=1:2) %>%
    subset(used1==1 & used2==0) %$%
    {
      palette(wesanderson::wes_palette('Zissou'))
      st_segment(x1,y1,x2,y2,st_crs(geometry)) %>%
        st_sf(stratum) %>%
        plot(main="",add=T,type='l',
             col=stratum)

      legend("topright",
             c("animal track",
               paste0("stratum",1:6),
               "..."),
             col=c('black',1:6,1), #seq_along(unique(stratum))
             lwd=c(2,
                   rep(1,6),1),
             bg = NA
      )
    }
  palette('default')
}

#' relative weight of regression coefficients
#'
#' @param fit
#' @param ...
#'
#' @return
#' @export
#'
#' @examples
#' states <- as.data.frame(state.x77[,c("Murder", "Population",
#' "Illiteracy", "Income", "Frost")])
#' zstates <- as.data.frame(scale(states))
#' zfit <- lm(Murder~Population + Income + Illiteracy + Frost, data=zstates)
#' coef(zfit)
#' relweights(zfit)
relweights <- function(fit,...){
  R <- cor(fit$model)
  nvar <- ncol(R)
  rxx <- R[2:nvar, 2:nvar]
  rxy <- R[2:nvar, 1]
  svd <- eigen(rxx)
  evec <- svd$vectors
  ev <- svd$values
  delta <- diag(sqrt(ev))
  lambda <- evec %*% delta %*% t(evec)
  lambdasq <- lambda ^ 2
  beta <- solve(lambda) %*% rxy
  rsquare <- colSums(beta ^ 2)
  rawwgt <- lambdasq %*% beta ^ 2
  import <- (rawwgt / rsquare) * 100
  import <- as.data.frame(import)
  row.names(import) <- names(fit$coefficients[2:nvar])
  names(import) <- "Weights"
  import <- import[order(import),1, drop=FALSE]
  dotchart(import$Weights, labels=row.names(import),
           xlab="% of R-Square", pch=19,
           main="Relative Importance of Predictor Variables",
           sub=paste("Total R-Square=", round(rsquare, digits=3)),
           ...)
  return(import)
}
#' extract rect x range
#'
#' @param x a vector of the x axis
#' @param where a logical vector indicating where to extract the x ranges
#' @param width the width of the length one bout i.e. a single T with F on its sides (..FTF...)
#'
#' @return
#' @export
#'
#' @examples
rect_xrange <- function(x, where,width=1){
  #mode2bout
  cumsum(c(1,diff(where) != 0) ) [where] %>% {
    data.frame(min_x=tapply(x[where],., min),
               max_x=tapply(x[where],., max))
  } -> x_range

  equal_idx <- x_range[,1] == x_range[,2]
  x_range[equal_idx, 1] <- x_range[equal_idx, 1] - width
  x_range[equal_idx, 2] <- x_range[equal_idx, 2] + width
  return(x_range)
}

#' interval intersection
#'
#' works with time and numeric vectors
#' @param ... timestamp vectors of type POSIXct and POSIXlt
#'
#' @return
#' @export
#'
#' @examples
interval_intersection <- function(...)
{
  times <- list(...)

  if (inherits(times[[1]],"POSIXlt"))
    times <- lapply(times,as.POSIXct)

  rn <- lapply(times, range)
  mat <- do.call(rbind,rn)
  rng <- c(max(mat[,1]),min(mat[,2]))

  if (inherits(times[[1]],"POSIXct"))
  {
    #they may have not this attribute
    tzo <- attr(times[[1]][1], "tzone") # works with both POSIXct and POSIXlt
    rng <- as.POSIXct(rng, origin='1970-01-01 00:00.00 UTC',tz = tzo)
  }
  return(rng)
}

#' animate tracks
#'
#' @param tracks a list of tracks
#' @param delay
#' @param skip_row
#' @param pdf_out
#'
#' @return
#' @export
#'
#' @examples
animate_move <- function(x1,
                         y1,
                         x2,
                         y2,
                         time,
                         name_vec = c("M1", "M2"),
                         delay = .2,
                         skip_row = 10) {
  xl <- range(union(x1,x2))
  yl <- range(union(y1,y2))
  for(i in seq(1, max(length(x1),length(x2)),skip_row)) {
    #make multi support (one or more than one)
    #add distance to this. make it sf aware
    #add basemap
    plot(NA,type='n',xlim=xl, ylim=yl , xlab="X",ylab="Y",asp=1)

    lines(x1[1:i], y1[1:i], col=alpha_col("blue",.5), cex=0.5, lwd=1.2)#,type='l', asp=1)
    points(x1[i], y1[i], pch=19, col=alpha_col("blue",.5),cex=1.2)

    lines(x2[1:i], y2[1:i], col=alpha_col("red",.5), cex=0.5, lwd=1.2)#,type='l', asp=1)
    points(x2[i], y2[i], pch=19, col=alpha_col("red",.5),cex=1.2)

    title(paste(time[i]))
    legend("topright",name_vec,col=c("blue","red"),lwd=1.2, pch=19)

    box()
    Sys.sleep(delay)
  }
}

#' rsf ready
#'
#' @param move_sfc object of type st_sfc (st_point)
#' @param n_sample number of samples
#' @param time is not used in RSF
#'
#' @description samples and adds the used field
#' @return
#' @export
#'
#' @examples
#' todo: make it work with temporal covariates
rsf_ready <- function(move_sfc, n_sample=length(move_sfc) * 2, time= NA_real_)
{
  used <- st_sf(geometry=move_sfc,used=rep(1,length(move_sfc))) #,time=time)

  samps <- st_sample(st_convex_hull(st_combine(used)),size=n_sample)
  avail <- st_sf(geometry=samps, used = rep(0,length(samps))) #, time=sample(time,length(samps),T))
  return(rbind(avail,used))
}

#' ssf ready
#' @description  adds used and stratum fields and calls caterpillar
#'
#' @return
#' @export
#'
#' @examples
ssf_ready <- function(move_sfc, time = NA_real_, sample_size = 10)
{
  #time ordering for caterpillar
  sf::st_sf(geometry=move_sfc, used=rep(1, length(move_sfc)),
            stratum = 1:length(move_sfc),
            time = time) -> used_sf

  cbind(caterpillar(move_sfc,time,sample_size),#returns sf
        used = rep(0, (nrow(used_sf)-2) * sample_size)) -> null_sf

  return(rbind(null_sf,used_sf))
}

#' report time interval
#'
#' @param time_dif
#' @param plot
#'
#' @return
#' @export
#'
#' @examples
#' report_time_interval(cheet$time,units = "hours")
report_time_interval <- function(time,units="hours",plot=T) {
  time_dif <- get_time_dif(time,units)

  tab <- table(time_dif)
  ints <- as.numeric(names(tab))
  #divide by the first non-zero interval
  ints <- ints/ints[which(ints!=0)[1]]
  total_obs_without_gaps <- sum( tab*ints ) + 1
  if(plot){
    barplot(tab,ylab="# Obs",xlab=paste0("Interval",units,collapse = ' '),
            col="grey",main="Gappy Observations")
    legend("topright",legend = paste("Num of regular obs =",total_obs_without_gaps))
  }
  return(total_obs_without_gaps)
}

#' Regularize Gappy Track
#' @description Gets the x,y,t of the track and fill in the gaps of the track. It is assumed that xyt is
#' regular (i.e. 4h interval) with gaps. uses the smallest interval for interpolation.
#'
#' @param x
#' @param y
#' @param datetime
#'
#' @return dataframe of x,y coordinates at regular interval
#' @export
#'
#' @examples
#'
regularize_gappy_track <- function(x,y,datetime,units="hours") {
  #assumption: time intervals are in hour.
  old_df <- data.frame(x,y,datetime)[order(datetime),]

  time_dif <- old_df$datetime %>%
    diff(units=units) %>%
    #as.numeric() %>%
    c(NA_integer_)

  tab <- table(time_dif)
  ints <- as.numeric(names(tab))
  #the smallest non-zero interval

  smallest <- min(ints[which(ints!=0)])
  time_dif <- time_dif/smallest

  #later: make more efficient
  #total_obs <- report_time_interval(time_dif,F)
  #xy <- matrix(0,nrow=total_obs ,ncol=2)

  new_df <- data.frame(x=double(),y=double()) #empty dataframe to fill

  for (i in 1:(nrow(old_df))) {
    if(is.na(time_dif[i])) #include the last observation
      new_df <- rbind(new_df,old_df[i,c('x','y')])
    else if (time_dif[i] == 0) #remove duplicate observations
      next
    else if (time_dif[i] > 1 )
    {
      xx <- old_df[i:(i+1),'x']
      yy <- old_df[i:(i+1),'y']
      tt = c(1,time_dif[i]+1)
      df <- data.frame(tt,xx,yy)

      mod <- lm(xx ~ tt,df)
      xs <- predict(mod,data.frame(tt=1:(time_dif[i]) ) )

      mod <- lm(yy ~ tt,df)
      ys <- predict(mod,data.frame(tt=1:(time_dif[i]) ) )

      new_df <- rbind(new_df,data.frame(x=xs,y=ys))
    }
    else #time_dif==1
      new_df <- rbind(new_df,old_df[i,c('x','y')])
  }
  #new_df$OBJECTID <- 1:nrow(new_df)
  new_df$time <- seq(old_df$datetime[1],
                     by = as.difftime(smallest,units = "hours"),
                     length.out = nrow(new_df))
  new_df
}

#' scan_track
#'
#' Plotting x-y, time-x, time-y scan of a track.  This function will take x, y, and time coordinates or a \code{track} class object
#'
#' @param time time (can be a \code{\link{POSIXt}})
#' @param x x Coordinate. x,y coordiantes an be two separate vectors OR a complex "x" OR a two-column matrix/date-frame.
#' @param y y coordinate.
#' @param col color vector t
#' @param alpha intensity of the color
#' @param cex character expansion of the points
#' @param ... options to be passed to plot functions
#'
#' @examples
#' scan_track(rnorm(100),rnorm(100),col=1:100)
#' @export
scan_track <- function(x, y, time=seq_along(x), main="", col=1, alpha=0.5, cex=0.5, ...)
{
  #TODO make generic with sf method
  layout(rbind(c(1,2), c(1,3)))
  op <- par(mar = c(0,4,0,0), oma = c(4,0,4,4), xpd = NA)
  #the resizing creates the problem however
  #plot(rnorm(100),rnorm(100),col=alpha_col(1,.1),type='o')  works fine even with resizing.
  #The problem may be from the layout or other params
  plot_track(x,y,col=alpha_col(col,alpha),...)
  plot(time,x, type="o", pch=19, col=alpha_col(col,alpha), xaxt="n", xlab="", cex=cex, ...)
  #lines(time,x)
  plot(time,y, type="o", pch=19, col=alpha_col(col,alpha), cex=cex, ...)
  #lines(time,y)
  par(op,mfrow=c(1,1))
  title(main,outer = F)
}

#' @export
plot_track <- function(x,...)
  UseMethod("plot_track")

#' @export
plot_track.sfc <- function(x, ...)
{
  xy <- do.call(rbind,x)
  plot_track(xy[,1],xy[,2],...)
}

#' plot track
#'
#' @param x
#' @param y
#' @param pch
#' @param col
#' @param cex
#' @param type
#' @param ...
#'
#' @return
#' @export
#'
#' @examples
#' plot_track(rnorm(100),rnorm(100))
plot_track.default <- function(x,y, pch=19, col =rgb(0,0,0,.2), cex = 0.5, type = "o", ...)
{
  plot(x, y, asp=1, type=type, pch=pch, cex=cex, col=col, ...)
  points(x[1], y[1], bg="green", pch=21)
  points(x[length(x)], y[length(x)], bg="red", pch=23)
  legend("topright",c("start","end"),pt.bg=c("green","red"),pch=c(21,23),bty='n',bg = "transparent")
}

#' aggr_track_by_bout
#'
#' @param x
#' @param y
#' @param time
#' @param move_bout
#' @param units
#'
#' @return
#' @export
#'
#' @examples
aggr_track_by_bout <- function(pt_sfc,time,mode_null, units = "hour"){
  #move %>% subset(ID='ch1') %$% cumsum(c(1,diff(mode_null %>% unclass) != 0 )) %>% unique
  #become 56. why?

  bout <- cumsum(c(1,diff(mode_null %>% unclass) != 0 ))

  start_time <- aggregate(time, list(bout), min)[,-1] #exclude the grouping. tapply %>% no simplify %>% do.call works
  end_time <- aggregate(time, list(bout), max)[,-1]
  duration <- difftime(end_time, start_time, units = units)  %>%
    as_units()

  #spatial
  geometry <- tapply(pt_sfc,
                     bout,
                     st_pt_line,
                     simplify = F) %>%
    do.call(c,.)

  displacement <- tapply(pt_sfc,
                         bout,
                         st_pt_line,
                         simplify = F) %>%
    lapply(st_length) %>%  #sapply removes units
    do.call(c,.) %>%
    set_units('km')

  beeline <- tapply(pt_sfc,
                    bout,
                    st_beeline,
                    simplify = F) %>%
    lapply(st_length) %>%  #sapply removes units
    do.call(c,.) %>%
    set_units('km')

  daily_displacement <- displacement/set_units(duration ,'day')

  #mode
  mode <- tapply(mode_null %>% unclass, list(bout), min) %>%
    factor(labels=c("Encamped","Moving"))

  data.frame(
    bout_id=seq_along(start_time),
    start_time,
    end_time,
    duration,
    displacement,
    beeline,
    daily_displacement,
    mode,
    geometry
  ) %>%
    st_sf
}

#' aggr_cov_by_bout
#'
#' @param covariates
#' @param move_bout
#'
#' @return
#' @export
#'
#' @examples
aggr_cov_by_bout <- function(covariates,mode_null){
  move_bout <- cumsum(c(1,diff(mode_null %>% unclass) != 0 ))

  mean_df <- aggregate(covariates, list(move_bout) , mean)[,-1] %>%
    `names<-`(sprintf("%s_mean",names(.)))

  sd_df <- aggregate(covariates, list(move_bout), sd)[,-1] %>%
    `names<-`(sprintf("%s_sd",names(.)))

  data.frame(mean_df,sd_df)
}