---
title: "Vertical Jump Data Processing Example"
author: "Colin Bond"
date: "3/22/2022"
output: html_document
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

New R users need to install the packages they require. Packages only need to be installed one time, unless a new version of the package is released and you wish to update to the new version. 

signal: Used to filter the ground reaction force data. https://cran.r-project.org/web/packages/signal/index.html
ggplot2: Used to plot results. https://cran.r-project.org/web/packages/ggplot2/index.html
All other functions used in this demonstration are available in Base R, which is what you get with installing R or RStudio.

```{r}
install.packages("signal") #Only run if you have not installed.
install.packages("ggplot2") #Only run if you have not installed.
```

For every R session, you need to "load" the packages that you want to use, assuming they are already installed. You only have to load packages one time per session. 
```{r}
library(signal) #Load package.
library(ggplot2) #Load package.
```

The next step is to import data. The below line creates a data frame called "data" by reading a csv file found at the given path with the first row in the csv file as a header row containing variable names. Note that for windows machines you need to use double "\\" between components in the path. The first 5 rows of data are then printed to visualize the structure of the data Fz is vertical ground reaction force and the variable of interest in this analysis. 

```{r}
data_files <- list.files("C:\\Users\\703354911\\Desktop\\NLACSM Presentation", pattern=".csv") #Identifies all files in the directory that end in .csv
data_files #Prints file names to the screen
all_data <- data.frame() #Creates an empty data frame that the vertical jump height data will be stored in.
for(i in 1:length(data_files)) {
data <- read.csv(paste0("C:\\Users\\703354911\\Desktop\\NLACSM Presentation\\",data_files[i]), header=TRUE) #Loads data from directory. Make sure to change directory to your own.
head(data) #Prints first 6 rows of data.
data$sample <- seq.int(nrow(data)) #Creates a sequential column from 1 to end called "sample".
Hz = 1/1000 #The sampling rate from the force plate used.
data$time <- data$sample*Hz #Creates a sequential column from 0.001 to end called "time".
bf <- butter(4, 0.1, type="low") #Sets up 4th order 10Hz lowpass Butterworth filter
data$Fz_filt <- filter(bf, data$Fz) #Applies filter.
print(BW <- mean(data[100:500,"Fz_filt"])) #BW is the mean of Fz_filt from 0.1 to 0.5 seconds.
data$Fz_filt_norm <- data$Fz_filt - BW #Subtracts BW from Fz_filt and creates new variable called Fz_filt_norm in data frame.
begin_negative_impulse = min(which(data[100:nrow(data),"Fz_filt_norm"] < -20))+100 #The beginning of the counter movement. Looks beginning at the 100th sample (0.1 s) and finds the first sample where Fz_filt_norm goes below -20 N. 20 or -20 N is commonly used threshold to distinguish an event, but this could be adjusted depending on the inherit noise in a given force plate.
begin_negative_impulse #Print sample to screen.
begin_positive_impulse = min(which(data[begin_negative_impulse:nrow(data),"Fz_filt_norm"] > 20))+begin_negative_impulse #The end of the counter movement. Looks beginning at the sample identified as begin_negative_impulse and finds the first sample where Fz_filt_norm goes above 20 N.
begin_positive_impulse #Print sample to screen.
end_positive_impulse = min(which(data[begin_positive_impulse:nrow(data),"Fz_filt_norm"] < 0))+begin_positive_impulse #The end of the positive impulse. Looks beginning at the sample identified as begin_positive_impulse and finds the first sample where Fz_filt_norm goes below 0 N.
end_positive_impulse #Print sample to screen.
Impulse = sum(data[begin_positive_impulse:end_positive_impulse,"Fz_filt_norm"]*Hz) #Calculates the impulse, or integral of force, from the beginning to end of the positive impulse.
Impulse #Prints impulse to screen.
Takeoff_Velocity = Impulse/(BW/9.81) #Calculates the take off velocity.
Takeoff_Velocity #Prints take off velocity to screen.
Jump_Height = 0.5 * (Takeoff_Velocity^2/9.81) #Calculated jump height in m.
Jump_Height #Prints jump height in m to screen.
Jump_Height_in = Jump_Height*39.37 #Converts jump height in m to jump height in inches.
Jump_Height_in #Prints jump height in inches to screen.
dat<- cbind(Jump_Height_in,data_files[i]) #Takes the vertical jump height in inches and the file name and puts them in two columns.
all_data <- rbind(all_data,dat) #Puts jump height in inches and data file name into our database that contains the results from all iterations of the loop.
}
colnames(all_data) = c("Jump_Height_in","Subject") #Renames the columns in our database.
all_data$Jump_Height_in <- as.numeric(all_data$Jump_Height_in) #Makes the vertical jump height a numeric variable.
all_data$Groups <- c("Female", "Female", "Female", "Male", "Male", "Male")
```

```{r}
plot1 <- ggplot(all_data, aes(x=Subject, y=Jump_Height_in)) + geom_point() + scale_y_continuous(limits=c(12,22)) #Plots our vertical jump heights for men and women.
print(plot1) #Prints plot to screen.
```

```{r}
plot2 <- ggplot(all_data, aes(x=Groups, y=Jump_Height_in)) + geom_point() + scale_y_continuous(limits=c(12,22)) + stat_summary(fun.data = "mean_cl_boot", colour = "red", size = .5) #Plots our vertical jump heights for men and women.
print(plot2) #Prints plot to screen
```

```{r}
compare <- t.test(Jump_Height_in ~ Groups, data = all_data, var.equal=TRUE)
compare
```


