#Installing Packages

# install.packages("readxl")
# install.packages("tidyverse")
# install.packages("janitor")
# install.packages("magrittr")
# install.packages("dplyr")
# install.packages("tidyr")
# install.packages("ggplot2")
# install.packages("ggsignif")

#This chunk exists to download all the packages needed to clean the data from the excel file and develop a figure. It does not need to be re-downloaded, unlike calling it from the library. Hashtags were added after installation. #Loading Packages

library("readxl")
library("tidyverse")
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.0     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library("janitor")
## 
## Attaching package: 'janitor'
## 
## The following objects are masked from 'package:stats':
## 
##     chisq.test, fisher.test
library("magrittr")
## 
## Attaching package: 'magrittr'
## 
## The following object is masked from 'package:purrr':
## 
##     set_names
## 
## The following object is masked from 'package:tidyr':
## 
##     extract
library("dplyr")
library("tidyr")
library("ggplot2")
library("ggsignif")

#These packages will need to be loaded in every time I reopen the file to perform more work on the data set. Doing so allows me to actually use the installed packages through the written code. #Importing dataset GAPDH

GAPDH_data <- readxl::read_excel("GAPDH_Amygdala.xls", sheet = "Results", range = cell_rows(42:106))
## New names:
## • `` -> `...1`
## • `` -> `...2`
## • `` -> `...3`
## • `` -> `...4`
## • `` -> `...5`
## • `` -> `...6`
## • `` -> `...7`
## • `` -> `...8`
## • `` -> `...9`
## • `` -> `...10`
## • `` -> `...11`
## • `` -> `...12`
## • `` -> `...13`
## • `` -> `...14`
## • `` -> `...15`
## • `` -> `...16`
## • `` -> `...17`
## • `` -> `...18`
## • `` -> `...19`
## • `` -> `...20`
## • `` -> `...21`
## • `` -> `...22`
## • `` -> `...23`
## • `` -> `...24`
## • `` -> `...25`
## • `` -> `...26`
## • `` -> `...27`
## • `` -> `...28`
## • `` -> `...29`
## • `` -> `...30`
## • `` -> `...31`
## • `` -> `...32`

#I was able to use the shortened code because both the R file and the excel file were in the same place (desktop). If it imports correctly, GAPDH should appear in the environment section of R. If you do not read in the excel file, none of the code will work. If there is an error in reading, change the 109 to 106. It did need to be changed to 106 due to a mismatch in number of rows between the file and the cleaning code. Before moving on to making the graph, the data must be cleaned first. #Cleaning and Organizing for Analysis GAPDH

GAPDH_data <- GAPDH_data %>% row_to_names(row_number = 1)%>% clean_names()

GAPDH_data_list <- data.frame(colnames(GAPDH_data)) %>% rowid_to_column(., "index")

GAPDH_data <- GAPDH_data[,c(1,2,4,5,15,16,17)] 

GAPDH_data <- GAPDH_data[!GAPDH_data$ct == "Undetermined", ]

GAPDH_data <- GAPDH_data %>% separate(sample_name, c('label','sample_num'))
GAPDH_data$sample_num <- as.numeric(GAPDH_data$sample_num)
GAPDH_data <- GAPDH_data[order(GAPDH_data$sample_num, decreasing = FALSE), ]

GAPDH_data$treatment <- c('M','M','M', 'S','S','S', 'M','M','M', 'S','S','S', 'S','S','S', 'M','M','M', 'M','M','M', 'M','M','M', 'S','S','S', 'S','S','S', 'M','M','M', 'M','M','M', 'M','M','M', 'S','S','S', 'S','S','S', 'S','S','S', 'S','S','S','S','S','S', 'M','M','M', 'M','M','M')

GAPDH_data <- select(GAPDH_data, "well", "well_position", "treatment", "sample_num", "target_name", "ct", "ct_mean", "ct_sd")

#This code is responsible for cleaning and organising data for creating a graph. Specifically, lines 43, 49, and 55 are creating and naming different columns that will exist within the dataset generated in the environment. Without this, the graph may contain null values or other errors that impact results. A dataset with 32 objects should appear in the environment after this step. #Removing Repeats from GAPDH data

GAPDH_data[ , 'ct_mean_ind'] = NA

for (row in 1:nrow(GAPDH_data)) {
        if(isTRUE(GAPDH_data$sample_num[row] == GAPDH_data$sample_num[row+1])==TRUE) {
           GAPDH_data$ct_mean_ind[row] = NA
        } else {
            GAPDH_data$ct_mean_ind[row] = GAPDH_data$ct_mean[row]
          }
        }

#This is a futher cleaning step that removes an repeating values in the dataset. Leaving in repeating values could drastically impact futher calculations and the eventual graph by introducing unneeded numerical data. #Data Frame Check

GAPDH_data
## # A tibble: 60 × 9
##    well  well_position treatment sample_num target_name ct         ct_mean ct_sd
##    <chr> <chr>         <chr>          <dbl> <chr>       <chr>      <chr>   <chr>
##  1 1     A1            M                  1 GAPDH       18.195129… 18.235… 0.03…
##  2 2     A2            M                  1 GAPDH       18.242349… 18.235… 0.03…
##  3 3     A3            M                  1 GAPDH       18.269241… 18.235… 0.03…
##  4 13    B1            S                  2 GAPDH       18.210189… 18.218… 0.01…
##  5 14    B2            S                  2 GAPDH       18.233800… 18.218… 0.01…
##  6 15    B3            S                  2 GAPDH       18.210977… 18.218… 0.01…
##  7 25    C1            M                  3 GAPDH       18.036970… 18.028… 0.01…
##  8 26    C2            M                  3 GAPDH       18.033191… 18.028… 0.01…
##  9 27    C3            M                  3 GAPDH       18.014335… 18.028… 0.01…
## 10 37    D1            S                  4 GAPDH       18.231857… 18.230… 0.01…
## # ℹ 50 more rows
## # ℹ 1 more variable: ct_mean_ind <chr>

#Running this code allows me to check the data set before moving on. It also allowed me to double check that the removing repeats step actually created a new column in the dataset since it ran with some plus signs in it. I verified that it indeed created a new column. #Calculating Values for Analysis involving CT

vect <- GAPDH_data$ct_mean_ind
vect <- as.numeric(vect)
vect <- na.omit(vect) 
pop_mean <- mean(vect)

GAPDH_data$ct_mean_ind <- as.numeric(GAPDH_data$ct_mean_ind)
GAPDH_data$dCT <- GAPDH_data$ct_mean_ind - pop_mean 

GAPDH_data$raised_dCT <- 2^-(GAPDH_data$dCT)

GAPDH_data$percent_change <- GAPDH_data$raised_dCT *100

#This step is one of the main goals of the overall analysis codes for GAPDH. First, the code develops the mean CT values of each sample of GAPDH. Next, the ΔΔCT is calculated, which normalizes the data and creates a baseline. After this, the code develops the fold change or relative expression of GAPDH between samples and controls using the ΔΔCT value dveloped in the previous step. Lastly, the code develops the percent change of GAPDH across the samples. #Running T-Test on Developed Values

GAPDH_Analysis <- t.test(raised_dCT ~ treatment, data = GAPDH_data)

#This part of the code runs a t-test comparing the percent change of the samples in the GAPDH dataset. Percent change was calculated in the previous chunk. Doing this allows for significance to be determined, which for GAPDH should be ideally not significant. This is because expression levels are relatively stable under normal conditions. Hvaing GAPDAH is needed for normalization steps when analyzing experiemntal data from BDNF/leptin trials.