expression() is an analysis function that computes the protein summary
statistics for a given tidyproteomics data object.
Usage
expression(
data = NULL,
...,
.pairs = NULL,
.method = stats::t.test,
.p.adjust = "BH"
)Arguments
- data
tidyproteomics data object
- ...
two sample comparison e.g. experimental/control
- .pairs
a list of vectors each containing two named sample groups
- .method
a two-distribution test function returning a p_value for the null hypothesis. Example functions include t.test, wilcox.test, stats::ks.test, additionally, the string "limma" can be used to select from the limma package to compute an empirical Bayesian estimation which performs better with non-linear distributions and uneven replicate balance between samples.
- .p.adjust
a stats::p.adjust string for multiple test correction, default is 'BH' (Benjamini & Hochberg, 1995)
Examples
library(dplyr, warn.conflicts = FALSE)
library(tidyproteomics)
# subset for fast execution
hela_proteins <- hela_proteins %>% subset(description %like% "Ribosome")
#>
#> ℹ Subsetting data: description %like% Ribosome
#> ✔ Subsetting data: description %like% Ribosome ... done
#>
# simple t.test expression analysis
hela_proteins %>%
expression(knockdown/control) %>%
export_analysis(knockdown/control, .analysis = "expression")
#> ℹ .. expression::t_test testing knockdown / control
#> ✔ .. expression::t_test testing knockdown / control [61ms]
#>
#> # A tibble: 17 × 27
#> protein imputed n imputed_control imputed_knockdown average_expression
#> <chr> <dbl> <int> <dbl> <dbl> <dbl>
#> 1 Q14146 0 6 0 0 667074.
#> 2 A0A087X020 0 6 0 0 27099436.
#> 3 Q969S9 0 6 0 0 1258238.
#> 4 Q9H7B2 0 6 0 0 17138597.
#> 5 Q15050 0 6 0 0 28765397.
#> 6 Q9Y221 0 6 0 0 9874377.
#> 7 O95478 0 6 0 0 10279868.
#> 8 Q14692 0 6 0 0 28455673.
#> 9 Q14137 0 6 0 0 41697675.
#> 10 Q8TDN6 0 6 0 0 35703697.
#> 11 Q9GZL7 0 6 0 0 46122319.
#> 12 Q9P2E9 0 6 0 0 146794734.
#> 13 Q9UHA3 0 6 0 0 2301616.
#> 14 Q96E11 0 6 0 0 11410697.
#> 15 Q9H9Y2 0 6 0 0 5643157.
#> 16 Q8N0V3 0.167 6 0.333 0 385672.
#> 17 Q9NZM5 0 6 0 0 2312164.
#> # ℹ 21 more variables: proportional_expression <dbl>, foldchange <dbl>,
#> # log2_foldchange <dbl>, p_value <dbl>, adj_p_value <dbl>,
#> # normalization <chr>, abundance_control_1 <dbl>, abundance_control_2 <dbl>,
#> # abundance_control_3 <dbl>, abundance_knockdown_1 <dbl>,
#> # abundance_knockdown_2 <dbl>, abundance_knockdown_3 <dbl>,
#> # description <chr>, biological_process <chr>, cellular_component <chr>,
#> # molecular_function <chr>, gene_id_entrez <chr>, gene_name <chr>, …
# a wilcox.test expression analysis
hela_proteins %>%
expression(knockdown/control, .method = stats::wilcox.test) %>%
export_analysis(knockdown/control, .analysis = "expression")
#> ℹ .. expression::wilcox_test testing knockdown / control
#> ✔ .. expression::wilcox_test testing knockdown / control [61ms]
#>
#> # A tibble: 17 × 27
#> protein imputed n imputed_control imputed_knockdown average_expression
#> <chr> <dbl> <int> <dbl> <dbl> <dbl>
#> 1 Q14146 0 6 0 0 667074.
#> 2 A0A087X020 0 6 0 0 27099436.
#> 3 Q969S9 0 6 0 0 1258238.
#> 4 Q9H7B2 0 6 0 0 17138597.
#> 5 Q15050 0 6 0 0 28765397.
#> 6 Q9Y221 0 6 0 0 9874377.
#> 7 O95478 0 6 0 0 10279868.
#> 8 Q14692 0 6 0 0 28455673.
#> 9 Q14137 0 6 0 0 41697675.
#> 10 Q8TDN6 0 6 0 0 35703697.
#> 11 Q9GZL7 0 6 0 0 46122319.
#> 12 Q9P2E9 0 6 0 0 146794734.
#> 13 Q9UHA3 0 6 0 0 2301616.
#> 14 Q96E11 0 6 0 0 11410697.
#> 15 Q9H9Y2 0 6 0 0 5643157.
#> 16 Q8N0V3 0.167 6 0.333 0 385672.
#> 17 Q9NZM5 0 6 0 0 2312164.
#> # ℹ 21 more variables: proportional_expression <dbl>, foldchange <dbl>,
#> # log2_foldchange <dbl>, p_value <dbl>, adj_p_value <dbl>,
#> # normalization <chr>, abundance_control_1 <dbl>, abundance_control_2 <dbl>,
#> # abundance_control_3 <dbl>, abundance_knockdown_1 <dbl>,
#> # abundance_knockdown_2 <dbl>, abundance_knockdown_3 <dbl>,
#> # description <chr>, biological_process <chr>, cellular_component <chr>,
#> # molecular_function <chr>, gene_id_entrez <chr>, gene_name <chr>, …
# a one-tailed wilcox.test expression analysis
wilcoxon_less <- function(x, y) {
stats::wilcox.test(x, y, alternative = "less")
}
hela_proteins <- hela_proteins %>%
expression(knockdown/control, .method = stats::wilcox.test)
#> ℹ .. expression::wilcox_test testing knockdown / control
#> ✔ .. expression::wilcox_test testing knockdown / control [61ms]
#>
hela_proteins %>% export_analysis(knockdown/control, .analysis = "expression")
#> # A tibble: 17 × 27
#> protein imputed n imputed_control imputed_knockdown average_expression
#> <chr> <dbl> <int> <dbl> <dbl> <dbl>
#> 1 Q14146 0 6 0 0 667074.
#> 2 A0A087X020 0 6 0 0 27099436.
#> 3 Q969S9 0 6 0 0 1258238.
#> 4 Q9H7B2 0 6 0 0 17138597.
#> 5 Q15050 0 6 0 0 28765397.
#> 6 Q9Y221 0 6 0 0 9874377.
#> 7 O95478 0 6 0 0 10279868.
#> 8 Q14692 0 6 0 0 28455673.
#> 9 Q14137 0 6 0 0 41697675.
#> 10 Q8TDN6 0 6 0 0 35703697.
#> 11 Q9GZL7 0 6 0 0 46122319.
#> 12 Q9P2E9 0 6 0 0 146794734.
#> 13 Q9UHA3 0 6 0 0 2301616.
#> 14 Q96E11 0 6 0 0 11410697.
#> 15 Q9H9Y2 0 6 0 0 5643157.
#> 16 Q8N0V3 0.167 6 0.333 0 385672.
#> 17 Q9NZM5 0 6 0 0 2312164.
#> # ℹ 21 more variables: proportional_expression <dbl>, foldchange <dbl>,
#> # log2_foldchange <dbl>, p_value <dbl>, adj_p_value <dbl>,
#> # normalization <chr>, abundance_control_1 <dbl>, abundance_control_2 <dbl>,
#> # abundance_control_3 <dbl>, abundance_knockdown_1 <dbl>,
#> # abundance_knockdown_2 <dbl>, abundance_knockdown_3 <dbl>,
#> # description <chr>, biological_process <chr>, cellular_component <chr>,
#> # molecular_function <chr>, gene_id_entrez <chr>, gene_name <chr>, …
# Note: the userdefined function is preserved in the operations tracking
hela_proteins %>% operations()
#> ℹ Data Transformations
#> • Data files (p97KD_HCT116_proteins.xlsx) were imported as proteins from
#> ProteomeDiscoverer
#> • Data subset `description` %like% `Ribosome`
#> • Analysis: expression difference wilcox_test knockdown/control, p.adjust =
#> BH
# limma expression analysis
hela_proteins %>%
expression(knockdown/control, .method = "limma") %>%
export_analysis(knockdown/control, .analysis = "expression")
#> ℹ .. expression::limma testing knockdown / control
#> ✔ .. expression::limma testing knockdown / control [129ms]
#>
#> # A tibble: 17 × 29
#> protein imputed n imputed_control imputed_knockdown log2_foldchange
#> <chr> <dbl> <int> <dbl> <dbl> <dbl>
#> 1 Q14146 0 6 0 0 1.52
#> 2 Q969S9 0 6 0 0 1.24
#> 3 A0A087X020 0 6 0 0 1.14
#> 4 Q9H7B2 0 6 0 0 0.778
#> 5 Q9Y221 0 6 0 0 0.608
#> 6 O95478 0 6 0 0 0.555
#> 7 Q15050 0 6 0 0 0.439
#> 8 Q14692 0 6 0 0 0.412
#> 9 Q8TDN6 0 6 0 0 0.385
#> 10 Q14137 0 6 0 0 0.380
#> 11 Q9P2E9 0 6 0 0 0.300
#> 12 Q9GZL7 0 6 0 0 0.260
#> 13 Q9UHA3 0 6 0 0 0.205
#> 14 Q96E11 0 6 0 0 0.167
#> 15 Q9H9Y2 0 6 0 0 0.130
#> 16 Q8N0V3 0.167 6 0.333 0 0.106
#> 17 Q9NZM5 0 6 0 0 -0.502
#> # ℹ 23 more variables: foldchange <dbl>, average_expression <dbl>,
#> # proportional_expression <dbl>, p_value <dbl>, adj_p_value <dbl>,
#> # limma_t_statistic <dbl>, limma_B_statistic <dbl>, normalization <chr>,
#> # abundance_control_1 <dbl>, abundance_control_2 <dbl>,
#> # abundance_control_3 <dbl>, abundance_knockdown_1 <dbl>,
#> # abundance_knockdown_2 <dbl>, abundance_knockdown_3 <dbl>,
#> # description <chr>, biological_process <chr>, cellular_component <chr>, …
# using the .pairs argument when multiple comparisons are needed
comps <- list(c("control","knockdown"),
c("knockdown","control"))
hela_proteins %>%
expression(.pairs = comps)
#> Expression Analysis - using the supplied 2 sample pairs ...
#> ℹ .. expression::t_test testing control / knockdown
#> ✔ .. expression::t_test testing control / knockdown [65ms]
#>
#> ℹ .. expression::t_test testing knockdown / control
#> ✔ .. expression::t_test testing knockdown / control [61ms]
#>
#>
#> ── Quantitative Proteomics Data Object ──
#>
#> Origin ProteomeDiscoverer
#> proteins (50.06 kB)
#> Composition 6 files
#> 2 samples (control, knockdown)
#> Quantitation 18 proteins
#> 2.7 log10 dynamic range
#> 7.41% missing values
#> *imputed
#> Accounting (4) num_peptides num_psms num_unique_peptides imputed
#> Annotations (9) description biological_process cellular_component molecular_function
#> gene_id_entrez gene_name wiki_pathway reactome_pathway
#> gene_id_ensemble
#> Analyses (2)
#> knockdown/control -> expression
#> control/knockdown -> expression
#>