Skip to content

Commit

Permalink
Merge pull request #154 from MoTrPAC/develop
Browse files Browse the repository at this point in the history
v0.6.9: metab density plots
  • Loading branch information
biodavidjm authored Jul 9, 2022
2 parents 0f8ffd4 + c3c2b09 commit 3c5d886
Show file tree
Hide file tree
Showing 84 changed files with 1,273 additions and 720 deletions.
8 changes: 4 additions & 4 deletions DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,8 +1,8 @@
Package: MotrpacBicQC
Type: Package
Title: QC/QA functions for the MoTrPAC community
Version: 0.6.8
Date: 2022-04-11
Version: 0.6.9
Date: 2022-07-09
Author: MoTrPAC Bioinformatics Center
Maintainer: David Jimenez-Morales <davidjm@stanford.edu>
Description: R Package for the analysis of MoTrPAC datasets.
Expand Down Expand Up @@ -41,10 +41,10 @@ Imports:
viridis
NeedsCompilation: no
VignetteBuilder: knitr
RoxygenNote: 7.1.2
RoxygenNote: 7.2.0
Roxygen: list(markdown = TRUE)
Suggests:
testthat,
markdown,
testthat,
rmarkdown,
rmdformats
1 change: 1 addition & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@ export(merge_metabolomics_metadata)
export(merge_phenotype_metabolomics)
export(open_file)
export(plot_basic_metabolomics_qc)
export(proteomics_plots_rii)
export(remove_empty_columns)
export(remove_empty_rows)
export(set_phase)
Expand Down
7 changes: 7 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
@@ -1,3 +1,10 @@
# MotrpacBicQC 0.6.9 (2022-07-09)

* Metabolomics: new density plots
* Code optimizations
* Bug fixes affecting file manifest checks


# MotrpacBicQC 0.6.8 (2022-04-11)

* Refactor the DMAQC validation. A new file will be required when:
Expand Down
45 changes: 40 additions & 5 deletions R/metabolomics_plots.R
Original file line number Diff line number Diff line change
Expand Up @@ -48,9 +48,37 @@ plot_basic_metabolomics_qc <- function(results,
labs(title = "MZ/RT density map",
subtitle = paste(output_prefix))
}


# Density plots------

if(verbose) message(" - (p) Density distributions")
mu <- results_long %>% group_by(sample_id) %>% summarise(grp.mean=mean(intensity))

den1 <- ggplot(data = results_long, aes(x = log2(intensity), color = sample_type)) +
geom_density(na.rm = TRUE) +
theme_light() +
labs(title = "Density distribution by sample type",
subtitle = paste(output_prefix))

den2 <- ggplot(data = results_long, aes(x = log2(intensity),
color = sample_id)) +
geom_density(na.rm = TRUE) +
theme_light() +
theme(legend.position="none") +
labs(title = "Density distribution by sample id",
subtitle = paste(output_prefix))

den3 <- ggplot(data = results_long, aes(x = log2(intensity),
color = sample_id)) +
geom_density(na.rm = TRUE) +
theme_minimal() +
theme(legend.position="none") +
facet_wrap(~sample_type) +
labs(title = "Density distributions",
subtitle = paste(output_prefix))

# Plot: sum of intensities------
if(verbose) message(" - (p) Plot sum of intensity/concentration values")
sum_int <- results_long %>%
group_by(sample_id, sample_type, sample_order) %>%
summarise(sum_quant = sum(intensity))
Expand Down Expand Up @@ -143,7 +171,6 @@ plot_basic_metabolomics_qc <- function(results,
piseio <- piseio + labs(title = "Sample Concentration distribution",
y = "log2(Concentration)")
}


# Plot number fo IDs by identification------
if(verbose) message(" - (p) Plot ID counts")
Expand Down Expand Up @@ -221,7 +248,9 @@ plot_basic_metabolomics_qc <- function(results,
tns <- unique(uid2[c("sample_id", "sample_type")]) %>% group_by(sample_type) %>% count(sample_type)
ptns <- ggplot(tns, aes(x = sample_type, y = n, fill = sample_type)) +
geom_bar(stat = "identity") +
labs(title = "Number of samples by sample type", y = "", x = "") +
labs(title = "Number of samples by sample type",
subtitle = paste(output_prefix),
y = "", x = "") +
theme_minimal() +
theme(legend.position = "none",
axis.text.x = element_text(size = 14),
Expand All @@ -242,7 +271,9 @@ plot_basic_metabolomics_qc <- function(results,
label = scales::percent((..count..)/sum(..count..))),
stat = "count", vjust = -0.25) +
scale_y_continuous(labels = percent) +
labs(title = "Proportion of Features Identified (named/unnamed)", y = "", x = "") +
labs(title = "Proportion of Features Identified (named/unnamed)",
subtitle = paste(output_prefix),
y = "", x = "") +
theme_light() +
theme(legend.position = "none",
axis.text.x = element_text(size = 18)) +
Expand All @@ -251,7 +282,8 @@ plot_basic_metabolomics_qc <- function(results,
pnids <- ggplot(results, aes(x = id_type, fill = id_type)) +
geom_bar() +
geom_text(stat = 'count',aes(label =..count.., vjust = -0.2)) +
labs(title = "Total Number of Features Identified (named/unnamed)", y = "", x = "") +
labs(title = "Total Number of Features Identified (named/unnamed)",
subtitle = paste(output_prefix), y = "", x = "") +
theme_light() +
theme(legend.position = "none",
axis.text.x = element_text(size = 18)) +
Expand Down Expand Up @@ -306,6 +338,9 @@ plot_basic_metabolomics_qc <- function(results,
if(printPDF) pdf(out_plot_summary, width = 12, height = 6)
print(ptns)
gridExtra::grid.arrange(ppids, pnids, ncol = 2)
print(den1)
print(den2)
print(den3)
if(!is.null(metametab)){print(pmzrt)}
if(printPDF) garbage <- dev.off()
}
Expand Down
171 changes: 171 additions & 0 deletions R/proteomics_plots.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,171 @@

# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
#' @title QC Plot of proteomics reporter ion intensity data
#'
#' @description check whether results file is following guidelines
#' @param all_samples (vector) all sample ids
#' @param all_vial_labels (vector) only vial_labels
#' @param peprii (char) Reporter Ion Intensity df
#' @param isPTM (logical) Is a PTM assay
#' @param v_m (df) sample metadata
#' @param out_qc_folder (char) if `f_proof is TRUE`, a folder path can be provided
#' (otherwise print in current working directory)
#' @param output_prefix (char) provide a prefix for the output name
#' @param printPDF (logical) if `TRUE` (default print plots to pdf)
#' @param verbose (logical) `TRUE` (default) shows messages
#' @return (int) number of issues identified
#' @examples {
#' check_results(r_m = results_named, m_s = metadata_sample_named, m_m = metadata_metabolites_named)
#' }
#' @export
proteomics_plots_rii <- function(all_vial_labels,
all_samples,
peprii,
isPTM,
v_m,
out_qc_folder = NULL,
output_prefix,
printPDF,
verbose){

if(verbose) message(" + (+) PLOTS RII------------------")

if( !is.null(all_vial_labels) ){
required_columns <- get_required_columns(isPTM = isPTM,
prot_file = "rii")
required_columns <- c(required_columns, all_samples)

if( all(required_columns %in% colnames(peprii)) ){
# Check distributions

if(isPTM){
r_c <- c("ptm_peptide", all_samples)
fpeprii <- subset(peprii, select = r_c)
peptides_long <- fpeprii %>% tidyr::pivot_longer(cols = -c(ptm_peptide),
names_to = "vial_label",
values_to = "ri_intensity")
}else{
r_c <- c("protein_id", "sequence", all_samples)
fpeprii <- subset(peprii, select = r_c)
peptides_long <- fpeprii %>% tidyr::pivot_longer(cols = -c(protein_id, sequence),
names_to = "vial_label",
values_to = "ri_intensity")
}

# Plot Intensity distribution-----
if(verbose) message(" - (p) Plot intensity distributions")

peptides_long <- merge(v_m, peptides_long, by = c("vial_label"))

peptides_long$vial_label <- as.character(peptides_long$vial_label)
peptides_long$vial_label <- as.factor(peptides_long$vial_label)
peptides_long$tmt_plex <- as.factor(peptides_long$tmt_plex)

pise <- ggplot2::ggplot(peptides_long,
aes(x = reorder(vial_label, log2(ri_intensity), FUN = median, na.rm = TRUE),
y = log2(ri_intensity),
fill = tmt_plex)) +
geom_boxplot(na.rm = TRUE) +
theme_linedraw() +
theme(axis.text.x = element_text(angle = 90,
hjust = 1,
vjust = 0.5,
size = 8)) + #legend.position = "none"
labs(x = "vial_label", y = "log2(ri_intensity)") +
labs(title = "Reporter Ion Intensity",
subtitle = output_prefix)


# Plot NA values------
if(verbose) message(" - (p) Plot NA values")

p_na_peprii <- peprii[required_columns] %>%
inspectdf::inspect_na() %>%
dplyr::arrange(match(col_name, colnames(peprii))) %>%
inspectdf::show_plot() +
ylim(0, 100) + theme_linedraw() +
theme(axis.text.x = element_text(angle = 90,
hjust = 1,
vjust = 0.5,
size = 8))

# Plot unique IDs-----
if(verbose) message(" - (p) Plot Unique IDs")

if(isPTM){
key_id <- "ptm_peptide"
}else{
key_id <- "protein_id"
}

uid <- peptides_long %>%
group_by(across(all_of(c(key_id, "vial_label", "tmt_plex")))) %>%
summarise(total_rii = ri_intensity, .groups = 'drop')
uid2 <- uid[which(!is.na(uid$total_rii)),]
uid3 <- unique(uid2[c(key_id, "vial_label", "tmt_plex")]) %>%
count(vial_label, tmt_plex)

puid1 <- ggplot(uid3, aes(x = reorder(vial_label, n), y = n, fill = tmt_plex)) +
geom_bar(stat = "identity") +
theme_linedraw() +
theme(
axis.text.x = element_text(
angle = 90,
hjust = 1,
vjust = 0.5,
size = 8
),
legend.position = "none"
) +
geom_text(
aes(label = n),
# vjust = -0.5,
hjust = 1,
size = 2.7,
angle = 90
) +
ggtitle("RII: Unique IDs in samples") +
facet_wrap(~ tmt_plex, scales = "free") +
xlab("Vial Labels")

puid2 <- ggplot(uid3, aes(x = reorder(vial_label, n), y = n, fill = tmt_plex)) +
geom_bar(stat = "identity",
na.rm = TRUE) +
theme_linedraw() +
theme(
axis.text.x = element_text(
angle = 90,
hjust = 1,
vjust = 0.5,
size = 8
)
) +
geom_text(
aes(label = n),
# vjust = -0.5,
hjust = 1,
size = 2.7,
angle = 90
) +
ggtitle("RII: Unique IDs in samples") +
xlab("Vial Labels")

if(is.null(out_qc_folder)){
out_plot_dist <- paste0(output_prefix,"-qc-rii-distribution.pdf")
}else{
out_plot_dist <- file.path(normalizePath(out_qc_folder), paste0(output_prefix,"-qc-rii-distribution.pdf"))
}

if(printPDF) pdf(out_plot_dist, width = 12, height = 8)
print(pise)
print(puid1)
print(puid2)
print(p_na_peprii)
if(printPDF) garbage <- dev.off()
}else{
if(verbose) message(" - (-) Not all required columns available in RII: print proof is not possible")
required_columns[!(required_columns %in% colnames(peprii))]
ic <- ic + 1
}
} # !is.null(all_vial_labels)
}
Loading

0 comments on commit 3c5d886

Please sign in to comment.