knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
library(pepdiff)

Load in data

d <- import_data(system.file("extdata", "anon.csv",package="pepdiff"),
                 gene_id = "gene_name",
                 peptide = "peptide_sequence",
                 treatment = "treatment_name",
                 )

Check number of missing data points

Includes at tech rep level

Missing data

missing_peptides_plot(d)

Replicate level

times_measured(d)
times_measured_plot(d)

Look at similarity of samples

plot_pca(d)

plot_kmeans(d)

Look at distribution of quantifications

plot_quant_distributions(d)
plot_quant_distributions(d, log = TRUE)
norm_qqplot(d)
norm_qqplot(d, log = TRUE)

Do a comparison

665e6428 0 seconds vs 665e6428 150 seconds

r <- compare(d, iters = 1000, tests=c("bootstrap_t", "wilcoxon", "rank_product"), 
             control = '665e6428', 
             c_seconds = '0', 
             treatment = '665e6428', 
             t_seconds = '150')

r

Look at distribution of calculated Fold Changes

plot_fc(r)
fc_qqplot(r)
plot_fc(r, log = TRUE)
fc_qqplot(r, log = TRUE)

Plot p-values

plot_result(r)

Test distribution of p-values for the different methods

p_value_hist(r )

Compare the significant peptides

compare_calls(r)
comparisons <- data.frame(
  control = c('665e6428','665e6428'),
  c_seconds = c(0,150),
  treatment = c('665e6428','665e6428'),
  t_seconds = c(150, 0)
)

many <- compare_many(d, comparisons, tests = c("bootstrap_t", "rank_product"))
many
plot_heatmap(many, metric = "bootstrap_t_p_val", log = TRUE,  
             col_order = c("665e6428_150-665e6428_0", "665e6428_0-665e6428_150"))

plot_heatmap(many, metric = "rank_prod_p1_p_val", log = TRUE)
#plot all peptides - even those that aren't significant
plot_heatmap(many, metric = "bootstrap_t_p_val", log = TRUE, sig_only = FALSE, 
             col_order = c("665e6428_150-665e6428_0", "665e6428_0-665e6428_150") )
volcano_plot(many, metric = "bootstrap_t")


TeamMacLean/pepdiff documentation built on April 24, 2023, 7:48 a.m.