PEIMAN2-vignette

3.1 Enrichment

runEnrichment() function can be used to run singular enrichment analysis for one list of protein. This function takes the following inputs:

• protein which is a character vector with protein UniProt accession codes.
• os.name which is a character vector of length one with exact taxonomy name of species.
• p.adj.method which is pvalue adjustment methos and optional. By default the value is set to ‘BH’. To see a possible list of values, type p.adjust.methods in R console.

As it was mentioned, the taxonomy name of species must be provided, e.g for a list of proteins belongs to human we pass os.name as ‘Homo sapiens (Human)’. The list is available at UniProt website. We also included a helper function named getTaxonomyName to help getting the exact taxonomy name. More on this function later.

The following lines of code illustrate the steps to run SEA on exmplData1. In runEnrichment function, we pass pl1 (a character vector of UniProt accession code) to perform SEA as follows and save the results in enrich1.

# Load PEIMAN2 package
library(PEIMAN2)
#> Loading required package: tidyverse
#> ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
#> ✔ dplyr     1.1.4     ✔ readr     2.1.5
#> ✔ forcats   1.0.0     ✔ stringr   1.5.1
#> ✔ ggplot2   3.5.1     ✔ tibble    3.2.1
#> ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
#> ✔ purrr     1.0.2     
#> ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
#> ✖ dplyr::filter() masks stats::filter()
#> ✖ dplyr::lag()    masks stats::lag()
#> ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

# Extract dataset and assign a variable name to it
pl1 <- exmplData1$pl1

# Run SEA on the list
enrich1 <- runEnrichment(protein = pl1, os.name = 'Homo sapiens (Human)')

The function returns a dataframe with the following columns:

• PTM: Post-translational modification (PTM).
• Freq in Uniprot: The total number of proteins with this PTM in UniProt.
• Freq in List: The total number of proteins with this PTM in the list.
• Sample: Number of proteins in the given list.
• Population: Total number of proteins in the current version of PEIMAN databse.
• pvalue: The p-value obtained from hypergeometric test (enrichment analysis).
• corrected pvalue: Adjusted p-value to correct for multiple testing.
• AC: Uniprot accession code (AC) of proteins with each PTM.

Note: As it was mentioned, the os.name is the exact taxonomy name of species that you are working with. The name should be exactly the same as UniProt definition. To facilitate searching for this name, you can pass your protein list with UniProt accession ID to getTaxonomyName function as follows. The result is the exact taxonomy name of protein list that you need to pass to runEnrichment. In the following example, the exact taxonomy name is printed:

getTaxonomyName(x = exmplData1$pl1)
#> [1] "Please use os.name = `Homo sapiens (Human)`"

Similarly, we can run SEA for the second list of proteins:

# Extract dataset and assign a variable name to it
pl2 <- exmplData1$pl2

# Run SEA on the list
enrich2 <- runEnrichment(protein = pl2, os.name = 'Homo sapiens (Human)')

3.2 Plotting SEA results

The plotEnrichment function can be used to visualize singular enrichment analysis for one set of proteins or match, analyse, and integrate results for two sets of proteins. To read more about this match and integration, please read details at (Nickchi, Jafari, and Kalantari 2015). We start by plotting the results for the firs list.

plotEnrichment(x = enrich1, sig.level = 0.05)

The results is a Lollipop plot which presents “Relative frequency” of each “PTM keywords” along with their corrected p-value measured in log scale. Note that only significant PTMs are shown. The default value for significance level is 5 percent. One can also visualize and match the results of two enrichment. For example, we can see the integrated results of enrich1 and enrich2 by the following line of code:

plotEnrichment(x = enrich1, y = enrich2, sig.level = 0.05)

The plot presents the ‘Relative frequency’ of common PTM terms among two enriched list (x and y). The coloring is the corrected p-value measured in log scale. By default a significance level of 5 percent is set to filter results. This can be modified by sig.level parameter.

4.1 PSEA

In order to run protein set enrichment analysis (PSEA), you can use runPSEA function. This function takes the following inputs:

• protein: A character vector with protein UniProt accession.
• os.name: A character of length one for the exact name of organism name.
• pexponent: Enrichment weighting exponent, p. The default value is 1. For values of p < 1, one can detect incoherent patterns in a set of protein. If one expects a small number of proteins to be coherent in a large set, then p > 1 is a good choice.
• nperm: Number of permutation to adjust for multiple testing in different pathways. Default is 1000.
• p.adj.method: The adjustment method to correct for multiple testing. Run p.adjust.methods to get a list of possible methods.
• sig.level: The significance level to filter pathways (applies to adjusted p-value), 0.05 is the default value.
• minSize: PTM pathways with a lower number of proteins than minSize are excluded. The default value is one.

psea_res <- runPSEA(protein = exmplData2, os.name = 'Rattus norvegicus (Rat)', nperm = 1000)

The result is a list with 6 elements. The first element of this list is important: A dataframe with protein set enrichment analysis (PSEA) results. Every row corresponds to a post-translational modification (PTM) pathway with the following columns:

• pval: p-value for singular enrichment analysis.
• pvaladj: Adjusted p-value
• ES: Enrichment score
• NES: Enrichmnt score normalized to mean enrichment of random samples of the same size.
• nMoreExtreme: Number of times the permuted sample resulted in profile with ES larger than abs(ES original)
• size: Number of proteins in the pathway
• Enrichment: Whether the proteins in the pathway have been enriched in the list.
• leadingEdge: UniProt accession code of leading edge proteins that drive the enrichment.

knitr::kable(psea_res[[1]], format = 'html')

4.2 Plotting

We now introduce the plotting features for protein set enrichment analysis. Two functions are included to visualize PSEA results returned from runPSEA function. The first plot is generated by plotPSEA function and shows Normalized Enrichment Score (NES) for each PTM pathway. User can restrict the number of pathways to draw based by adjusting sig.level parameter (default value is 0.05). The coloring of the plot indicates if the pathway is enriched or not.

plotPSEA(x = psea_res)

The second plot is generated by plotRunningScore function. A running enrichment score plot for each PTM can be plotted.