---
title: "Entropy Analysis"
output: html_notebook
---

This is an [R Markdown](http://rmarkdown.rstudio.com) Notebook. When you execute code within the notebook, the results appear beneath the code. 

Try executing this chunk by clicking the *Run* button within the chunk or by placing your cursor inside it and pressing *Cmd+Shift+Enter*. 

```{r Load packages}
library(entropy)
library(ccber)
library(entropy)
library(dplyr)
library(lubridate)

```

###1. Simulated data
```{r Set parameters}

#actual minimum and maximum NDVI
ndvi.max <- 0.6
ndvi.min <- 0.25

#generate some values in between
range.ndvi <- sin(seq(0,2*pi,length.out = 100))
range.min <- min(range.ndvi)
range.max <- max(range.ndvi)
sim.ndvi <- ((range.ndvi - min(range.ndvi))*(ndvi.max-ndvi.min)/(range.max-range.min)) +ndvi.min

#5 Behaviours simulated: rest, forage, move, groom and autogroom 
possible_behaviours = c("R","F","L","G","U")

#Probability of behaviours at maximum NDVI based on previous work by Young et al. (2019)
NDVI_max = c(0.326,0.373,0.198,0.090,0.013)

sequence_size = 20 #set seuqence length according to focal length and chosen sampling frequency

number_of_focals = 60 #number of focals per inidividual to be simulated

number_of_individuals = 27 #number of individuals

```

```{r Simulate data}

#make a dataframe to store the data
sim.df.ndvi <- data.frame(ID=(1), focal_seq_max=-1, ndvi=-1)

#make a loop for individuals
for(j in 1:number_of_individuals){
  
  #Make a loop and generate a line of data each time
  for(i in 1:number_of_focals){
    
    #get a focal for each ndvi value
    for(n in sim.ndvi){
      
      #use NDVI values to adjust behaviour frequency
      NDVI.now <- NDVI_max
      NDVI.now[1] <- NDVI_max[1] + 0.05*n 
      NDVI.now[2] <- NDVI_max[2] - 0.05*n 
      
      focal_seq_max <-  sample(possible_behaviours, size = sequence_size, prob = NDVI.now, replace = TRUE)
      
      #add the sequence to the growing dataset
      sim.df.ndvi <- rbind(sim.df.ndvi, data.frame(ID=j,focal_seq_max=paste0(focal_seq_max, collapse = ''), ndvi = n))
    }
    
  }
  
}

#remove the first row of the dataframe (just used to initialize it with the right type of variable)
sim.df.ndvi<-sim.df.ndvi[-1,]

#subset for each id (60)
sim.ndvi.final <- sim.df.ndvi %>% group_by(ID) %>% sample_n(60)

sim.ndvi.final

```

```{r Calculate entropy rate}

#function to help calculate entropy for each row (CalcEntropyRate from ccber)

get.entropy <- function(x){
  x <- as.numeric(as.factor(unlist(strsplit(as.character(x[2]), split=NULL))))
  
  quicker_estimate <- CalcEntropyRate(x, 
                                      state_space = 1:max(x), 
                                      stat_method = "Empirical")
  return(quicker_estimate)
}


#calculate for each row in the dataframe
ent_rate<-apply(sim.ndvi.final, 1, get.entropy)

#put it back in the dataframe
sim.ndvi.final<-cbind(sim.ndvi.final,ent_rate)

#Calculate sequence length

get.length <- function(x){
  x <- as.numeric(as.factor(unlist(strsplit(as.character(x[2]), split=NULL))))
  
  length.seq <- length(x)
  return(length.seq)
}

#calculate for each row in the dataframe
seq.length<-apply(sim.ndvi.final, 1, get.length)

#put it back in the dataframe
sim.ndvi.final<-cbind(sim.ndvi.final,seq.length)

```
##2. Focal data
```{r Converting focal samples to sequences}

#Import data
df.focal <- read.csv("focal_data_example.csv")

#set the sampling frequency resolution: how often a behaviour is sampled
res = seconds("30")

#make the time column a real time column
df.focal$Time <- hms(df.focal$Time)
DateTime <- ymd_hms(paste0(as.Date(dmy_hms(df.focal$TimeStamp) ), " ", df.focal$Time))
df.focal <- cbind(df.focal, DateTime)

#make the date/time column a real date/time column
df.focal$date_time <- as_datetime(dmy_hm(df.focal$date_time))

#make sure data are in time order
df.focal <- df.focal %>% arrange(DateTime)

#create unique id for each focal: ID_Date  (this only works if an individual is never sampled twice durring a day...)
df.focal$focal_ID <- paste0(df.focal$Subject,"_", df.focal$date_time)

#get each focal
unique.focals <- unique(df.focal$focal_ID)

#run a loop to extract all sequences
focal.seq <- data.frame(Subject="Lucy",focal_ID="Lucy",sequence="TTTTTTTEEEESSSTTTTT", stringsAsFactors = F)

for(i in unique.focals){
  
  #get data from only that focal
  df.sub <- df.focal[df.focal$focal_ID==i,]
  
  #get all behavioural intervals
  df.sub$timeDate <- ymd_hms(paste0("2019-03-12"," ", df.sub$Time)) #the date here is just so we can calculate an interval
  df.sub <- df.sub %>% mutate(leadT = lead(DateTime))
  df.sub$intervals <- interval(df.sub$DateTime,df.sub$leadT-seconds(1)) 
  
  #sample from intervals using a moving time slot
  startTime=df.sub[1,]$DateTime 
  endTime=df.sub[nrow(df.sub),]$DateTime
  
  j=startTime
  sequence <- vector()
  while(j<endTime){
    
    #get which behaviour is sampled at that time
    sequence[length(sequence)+1]<-as.character(df.sub[(j %within% df.sub$intervals),]$Behaviour)
    
       #jump in time
    j<- j+res
    
  }
  
  #store the data for each focal 
  focal.seq <- rbind(focal.seq, data.frame(Subject =df.sub$Subject[1] , focal_ID=i,sequence=paste(sequence,collapse='') ))
}

#drop the first row (used to format the dataframe)
focal.seq<-focal.seq[-1,]

#could then bind it back to the original data frame if we'd like to get Subject and group
#df.groupID<-df.focal%>%group_by(Subject,date_time)%>%slice(1)%>%dplyr::select(Troop,focal_ID,Date)

#get every individuals group
#df.out <- left_join(focal.seq,df.focal, by="focal_ID")

#take a look
focal.seq

```

```{r Entropy rate of focals}

#function to help calculate entropy for each row
get.entropy <- function(x){
  x <- as.numeric(as.factor(unlist(strsplit(as.character(x[3]), split=NULL)))) #3=column with sequence
  
  quicker_estimate <- CalcEntropyRate(x, 
                                      state_space = 1:max(x), 
                                      stat_method = "Empirical")
  return(quicker_estimate)
}


#calculate for each row in the dataframe
ent_rate<-apply(df.out, 1, get.entropy)

#put it back in the dataframe
df.ent.rate<-cbind(df.out,ent_rate)

#sequence length#
get.length <- function(x){
  x <- as.numeric(as.factor(unlist(strsplit(as.character(x[3]), split=NULL))))
  
  length.seq <- length(x)
  return(length.seq)
}

#calculate for each row in the dataframe
seq.length<-apply(df.out, 1, get.length)

df.ent.rate<-cbind(seq.length,df.ent.rate)

```


Add a new chunk by clicking the *Insert Chunk* button on the toolbar or by pressing *Cmd+Option+I*.

When you save the notebook, an HTML file containing the code and output will be saved alongside it (click the *Preview* button or press *Cmd+Shift+K* to preview the HTML file). 

The preview shows you a rendered HTML copy of the contents of the editor. Consequently, unlike *Knit*, *Preview* does not run any R code chunks. Instead, the output of the chunk when it was last run in the editor is displayed.