plot contrasts
plot contrasts
Other modelling:
Contrasts,
ContrastsMissing,
ContrastsModerated,
ContrastsProDA,
ContrastsROPECA,
ContrastsTable,
INTERNAL_FUNCTIONS_BY_FAMILY,
LR_test(),
Model,
build_model(),
contrasts_fisher_exact(),
get_anova_df(),
get_complete_model_fit(),
get_p_values_pbeta(),
isSingular_lm(),
linfct_all_possible_contrasts(),
linfct_factors_contrasts(),
linfct_from_model(),
linfct_matrix_contrasts(),
merge_contrasts_results(),
model_analyse(),
model_summary(),
moderated_p_limma(),
moderated_p_limma_long(),
my_contest(),
my_contrast(),
my_contrast_V1(),
my_contrast_V2(),
my_glht(),
pivot_model_contrasts_2_Wide(),
plot_lmer_peptide_predictions(),
sim_build_models_lm(),
sim_build_models_lmer(),
sim_make_model_lm(),
sim_make_model_lmer(),
strategy_lmer(),
summary_ROPECA_median_p.scaled()
Other plotting:
INTERNAL_FUNCTIONS_BY_FAMILY,
UpSet_interaction_missing_stats(),
UpSet_missing_stats(),
medpolish_estimate_df(),
missigness_histogram(),
missingness_per_condition(),
missingness_per_condition_cumsum(),
plot_NA_heatmap(),
plot_estimate(),
plot_heatmap(),
plot_heatmap_cor(),
plot_heatmap_cor_iheatmap(),
plot_hierarchies_add_quantline(),
plot_hierarchies_boxplot_df(),
plot_hierarchies_line(),
plot_hierarchies_line_df(),
plot_intensity_distribution_violin(),
plot_pca(),
plot_raster(),
plot_sample_correlation(),
plot_screeplot()
contrastDFdata frame with contrasts
modelNameof column with model name
subject_Idhierarchy key columns
prefixdefault Contrasts - used to generate file names
diffcolumn with fold change differences
contrastcolumn with contrasts names, default "contrast"
volcano_specvolcano plot specification
score_specscore plot specification
histogram_specplot specification
fcthreshfold change threshold
avg.abundancename of column containing avg abundance values.
protein_annotprotein annotation
new()create Crontrast_Plotter
ContrastsPlotter$new(
contrastDF,
subject_Id,
volcano = list(list(score = "FDR", thresh = 0.1)),
histogram = list(list(score = "p.value", xlim = c(0, 1, 0.05)), list(score = "FDR",
xlim = c(0, 1, 0.05))),
score = list(list(score = "statistic", thresh = NULL)),
fcthresh = 1,
modelName = "modelName",
diff = "diff",
contrast = "contrast",
avg.abundance = "avgAbd",
protein_annot = NULL
)contrastDFframe with contrast data
subject_Idcolumns containing subject Identifier
volcanowhich score to plot and which ablines to add.
histogramwhich scores to plot and which range (x) should be shown.
scorescore parameters
fcthreshdefault 1 (log2 FC threshold)
modelNamename of column with model names
difffold change (difference) diff column
contrastcontrast column
avg.abundancename of column with average abundance
protein_annotadd protein annotation (optional)
volcano()volcano plots (fold change vs FDR)
ContrastsPlotter$volcano(
colour,
legend = TRUE,
scales = c("fixed", "free", "free_x", "free_y"),
min_score = 1e-04
)volcano_plotly()plotly volcano plots
ContrastsPlotter$volcano_plotly(
colour,
legend = TRUE,
scales = c("fixed", "free", "free_x", "free_y"),
min_score = 1e-04
)ma_plot()ma plot
MA plot displays the effect size estimate as a function of the mean protein intensity across groups. Each dot represents an observed protein. Red horizontal lines represent the fold-change threshold.
Sometimes measured effects sizes (differences between samples groups) are biased by the signal intensity (here protein abundance). Such systematic effects can be explored using MA-plots.
ma_plotly()ma plotly
score_plot()plot a score against the log2 fc e.g. t-statistic
score_plotly()plot a score against the log2 fc e.g. t-statistic
istar <- sim_lfq_data_peptide_config(Nprot = 100)
#> creating sampleName from fileName column
#> completing cases
#> completing cases done
#> setup done
modelName <- "Model"
modelFunction <-
strategy_lmer("abundance ~ group_ + (1 | peptide_Id) + (1|sample)",
model_name = modelName)
pepIntensity <- istar$data
config <- istar$config
mod <- build_model(
pepIntensity,
modelFunction,
modelName = modelName,
subject_Id = config$table$hierarchy_keys_depth())
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> boundary (singular) fit: see help('isSingular')
#> Warning: There were 22 warnings in `dplyr::mutate()`.
#> The first warning was:
#> ℹ In argument: `linear_model = purrr::map(data, model_strategy$model_fun, pb =
#> pb)`.
#> ℹ In group 4: `protein_Id = "3QLHfm~8938"`.
#> Caused by warning in `value[[3L]]()`:
#> ! WARN :Error: grouping factors must have > 1 sampled level
#> ℹ Run `dplyr::last_dplyr_warnings()` to see the 21 remaining warnings.
#> Joining with `by = join_by(protein_Id)`
Contr <- c("group_A_vs_Ctrl" = "group_A - group_Ctrl",
"group_B_vs_Ctrl" = "group_B - group_Ctrl"
)
contrast <- prolfqua::Contrasts$new(mod,
Contr)
contr <- contrast$get_contrasts()
#> determine linear functions:
#> compute contrasts:
#> computing contrasts.
#> Joining with `by = join_by(protein_Id, contrast)`
cp <- ContrastsPlotter$new(contr,
contrast$subject_Id,
volcano = list(list(score = "FDR", thresh = 0.1)),
histogram = list(list(score = "p.value", xlim = c(0,1,0.05)),
list(score = "FDR", xlim = c(0,1,0.05))),
score =list(list(score = "statistic", thresh = 5))
)
stopifnot("plotly" %in% class(cp$volcano_plotly()$FDR))
stopifnot("ggplot" %in% class(cp$score_plot(legend=FALSE)$statistic))
p <- cp$histogram()
stopifnot("ggplot" %in% class(p$FDR))
stopifnot("ggplot" %in% class(p$p.value))
p <- cp$histogram_estimate()
stopifnot("ggplot" %in% class(p))
res <- cp$volcano()
stopifnot("ggplot" %in% class(res$FDR))
respltly <- cp$volcano_plotly()
stopifnot("plotly" %in% class(respltly$FDR))
stopifnot("ggplot" %in% class(cp$ma_plot()))
stopifnot("plotly" %in% class(cp$ma_plotly(rank=TRUE)))
res <- cp$barplot_threshold()
stopifnot("ggplot" %in% class(res$FDR$plot))
stopifnot("ggplot" %in% class(cp$histogram_diff()))