library(openxlsx)
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──
## ✓ ggplot2 3.3.5 ✓ purrr 0.3.4
## ✓ tibble 3.1.6 ✓ dplyr 1.0.8
## ✓ tidyr 1.2.0 ✓ stringr 1.4.0
## ✓ readr 2.1.2 ✓ forcats 0.5.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(lubridate)
##
## Attaching package: 'lubridate'
## The following objects are masked from 'package:base':
##
## date, intersect, setdiff, union
library(dplyr)
library(tidyr)
library(ggplot2)
library("Rmisc")
## Loading required package: lattice
## Loading required package: plyr
## ------------------------------------------------------------------------------
## You have loaded plyr after dplyr - this is likely to cause problems.
## If you need functions from both plyr and dplyr, please load plyr first, then dplyr:
## library(plyr); library(dplyr)
## ------------------------------------------------------------------------------
##
## Attaching package: 'plyr'
## The following objects are masked from 'package:dplyr':
##
## arrange, count, desc, failwith, id, mutate, rename, summarise,
## summarize
## The following object is masked from 'package:purrr':
##
## compact
library("plyr")
library("parallel")
library(dabestr)
## Loading required package: magrittr
##
## Attaching package: 'magrittr'
## The following object is masked from 'package:purrr':
##
## set_names
## The following object is masked from 'package:tidyr':
##
## extract
library(RColorBrewer)
library(ggridges)
library(ggrepel)
t5<-Sys.time()
# The genotypes will be arranged by? Choose "default" or "customize".
analysis_method <- "default"
####################################### ↓↓↓↓ FUNCTION ↓↓↓↓ #####################################
# Define the funcion of extract the information of SEXUAL, GENO, expgroup and genptype.
extract_f <- function(x){
# Extract the information of SEXUAL and GENO
data_sex <- as.character( x[,'genotype']) %>% strsplit( .,'_') %>% as.data.frame(head=FALSE) %>%t() %>% .[,2]%>% as.character()
data_geno <- as.character( x[,'genotype']) %>% strsplit( .,'_') %>% as.data.frame(head=FALSE) %>%t() %>% .[,1]%>% as.character()
# Add the information of SEXUAL and GENO, differ the group of CTRL and EXP, filter the group of UAS
df <- x%>%
mutate(genotype.s=data_geno)%>%
mutate(sex=data_sex) %>%
mutate(group.c= ifelse(grepl( '/+', genotype.s),'ctrl', ifelse(grepl( '>', genotype.s),'exp',ifelse(grepl( 'UAS', genotype.s),'UAS','NA') ))) %>%
filter(group.c!="UAS")
# Then extract the information of genotype,
genotype_t <- as.character( df[,'genotype.s']) %>% strsplit( .,'>')%>% as.data.frame(head=FALSE)%>% t() %>% .[,1]%>% as.character()%>% strsplit( .,'\\-') %>% as.data.frame(head=FALSE)%>%t() %>% .[,1]%>% as.character()
df <- df%>%
mutate(genotype.t=genotype_t)
df
}
####################################### ↑↑↑↑ FUNCTION ↑↑↑↑ #####################################
####################################### ↓↓↓↓ FUNCTION ↓↓↓↓ #####################################
# Define the funcion of read csv in the folder.
read_d_f <- function(x){
df <- read.csv(paste0( getwd(),"/",x))
for (i in 1:files_number) {
df_path <- paste0(getwd(),"/DAM data/",file_name[i],"/",x)
df_raw <- read.csv(df_path,header = TRUE)
df <- rbind(df,df_raw)
}
df
}
####################################### ↑↑↑↑ FUNCTION ↑↑↑↑ #####################################
#get the list of files under /DAM data
file_name <- list.files(paste0( getwd(),"/DAM data/"))
files_number <- length(file_name)
############### Read the sleep data frame from "DAM data " folder
s_total_data <- read_d_f("sleep_dn_total.csv")
sbd_data <- read_d_f("sbd_period_bar.csv")
sb_data <- read_d_f("sleep_bout_period_bar_data.csv")
############### ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
s_data_t <- extract_f(s_total_data)
sbd_data_t <- extract_f(sbd_data)
sb_data_t <- extract_f(sb_data)
#Get the genotype list
f_list <- s_data_t$genotype.s %>% unlist() %>% unique()
# create auto_geno_list
f_list_m2 <- s_data_t$genotype.t %>% unlist() %>% unique()
l_2 <-s_data_t %>% dplyr::arrange(.,group.c)%>% dplyr::select(.,genotype.t,genotype.s) %>% unique()
l_2 <- split(l_2, l_2$genotype.t)
auto_geno_list <- list()
for (i in 1:length(l_2)) {
auto_geno_list[i] <- l_2[[i]][2]
}
genotypecsv <- data.frame(genotype= character())
for (n in 1:length(f_list_m2)) {
genotypecsv[n,'Genotype'] <- f_list_m2[n]
n=n+1
}
#create input_geno_list
input_geno_list <- read.csv(file = paste0(getwd(),"/genotypes.csv"),header = FALSE) %>% t() %>% as.data.frame() %>% as.list()
#choose the type of arrange AUTO or USER INPUT
print(paste0("Now the Analysis method is ", analysis_method, "."))
## [1] "Now the Analysis method is default."
write_geno_list <- auto_geno_list %>% as.data.frame() %>% t()
write.csv(write_geno_list, paste0( getwd(),"/genotypes.csv"), na = '',col.names=TRUE,row.names = FALSE,)
ifelse(analysis_method=="customize",geno_list <- input_geno_list, geno_list <- auto_geno_list)
## [[1]]
## [1] "104y-Gal4/+" "104y>DmGluR-RNAi34872"
#Output graph of sleep
g_d_sleep_total <- s_data_t %>% group_by(genotype, F_NO,genotype.t,genotype.s,group.c) %>% dplyr::summarise(sleep_sum_hr_t=sum(sleep_sum_hr)) %>%
dabest(genotype.s, sleep_sum_hr_t,
# The idx below passes "Control" as the control group,
# and "Group1" as the test group. The mean difference
# will be computed as mean(Group1) - mean(Control1).
idx = geno_list,
paired = FALSE) %>% mean_diff() %>% plot(
.,
palette = c("#999999", "#D55E00"),
color.column = group.c,
show.legend=FALSE,
rawplot.ylabel = "Total sleep (hr)",
rawplot.ylim = c(0, 24),
effsize.ylim = c(-6, 6)
)
g_d_sleep_light <- s_data_t %>% filter(light==1) %>%
dabest(genotype.s, sleep_sum_hr,
# The idx below passes "Control" as the control group,
# and "Group1" as the test group. The mean difference
# will be computed as mean(Group1) - mean(Control1).
idx = geno_list,
paired = FALSE) %>% mean_diff() %>%
plot(
.,
palette = c("#999999", "#D55E00"),
color.column = group.c,
show.legend=FALSE,
rawplot.ylabel = "Sleep of Daytime (hr)",
rawplot.ylim = c(0, 12),
effsize.ylim = c(-4, 4)
)
g_d_sleep_dark <- s_data_t %>% filter(light==0) %>%
dabest(genotype.s, sleep_sum_hr,
# The idx below passes "Control" as the control group,
# and "Group1" as the test group. The mean difference
# will be computed as mean(Group1) - mean(Control1).
idx = geno_list,
paired = FALSE) %>% mean_diff() %>%
plot(
.,
palette = c("#999999", "#D55E00"),
color.column = group.c,
show.legend=FALSE,
rawplot.ylabel = "Sleep of Nighttime (hr)",
rawplot.ylim = c(0, 12),
effsize.ylim = c(-4, 4)
)
g_d_sleep_total
g_d_sleep_light
g_d_sleep_dark
#Output graph of sleep bout duration
sum_sbd_data_t <- sbd_data_t %>% group_by(genotype, F_NO,genotype.t,genotype.s,group.c) %>%
dplyr::summarise(sbd_t=sum(sbd_mean))
max_sbd_t <- max(sum_sbd_data_t$sbd_t)
max_sbd_p <- max(sbd_data_t$sbd_mean)
g_d_sbd_total <- sum_sbd_data_t%>%
dabest(genotype.s, sbd_t,
# The idx below passes "Control" as the control group,
# and "Group1" as the test group. The mean difference
# will be computed as mean(Group1) - mean(Control1).
idx = geno_list,
paired = FALSE) %>% mean_diff() %>% plot(
.,
palette = c("#999999", "#D55E00"),
color.column = group.c,
show.legend=FALSE,
rawplot.ylabel = "Sbd of Total",
rawplot.ylim = c(0, max_sbd_t),
effsize.ylim = c(-160, 160)
)
g_d_sbd_light <- sbd_data_t %>% filter(exp_period=="light") %>%
dabest(genotype.s, sbd_mean,
# The idx below passes "Control" as the control group,
# and "Group1" as the test group. The mean difference
# will be computed as mean(Group1) - mean(Control1).
idx = geno_list,
paired = FALSE) %>% mean_diff() %>%
plot(
.,
palette = c("#999999", "#D55E00"),
color.column = group.c,
show.legend=FALSE,
rawplot.ylabel = "Sbd of Daytime",
rawplot.ylim = c(0, max_sbd_p),
effsize.ylim = c(-80, 80)
)
g_d_sbd_dark <- sbd_data_t %>% filter(exp_period=="dark") %>%
dabest(genotype.s, sbd_mean,
# The idx below passes "Control" as the control group,
# and "Group1" as the test group. The mean difference
# will be computed as mean(Group1) - mean(Control1).
idx = geno_list,
paired = FALSE) %>% mean_diff() %>%
plot(
.,
palette = c("#999999", "#D55E00"),
color.column = group.c,
show.legend=FALSE,
rawplot.ylabel = "Sbd of Nighttime",
rawplot.ylim = c(0, max_sbd_p),
effsize.ylim = c(-80, 80)
)
g_d_sbd_total
g_d_sbd_light
g_d_sbd_dark
#Output graph of sleep bout
sum_sb_data_t <- sb_data_t %>% group_by(genotype, F_NO,genotype.t,genotype.s,group.c) %>%
dplyr::summarise(sb_t=sum(n_perday))
max_sb_t <- max(sum_sb_data_t$sb_t)
max_sb_p <- max(sb_data_t$n_perday)
g_d_sb_total <- sum_sb_data_t%>%
dabest(genotype.s, sb_t,
# The idx below passes "Control" as the control group,
# and "Group1" as the test group. The mean difference
# will be computed as mean(Group1) - mean(Control1).
idx = geno_list,
paired = FALSE) %>% mean_diff() %>% plot(
.,
palette = c("#999999", "#D55E00"),
color.column = group.c,
show.legend=FALSE,
rawplot.ylabel = "Sleep bout of Total",
rawplot.ylim = c(0, max_sb_t),
effsize.ylim = c(-60, 60)
)
g_d_sb_light <- sb_data_t %>% filter(exp_period=="light") %>%
dabest(genotype.s, n_perday,
# The idx below passes "Control" as the control group,
# and "Group1" as the test group. The mean difference
# will be computed as mean(Group1) - mean(Control1).
idx = geno_list,
paired = FALSE) %>% mean_diff() %>%
plot(
.,
palette = c("#999999", "#D55E00"),
color.column = group.c,
show.legend=FALSE,
rawplot.ylabel = "Sleep bout of Daytime",
rawplot.ylim = c(0, max_sb_p),
effsize.ylim = c(-30, 30)
)
g_d_sb_dark <- sb_data_t %>% filter(exp_period=="dark") %>%
dabest(genotype.s, n_perday,
# The idx below passes "Control" as the control group,
# and "Group1" as the test group. The mean difference
# will be computed as mean(Group1) - mean(Control1).
idx = geno_list,
paired = FALSE) %>% mean_diff() %>%
plot(
.,
palette = c("#999999", "#D55E00"),
color.column = group.c,
show.legend=FALSE,
rawplot.ylabel = "Sleep bout of Nighttime",
rawplot.ylim = c(0, max_sb_p),
effsize.ylim = c(-30, 30)
)
g_d_sb_total
g_d_sb_light
g_d_sb_dark
now_time <- Sys.time() %>% as.POSIXct()
today_date <- paste0(year(now_time),month(now_time),day(now_time),hour(now_time),minute(now_time))
output_path <- paste0(getwd(),"/Output/",today_date)
dir.create(output_path)
#DAM data dabest analysis -sleep
ggsave(file = paste0(output_path,"/","g_d_sleep_total.png"), plot = g_d_sleep_total, dpi = 600, width =12, height = 6) #total sleep
ggsave(file = paste0(output_path,"/","g_d_sleep_light.png"), plot = g_d_sleep_light, dpi = 600, width =12, height = 6) #daytime sleep
ggsave(file = paste0(output_path,"/","g_d_sleep_dark.png"), plot = g_d_sleep_dark, dpi = 600, width =12, height = 6) #nighttime sleep
#DAM data dabest analysis -sleep bout duration
ggsave(file = paste0(output_path,"/","g_d_sbd_total.png"), plot = g_d_sbd_total, dpi = 600, width =12, height = 6) #total sleep bout duration
ggsave(file = paste0(output_path,"/","g_d_sbd_light.png"), plot = g_d_sbd_light, dpi = 600, width =12, height = 6) #daytime sleep bout duration
ggsave(file = paste0(output_path,"/","g_d_sbd_dark.png"), plot = g_d_sbd_dark, dpi = 600, width =12, height = 6) #nighttime sleep bout duration
#DAM data dabest analysis -sleep bout
ggsave(file = paste0(output_path,"/","g_d_sb_total.png"), plot = g_d_sb_total, dpi = 600, width =12, height = 6) #total sleep bout duration
ggsave(file = paste0(output_path,"/","g_d_sb_light.png"), plot = g_d_sb_light, dpi = 600, width =12, height = 6) #daytime sleep bout duration
ggsave(file = paste0(output_path,"/","g_d_sb_dark.png"), plot = g_d_sb_dark, dpi = 600, width =12, height = 6) #nighttime sleep bout duration
t6<-Sys.time()
print(paste0("Complete! Analysised data save directory: ",output_path))
## [1] "Complete! Analysised data save directory: /Users/leelv/Documents/Programing/R/2s_merge sleep v1.2/Output/20222252220"
print(paste0("Analysis took ",round(difftime(t6,t5),1)," seconds"))
## [1] "Analysis took 50.2 seconds"