diff --git a/.Rhistory b/.Rhistory
index eafeceb68ebcb3a4f99f5521ef509431f3f24baf..97ff8ae14324425df3fa0294343461689c85ea38 100644
--- a/.Rhistory
+++ b/.Rhistory
@@ -1,170 +1,3 @@
-library(here)
-load(here("./runs/GC-MS normalization/Image.RData"))
-library(openxlsx)
-library(tidyverse)
-library(car)
-library(pheatmap)
-library(broom)
-library(ggpubr)
-library(viridisLite)
-library(modelr)
-#library(dlookr)
-#library(imputeLCMD)
-library(ggrepel)
-sam_dat_dic <- colnames(sam_dat) %>%
-as_tibble()
-orig <- colnames(sam_dat1) %>%
-as_tibble()
-break_lines <- function(x){
-group = vector(mode = "list", length = length(x))
-#id <- vector(mode = "character", length = length(x))
-#isa_unit <-
-#term_source <- vector(mode = "character", length = length(x))
-#term_acc <- vector(mode = "character", length = length(x))
-for(i in seq_along(x)){
-if(str_detect(x[[i]], "Characteristic|Parameter|Factor|Component")){
-group[[i]] <- tibble(id = x[[i]])
-if(str_detect(x[[i+1]], "Unit")) {
-group[[i]] <- group[[i]] %>% mutate(isa_unit = x[[i+1]])
-if(str_detect(x[[i+2]], "Term Source REF")) {
-group[[i]] <- group[[i]] %>% mutate(term_source = x[[i+2]])
-if(str_detect(x[[i+3]], "Term Accession Number")){
-group[[i]] <- group[[i]] %>% mutate(term_acc = x[[i+3]])
-}
-}
-} else if(str_detect(x[[i+1]], "Term Source REF")) {
-group[[i]] <- group[[i]] %>% mutate(term_source = x[[i+1]])
-if(str_detect(x[[i+2]], "Term Accession Number")){
-group[[i]] <- group[[i]] %>% mutate(term_acc = x[[i+2]])
-}
-}
-} else {
-#id[[i]] <- x[[i]]
-}
-}
-group_tib <- group %>%
-discard(.p = is.null) %>%
-reduce(.f = bind_rows)
-#test_tib <- tibble(id = id,
-# unit = isa_unit,
-# term_source = term_source,
-# term_acc = term_acc) %>%
-# filter()
-assign("output",group, pos = 1)
-}
-break_lines(orig$value)
-break_lines <- function(x){
-group = vector(mode = "list", length = length(x))
-#id <- vector(mode = "character", length = length(x))
-#isa_unit <-
-#term_source <- vector(mode = "character", length = length(x))
-#term_acc <- vector(mode = "character", length = length(x))
-for(i in seq_along(x)){
-if(str_detect(x[[i]], "Characteristic|Parameter|Factor|Component")){
-group[[i]] <- tibble(id = x[[i]])
-if(str_detect(x[[i+1]], "Unit")) {
-group[[i]] <- group[[i]] %>% mutate(isa_unit = x[[i+1]])
-if(str_detect(x[[i+2]], "Term Source REF")) {
-group[[i]] <- group[[i]] %>% mutate(term_source = x[[i+2]])
-if(str_detect(x[[i+3]], "Term Accession Number")){
-group[[i]] <- group[[i]] %>% mutate(term_acc = x[[i+3]])
-}
-}
-} else if(str_detect(x[[i+1]], "Term Source REF")) {
-group[[i]] <- group[[i]] %>% mutate(term_source = x[[i+1]])
-if(str_detect(x[[i+2]], "Term Accession Number")){
-group[[i]] <- group[[i]] %>% mutate(term_acc = x[[i+2]])
-}
-}
-} else {
-#id[[i]] <- x[[i]]
-}
-}
-group_tib <- group %>%
-discard(.p = is.null) %>%
-reduce(.f = bind_rows)
-#test_tib <- tibble(id = id,
-# unit = isa_unit,
-# term_source = term_source,
-# term_acc = term_acc) %>%
-# filter()
-assign("output",group_tib, pos = 1)
-}
-break_lines(orig$value)
-View(output)
-str_detect(c("test", "different", "capture", "groups"), "test")
-str_detect(c("test", "different", "capture", "groups"), "test|different")
-break_lines <- function(x){
-group = vector(mode = "list", length = length(x))
-for(i in seq_along(x)){
-if(str_detect(x[[i]], "Characteristic|Parameter|Factor|Component")){
-group[[i]] <- tibble(id = x[[i]])
-if(str_detect(x[[i+1]], "Unit")) {
-group[[i]] <- group[[i]] %>% mutate(isa_unit = x[[i+1]])
-if(str_detect(x[[i+2]], "Term Source REF")) {
-group[[i]] <- group[[i]] %>% mutate(term_source = x[[i+2]])
-if(str_detect(x[[i+3]], "Term Accession Number")){
-group[[i]] <- group[[i]] %>% mutate(term_acc = x[[i+3]])
-}
-}
-} else if(str_detect(x[[i+1]], "Term Source REF")) {
-group[[i]] <- group[[i]] %>% mutate(term_source = x[[i+1]])
-if(str_detect(x[[i+2]], "Term Accession Number")){
-group[[i]] <- group[[i]] %>% mutate(term_acc = x[[i+2]])
-}
-}
-} else if(str_detect(x[[i]], "Source REF|Accession Number|Unit",negate = T)){
-group[[i]] <- tibble(id = x[[i]])
-}
-}
-group_tib <- group %>%
-discard(.p = is.null) %>%
-reduce(.f = bind_rows)
-assign("output",group_tib, pos = 1)
-}
-break_lines(orig$value)
-View(output)
-here()
-View(sam_dat_dic)
-View(sam_dat)
-View(isa_ext)
-View(met_dat)
-break_lines <- function(x){
-x <- colnames(x) %>%
-as_tibble()
-x <- x$value
-group = vector(mode = "list", length = length(x))
-for(i in seq_along(x)){
-if(str_detect(x[[i]], "Characteristic|Parameter|Factor|Component")){
-group[[i]] <- tibble(id = x[[i]])
-if(str_detect(x[[i+1]], "Unit")) {
-group[[i]] <- group[[i]] %>% mutate(isa_unit = x[[i+1]])
-if(str_detect(x[[i+2]], "Term Source REF")) {
-group[[i]] <- group[[i]] %>% mutate(term_source = x[[i+2]])
-if(str_detect(x[[i+3]], "Term Accession Number")){
-group[[i]] <- group[[i]] %>% mutate(term_acc = x[[i+3]])
-}
-}
-} else if(str_detect(x[[i+1]], "Term Source REF")) {
-group[[i]] <- group[[i]] %>% mutate(term_source = x[[i+1]])
-if(str_detect(x[[i+2]], "Term Accession Number")){
-group[[i]] <- group[[i]] %>% mutate(term_acc = x[[i+2]])
-}
-}
-} else if(str_detect(x[[i]], "Source REF|Accession Number|Unit",negate = T)){
-group[[i]] <- tibble(id = x[[i]])
-}
-}
-group_tib <- group %>%
-discard(.p = is.null) %>%
-reduce(.f = bind_rows)
-assign("output",group_tib, pos = 1)
-}
-break_lines(sam_dat1)
-break_lines <- function(x, output_name){
-x <- colnames(x) %>%
-as_tibble()
-x <- x$value
group = vector(mode = "list", length = length(x))
for(i in seq_along(x)){
if(str_detect(x[[i]], "Characteristic|Parameter|Factor|Component")){
@@ -220,3 +53,460 @@ here(out,
View(isa_ms)
View(isa_ms_tidy)
1176-1145
+source("C:/Users/Micha/Nextcloud/Data_Plant/cmqtl_val1_arc/workflows/GC_MS_normalization/210927_primary_normalization_with_split.R", echo=TRUE)
+#PCA testing chunk####
+##https://www.intechopen.com/books/metabolomics-fundamentals-and-applications/processing-and-visualization-of-metabolomics-data-using-r
+#exp1
+pca_base <- features_all %>%
+ungroup() %>%
+filter(exp == 1) %>%
+arrange(batch_pca, daily_num) %>%
+mutate(log_imp = log2(imp)) %>%
+pivot_wider(id_cols = c(tissue, treatment, batch_GC, daily_num, class, machine_num_GC),
+names_from = met,
+values_from = linnorm_fw_log) %>%
+filter(!is.na(m_1))
+pca <- summary(prcomp(
+pca_base%>%
+select(starts_with("m_"))))
+sam_vars <- pca_base %>% select(-starts_with("m_")) %>% colnames()
+pca_plot <- as_tibble(pca$x) %>%
+mutate(join_num = 1:nrow(pca$x)) %>%
+left_join(pca_base %>%
+mutate(join_num=1:nrow(pca$x))) %>%
+#left_join(sam_dat) %>%
+mutate(batch_GC= as_factor(batch_GC),
+class= as_factor(class)) %>%
+select(all_of(sam_vars), everything())
+exp_var <- as_tibble(pca[["importance"]])
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = tissue)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = tissue)) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+ggsave(here(out,
+"cmQTL_val1_tissue_lin_PCA.jpg"
+),
+width = 183,
+height = 100,
+units = "mm",
+dpi = 300)
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = treatment)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = treatment)) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+ggsave(here(out,
+"cmQTL_val1_treatment_PCA.jpg"
+),
+width = 183,
+height = 100,
+units = "mm",
+dpi = 300)
+dev.off()
+pca_base <- features_all %>%
+ungroup() %>%
+filter(exp == 1) %>%
+arrange(batch_pca, daily_num) %>%
+mutate(log_imp = log2(imp)) %>%
+pivot_wider(id_cols = c(tissue, treatment, batch_GC, daily_num, class, machine_num_GC),
+names_from = met,
+values_from = loess_norm_fw_log) %>%
+filter(!is.na(m_1))
+pca <- summary(prcomp(
+pca_base%>%
+select(starts_with("m_"))))
+sam_vars <- pca_base %>% select(-starts_with("m_")) %>% colnames()
+pca_plot <- as_tibble(pca$x) %>%
+mutate(join_num = 1:nrow(pca$x)) %>%
+full_join(pca_base %>%
+mutate(join_num=1:nrow(pca$x))) %>%
+#left_join(sam_dat) %>%
+mutate(batch_pca= as_factor(batch_GC),
+class= as_factor(class)) %>%
+select(all_of(sam_vars), everything())
+exp_var <- as_tibble(pca[["importance"]])
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = tissue)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = tissue)) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+ggsave(here(out,
+"cmQTL_val1_tissue_loess_PCA.jpg"
+),
+width = 183,
+height = 100,
+units = "mm",
+dpi = 300)
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = treatment)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = treatment)) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+ggsave(here(out,
+"cmQTL_val1_treatment_loess_PCA.jpg"
+),
+width = 183,
+height = 100,
+units = "mm",
+dpi = 300)
+pca_base <- features_all %>%
+ungroup() %>%
+filter(exp == 2) %>%
+arrange(batch_pca, daily_num) %>%
+mutate(log_imp = log2(isnorm)) %>%
+pivot_wider(id_cols = c(tissue, treatment, batch_GC, daily_num, class, machine_num_GC),
+names_from = met,
+values_from = log_imp)
+pca <- summary(prcomp(
+pca_base%>%
+select(starts_with("m_"))))
+sam_vars <- pca_base %>% select(-starts_with("m_")) %>% colnames()
+pca_plot <- as_tibble(pca$x) %>%
+mutate(join_num = 1:nrow(pca$x)) %>%
+full_join(pca_base %>%
+mutate(join_num=1:nrow(pca$x))) %>%
+#left_join(sam_dat) %>%
+mutate(batch_pca= as_factor(batch_GC),
+class= as_factor(class)) %>%
+select(all_of(sam_vars), everything())
+exp_var <- as_tibble(pca[["importance"]])
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = class)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = class)) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = batch_GC)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = batch_GC)) +
+geom_text(aes(x = PC1, y = PC2,
+label = if_else(machine_num_GC %in% exclude_samples, machine_num_GC, ""))) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+ggsave(here(out,
+"cmQTL_val1_treatment_loess_outlier_exp2_loess_norm_PCA.jpg"
+),
+width = 183,
+height = 100,
+units = "mm",
+dpi = 300)
+for(tiss in seq_along(tissues)){
+pca_base <- features_all %>%
+ungroup() %>%
+filter(exp == 1, tissue == tissues[[tiss]]) %>%
+arrange(batch_pca, daily_num) %>%
+mutate(log_imp = log2(imp)) %>%
+pivot_wider(id_cols = c(tissue, treatment,batch_GC, daily_num, class, machine_num_GC),
+names_from = met,
+values_from = loess_norm_fw_log) %>%
+filter(!is.na(m_1))
+pca <- summary(prcomp(
+pca_base%>%
+select(starts_with("m_"))))
+sam_vars <- pca_base %>% select(-starts_with("m_")) %>% colnames()
+pca_plot <- as_tibble(pca$x) %>%
+mutate(join_num = 1:nrow(pca$x)) %>%
+full_join(pca_base %>%
+mutate(join_num=1:nrow(pca$x))) %>%
+left_join(sam_dat) %>%
+mutate(batch_pca= as_factor(batch_GC),
+class= as_factor(class)) %>%
+select(all_of(sam_vars), everything())
+exp_var <- as_tibble(pca[["importance"]])
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = treatment)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = treatment)) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+ggsave(here(out,
+"cmQTL_val1", tissues[[tiss]], "treatment_loess_PCA.jpg"
+),
+width = 183,
+height = 100,
+units = "mm",
+dpi = 300)
+}
+sam_vars <- pca_base %>% select(-starts_with("m_")) %>% colnames()
+pca_plot <- as_tibble(pca$x) %>%
+mutate(join_num = 1:nrow(pca$x)) %>%
+full_join(pca_base %>%
+mutate(join_num=1:nrow(pca$x))) %>%
+left_join(sam_dat) %>%
+mutate(batch_pca= as_factor(batch_pca),
+class= as_factor(class)) %>%
+select(all_of(sam_vars), everything())
+exp_var <- as_tibble(pca[["importance"]])
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = tissue)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = tissue)) +
+geom_text(aes(x = PC1, y = PC2,
+label = if_else(machine_num_GC == "21106rA_31", machine_num_GC, "")),
+nudge_y = -2,
+nudge_x = 4) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+ggsave(here(out,
+"cmQTL_val1_tissue_imp_PCA.jpg"
+),
+width = 183,
+height = 100,
+units = "mm",
+dpi = 300)
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = batch_GC)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = batch_GC)) +
+geom_text(aes(x = PC1, y = PC2,
+label = if_else(machine_num_GC == "21106rA_31", machine_num_GC, "")),
+nudge_y = -2,
+nudge_x = 4) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+ggsave(here(out,
+"cmQTL_val1_batch_imp_PCA.jpg"
+),
+width = 183,
+height = 100,
+units = "mm",
+dpi = 300)
+pca_plot %>%
+# filter(treatment == "HL") %>%
+ggplot()+
+geom_jitter(aes(x=PC1, y=PC2, color = class)) +
+stat_ellipse(aes(x=PC1, y=PC2, color = class)) +
+geom_text(aes(x = PC1, y = PC2,
+label = if_else(machine_num_GC == "21106rA_31", machine_num_GC, "")),
+nudge_y = -2,
+nudge_x = 4) +
+ylab(str_c("PC2 ", "(", exp_var$PC2[[2]]*100, "%)")) +
+xlab(str_c("PC1 ", "(", exp_var$PC1[[2]]*100, "%)"))
+ggsave(here(out,
+"cmQTL_val1_class_imp_PCA.jpg"
+),
+width = 183,
+height = 100,
+units = "mm",
+dpi = 300)
+#cmQTLval2
+pca_base <- features_all %>%
+ungroup() %>%
+filter(exp == 2, batch_GC == "3_split") %>%
+arrange(batch_GC, daily_num) %>%
+mutate(log_imp = log2(isnorm)) %>%
+pivot_wider(id_cols = c(tissue, treatment, batch_GC, daily_num, class, machine_num_GC),
+names_from = met,
+values_from = log_imp)%>%
+filter(!is.na(m_1))
+pca <- summary(prcomp(
+pca_base%>%
+select(starts_with("m_"))))
+sam_vars <- pca_base %>% select(-starts_with("m_")) %>% colnames()
+skip <- T
+if(skip == T) {
+} else{
+# Plot all metabolites ----------------------------------------------------
+plotmets <- features_all %>% distinct(met) %>% as_vector
+plotmet_labs <- plotmets %>% as_tibble() %>%
+left_join(met_dat, by= c("value" = "met")) %>%
+#mutate(peak_num = base::rank(HMDB_clear_name, ties.method = "first"),
+# dup = sum(peak_num),
+# HMDB_clear_name_unique = if_else(dup>1, str_c(HMDB_clear_name, peak_num), HMDB_clear_name)) %>%
+#ungroup() %>%
+select(Compound_Name) %>% as_vector()
+plot_out <- vector("list", length = length(plotmets))
+for (meta in seq_along(plotmets)) {
+plot_out [[meta]] <- features_all %>%
+mutate(xint = if_else(daily_num == 4, run_num_GC-3.5, max(run_num_GC)),
+is_miss = as_factor(if_else(is.na(area), T, F))) %>%
+filter(met == plotmets[[meta]]) %>%
+filter(class!="blank") %>%
+ggplot(aes(x=run_num_GC, y=isnorm)) +
+geom_point(aes(color = class, shape = is_miss)) +
+geom_point(aes(y=predict), color="black", size=0.1) +
+geom_vline(aes(xintercept = xint))+
+# facet_grid(rows = vars(treatment), cols = vars(rep), scales = "free") +
+ggtitle(label = plotmet_labs[[meta]])
+}
+pdf(file = here(out,"_loess_fit.pdf"))
+for (meta in seq_along(plotmets)) {
+print(plot_out[[meta]])
+}
+dev.off()
+plot_out <- vector("list", length = length(plotmets))
+for (meta in seq_along(plotmets)) {
+plot_out [[meta]] <- features_all %>%
+mutate(xint = if_else(daily_num == 4, run_num_GC-3.5, max(run_num_GC))) %>%
+filter(met == plotmets[[meta]]) %>%
+filter(class!="blank") %>%
+ggplot(aes(x=run_num_GC, y=loess_norm_fw_log)) +
+geom_point(aes(color=class)) +
+#geom_point(aes(y=predp), color="black", size=0.1) +
+geom_vline(aes(xintercept = xint))+
+# facet_grid(rows = vars(treatment), cols = vars(rep), scales = "free") +
+ggtitle(label = plotmet_labs[[meta]])
+}
+pdf(file = here(out,"_loess_norm_fw_log.pdf"))
+for (meta in seq_along(plotmets)) {
+print(plot_out[[meta]])
+}
+dev.off()
+for (meta in seq_along(plotmets)) {
+plot_out [[meta]] <- features_all %>%
+mutate(is_miss = as_factor(if_else(is.na(area), T, F))) %>%
+filter(met == plotmets[[meta]]) %>%
+filter(class!="blank") %>%
+ggplot(aes(x=daily_num, y=isnorm)) +
+geom_point(aes(color = class, shape = is_miss)) +
+geom_point(aes(y=predict), color="black", size=0.1) +
+facet_wrap(facets = vars(batch_GC), nrow = 1, scales = "free") +
+ggtitle(label = plotmet_labs[[meta]]) +
+theme(axis.ticks.y = element_blank(),
+axis.text.y = element_blank(),
+panel.spacing = unit(-0.1,"lines"))
+}
+pdf(file = here(out,"_loess_norm_no_scale.pdf"),
+height = 3.5)
+for (meta in seq_along(plotmets)) {
+print(plot_out[[meta]])
+}
+dev.off()
+for (meta in seq_along(plotmets)) {
+plot_out [[meta]] <- features_all %>%
+mutate(xint = if_else(daily_num == 4, run_num_GC-3.5, max(run_num_GC)),
+is_miss = as_factor(if_else(is.na(area), T, F))) %>%
+filter(met == plotmets[[meta]]) %>%
+filter(class!="blank") %>%
+ggplot(aes(x=run_num_GC, y=isnorm)) +
+geom_point(aes(color = class, shape = is_miss)) +
+geom_point(aes(y=predp), color="black", size=0.1) +
+geom_vline(aes(xintercept = xint))+
+# facet_grid(rows = vars(treatment), cols = vars(rep), scales = "free") +
+ggtitle(label = plotmet_labs[[meta]])
+}
+pdf(file = here(out,"_linnorm.pdf"))
+for (meta in seq_along(plotmets)) {
+print(plot_out[[meta]])
+}
+dev.off()
+rescale <- features_all %>%
+group_by(met, tissue, exp) %>%
+summarise(rescale = median(imp))
+sam_vars <- c(#"plantline","extraction_num", "sample_num", #in case I add them back
+"alias", "LIMS_ID",
+"treatment", "tissue", "batch_GC", "run_date_GC",
+"machine_num_GC",
+"class", "run_num_GC", "sample_weight", "exp", "genotype")
+features_out <- features_all %>%
+left_join(rescale) %>%
+filter(class != "blank",
+(component %in% take_split & str_detect(batch_GC, "split"))|
+(!component %in% take_split & str_detect(batch_GC, "split", negate = T))|
+tissue == "leaves") %>%
+group_by(met, tissue, exp) %>%
+mutate(loess_norm_med = loess_norm_fw/median(loess_norm_fw),
+rescaled = loess_norm_med*rescale) %>%
+select(all_of(sam_vars), met, Compound_Name, loess_norm_fw, loess_norm_med, area, rescaled)
+#features_out %>%
+# filter(met == "m_02177", tissue == "fruits", exp == 1) %>%
+# ggplot(aes(x = run_num_GC, y = loess_norm_fw)) +
+# geom_point()
+}
+write_csv(features_out,
+here(out,
+"cmQTL_val_1_2_feat_dat_GC.csv"
+))
+here()
+write_csv(features_out,
+here(out,
+"cmQTL_val_1_2_feat_dat_GC.csv"
+))
+features_out
+write_csv(features_out,
+here(out,
+"cmQTL_val_1_2_feat_dat_GC.csv"
+))
+getwd()
+here()
+write_csv(features_out,
+here(out,
+"cmQTL_val_1_2_feat_dat_GC.csv"
+))
+here::dr_here()
+here::i_am("workflows/GC_MS_normalization/210927_primary_normalization_with_split.R")
+write_csv(features_out,
+here(out,
+"cmQTL_val_1_2_feat_dat_GC.csv"
+))
+here(out,
+"cmQTL_val_1_2_feat_dat_GC.csv"
+)
+out
+out <- here("runs/GC-MS normalization")
+out
+write_csv(features_out,
+here(out,
+"cmQTL_val_1_2_feat_dat_GC.csv"
+))
+write_csv(met_dat,
+here(out,
+"cmQTL_val_1_2_met_dat_GC.csv"
+))
+write_csv(sam_dat,
+here(out,
+"cmQTL_val_1_2_sam_dat_GC.csv"
+))
+write_csv(missingness,
+here(out,
+"cmQTL_val_1_2_miss_GC.csv"
+))
+View(features_out)
+source("C:/Users/Micha/Nextcloud/Data_Plant/cmqtl_val1_arc/workflows/GC_MS_normalization/210927_primary_normalization_with_split.R", echo=TRUE)
+View(features_out)
+View(features_all)
+features_out <- features_all %>%
+left_join(rescale) %>%
+filter(class != "blank",
+(component %in% take_split & str_detect(batch_GC, "split"))|
+(!component %in% take_split & str_detect(batch_GC, "split", negate = T))|
+tissue == "leaves") %>%
+group_by(met, tissue, exp) %>%
+mutate(loess_norm_med = loess_norm_fw/median(loess_norm_fw),
+rescaled = loess_norm_med*rescale) %>