Gene Set Analysis (GSA) via self-contained test for group
effect on signaling (directed) pathways based on SEM. The core of the
methodology is implemented in the RICF algorithm of SEMrun(),
recovering from RICF output node-specific group effect p-values, and
Brown’s combined permutation p-values of node activation and inhibition.
SEMgsa(g = list(), data, group, method = "BH", alpha = 0.05, n_rep = 1000, ...)A list of pathways to be tested.
A matrix or data.frame. Rows correspond to subjects, and columns to graph nodes (variables).
A binary vector. This vector must be as long as the number of subjects. Each vector element must be 1 for cases and 0 for control subjects.
Multiple testing correction method. One of the values
available in p.adjust. By default, method is set
to "BH" (i.e., Benjamini-Hochberg correction).
Gene set test significance level (default = 0.05).
Number of randomization replicates (default = 1000).
Currently ignored.
A list of 2 objects:
"gsa", A data.frame reporting the following information for each pathway in the input list:
"No.nodes", pathway size (number of nodes);
"No.DEGs", number of differential espression genes (DEGs) within
the pathway, after multiple test correction with one of the methods
available in p.adjust;
"pert", pathway perturbation status (see details);
"pNA", Brown's combined P-value of pathway node activation;
"pNI", Brown's combined P-value of pathway node inhibition;
"PVAL", Bonferroni combined P-value of pNA, and pNI; i.e., 2* min(pNA, PNI);
"ADJP", Adjusted Bonferroni P-value of pathway perturbation; i.e., min(No.pathways * PVAL; 1).
"DEG", a list with DEGs names per pathways.
For gaining more biological insights into the functional roles of pre-defined subsets of genes, node perturbation obtained from RICF fitting has been combined with up- or down-regulation of genes from KEGG to obtain overall pathway perturbation as follows:
The node perturbation is defined as activated when the minimum among the p-values is positive; if negative, the status is inhibited.
Up- or down- regulation of genes (derived from KEGG database) has been obtained from the weighted adjacency matrix of each pathway as column sum of weights over each source node. If the overall sum of node weights is below 1, the pathway is flagged as down-regulated otherwise as up-regulated.
The combination between these two quantities allows to define the direction (up or down) of gene perturbation. Up- or down regulated gene status, associated with node inhibition, indicates a decrease in activation (or increase in inhibition) in cases with respect to control group. Conversely, up- or down regulated gene status, associated with node activation, indicates an increase in activation (or decrease in inhibition) in cases with respect to control group.
Grassi, M., Tarantino, B. SEMgsa: topology-based pathway enrichment analysis with structural equation models. BMC Bioinformatics 23, 344 (2022). https://doi.org/10.1186/s12859-022-04884-8
# \dontrun{
# Nonparanormal(npn) transformation
als.npn <- transformData(alsData$exprs)$data
#> Conducting the nonparanormal transformation via shrunkun ECDF...done.
# Selection of FTD-ALS pathways from kegg.pathways.Rdata
paths.name <- c("MAPK signaling pathway",
"Protein processing in endoplasmic reticulum",
"Endocytosis",
"Wnt signaling pathway",
"Neurotrophin signaling pathway",
"Amyotrophic lateral sclerosis")
j <- which(names(kegg.pathways) %in% paths.name)
GSA <- SEMgsa(kegg.pathways[j], als.npn, alsData$group,
method = "bonferroni", alpha = 0.05,
n_rep = 1000)
#> k = 1 Protein processing in endoplasmic reticulum
#> k = 2 MAPK signaling pathway
#> k = 3 Wnt signaling pathway
#> k = 4 Endocytosis
#> k = 5 Neurotrophin signaling pathway
#> k = 6 Amyotrophic lateral sclerosis
GSA$gsa
#> No.nodes No.DEGs pert
#> Neurotrophin signaling pathway 119 17 up act
#> Protein processing in endoplasmic reticulum 170 9 down act
#> Wnt signaling pathway 174 8 down inh
#> MAPK signaling pathway 301 47 up act
#> Amyotrophic lateral sclerosis 364 9 down inh
#> Endocytosis 250 2 <NA>
#> pNa pNi
#> Neurotrophin signaling pathway 2.353673e-14 9.359608e-07
#> Protein processing in endoplasmic reticulum 4.187547e-09 1.130873e-12
#> Wnt signaling pathway 4.591216e-12 2.752644e-01
#> MAPK signaling pathway 6.054046e-12 4.914468e-09
#> Amyotrophic lateral sclerosis 8.383960e-12 1.452882e-03
#> Endocytosis 1.564526e-01 3.142695e-04
#> PVAL ADJP
#> Neurotrophin signaling pathway 4.707346e-14 2.824407e-13
#> Protein processing in endoplasmic reticulum 2.261746e-12 1.357048e-11
#> Wnt signaling pathway 9.182433e-12 5.509460e-11
#> MAPK signaling pathway 1.210809e-11 7.264855e-11
#> Amyotrophic lateral sclerosis 1.676792e-11 1.006075e-10
#> Endocytosis 6.285389e-04 3.771234e-03
GSA$DEG
#> $`Neurotrophin signaling pathway`
#> [1] "207" "208" "1398" "1399" "3551" "3845" "4215" "4217" "5602" "5600"
#> [11] "5599" "5601" "4893" "5879" "5906" "9252" "6654"
#>
#> $`Protein processing in endoplasmic reticulum`
#> [1] "1649" "5602" "5599" "5601" "578" "2081" "4780" "23645" "7466"
#>
#> $`Wnt signaling pathway`
#> [1] "6885" "4772" "4775" "5532" "5534" "5567" "5879" "5881"
#>
#> $`MAPK signaling pathway`
#> [1] "207" "208" "774" "775" "776" "9254"
#> [7] "836" "1398" "1436" "1946" "8822" "2885"
#> [13] "3480" "8517" "8681" "4254" "4296" "7786"
#> [19] "4215" "9064" "8491" "1432" "5599" "23162"
#> [25] "5601" "2122" "8569" "4763" "4772" "4775"
#> [31] "56034" "5228" "100137049" "5494" "5532" "5534"
#> [37] "5567" "84867" "5879" "5881" "5906" "5921"
#> [43] "5922" "9252" "6654" "23118" "57551"
#>
#> $`Amyotrophic lateral sclerosis`
#> [1] "836" "5608" "5600" "4747" "4741" "5532" "5534" "5879" "2081"
#>
#> $Endocytosis
#> [1] "2261" "2264"
#>
# }