Gene Ontology (GO)¶
Pathways represent interconnected molecular networks of signaling cascades that govern critical cellular processes. They provide understandings cellular behavior mechanisms, insights of disease progression and treatment responses. In an R&D organization, managing pathways across different datasets are crucial for gaining insights of potential therapeutic targets and intervention strategies.
In this notebook we manage a pathway registry based on “2023 GO Biological Process” ontology. We’ll walk you through the steps of registering pathways and link them to genes.
In the following Standardize metadata on-the-fly notebook, we’ll demonstrate how to perform a pathway enrichment analysis and track the dataset with LaminDB.
# !pip install 'lamindb[jupyter,bionty]'
!lamin init --storage ./use-cases-registries --schema bionty
import lamindb as ln
import bionty as bt
import gseapy as gp
bt.settings.organism = "human" # globally set organism
→ connected lamindb: testuser1/use-cases-registries
Fetch GO pathways annotated with human genes using Enrichr¶
First we fetch the “GO_Biological_Process_2023” pathways for humans using GSEApy which wraps GSEA and Enrichr.
go_bp = gp.get_library(name="GO_Biological_Process_2023", organism="Human")
print(f"Number of pathways {len(go_bp)}")
Number of pathways 5406
go_bp["ATF6-mediated Unfolded Protein Response (GO:0036500)"]
['MBTPS1', 'MBTPS2', 'XBP1', 'ATF6B', 'DDIT3', 'CREBZF']
Parse out the ontology_id from keys, convert into the format of {ontology_id: (name, genes)}
def parse_ontology_id_from_keys(key):
"""Parse out the ontology id.
"ATF6-mediated Unfolded Protein Response (GO:0036500)" -> ("GO:0036500", "ATF6-mediated Unfolded Protein Response")
"""
name, id = key.rsplit(" (", 1)
id = id.rstrip(")")
return id, name
go_bp_parsed = {}
for key, genes in go_bp.items():
id, name = parse_ontology_id_from_keys(key)
go_bp_parsed[id] = (name, genes)
go_bp_parsed["GO:0036500"]
('ATF6-mediated Unfolded Protein Response',
['MBTPS1', 'MBTPS2', 'XBP1', 'ATF6B', 'DDIT3', 'CREBZF'])
Register pathway ontology in LaminDB¶
bionty = bt.Pathway.public()
bionty
Show code cell output
PublicOntology
Entity: Pathway
Organism: all
Source: go, 2024-06-17
#terms: 47856
Next, we register all the pathways and genes in LaminDB to finally link pathways to genes.
Register pathway terms¶
To register the pathways we make use of .from_values
to directly parse the annotated GO pathway ontology IDs into LaminDB.
pathway_records = bt.Pathway.from_values(go_bp_parsed.keys(), bt.Pathway.ontology_id)
ln.save(pathway_records)
Register gene symbols¶
Similarly, we use .from_values
for all Pathway associated genes to register them with LaminDB.
all_genes = bt.Gene.standardize(list({g for genes in go_bp.values() for g in genes}))
gene_records = bt.Gene.from_values(all_genes, bt.Gene.symbol)
ln.save(gene_records);
! found 56 synonyms in Bionty: ['RNF165', 'C11ORF65', 'C1ORF68', 'C1ORF43', 'SPATA5', 'HSPB11', 'DUSP13', 'C18ORF32', 'C2ORF69', 'C1ORF112', 'SPATA5L1', 'TDGF1', 'C9ORF78', 'C12ORF4', 'C3ORF70', 'FAM172A', 'TRB', 'C6ORF89', 'C8ORF88', 'TTC26', 'C9ORF72', 'PDZD3', 'FAM172BP', 'C12ORF29', 'C21ORF91', 'TTC30A', 'C1ORF56', 'C12ORF57', 'C3ORF38', 'SLC9A3R1', 'C18ORF25', 'C1ORF109', 'PIFO', 'ANP32C', 'C15ORF62', 'C10ORF90', 'C20ORF173', 'SLC9A3R2', 'C6ORF15', 'GBAP1', 'C8ORF17', 'C1ORF146', 'C11ORF80', 'C17ORF97', 'C12ORF50', 'SOGA1', 'C10ORF71', 'C17ORF75', 'C3ORF33', 'C19ORF12', 'C1ORF131', 'C2ORF49', 'THEG', 'ADGRF2', 'C17ORF99', 'USP41']
please add corresponding Gene records via `.from_values(['IFT70A', 'C1orf131', 'C2orf49', 'C17orf75', 'C11orf65', 'C1orf146', 'SPMAP2', 'NHERF1', 'C18orf32', 'ADGRF2P', 'C20orf173', 'C8orf17', 'AFG2B', 'ANP32CP', 'C12orf50', 'C3orf33', 'ARB2A', 'C19orf12', 'AIRIM', 'C1orf43', 'ARK2C', 'C9orf78', 'C6orf89', 'C12orf57', 'FIRRM', 'C15orf62', 'C6orf15', 'C3orf38', 'MTCL2', 'KPLCE', 'C1orf56', 'CIMAP3', 'C8orf88', 'ARB2BP', 'C21orf91', 'GBA1LP', 'C3orf70', 'NHERF2', 'TOP6BL', 'IFT56', 'AFG2A', 'C9orf72', 'LIAT1', 'ARK2N', 'C17orf99', 'NHERF4', 'C12orf4', 'IFT25', 'C10orf90', 'USP41P', 'DUSP13B', 'RLIG1', 'C10orf71', 'THRB', 'CRIPTO', 'C2orf69'])`
! ambiguous validation in Bionty for 1104 records: 'AGER', 'NPEPPSP1', 'YTHDC1', 'PCDHB4', 'CUTA', 'IFI27', 'KCND1', 'DDX39B', 'RAB11FIP3', 'RPP21', 'SLC12A7', 'GNL1', 'AIF1', 'GOLGA8Q', 'CCDC120', 'DCAF11', 'CASTOR3P', 'WDR73', 'BRD2', 'PCP4', ...
! did not create Gene records for 38 non-validated symbols: 'AFD1', 'AZF1', 'CCL3L1', 'CCL4L1', 'DGS2', 'DUX3', 'DUX5', 'FOXL3-OT1', 'IGL', 'LOC100653049', 'LOC102723475', 'LOC102723996', 'LOC102724159', 'LOC107984156', 'LOC112268384', 'LOC122319436', 'LOC122513141', 'LOC122539214', 'LOC344967', 'MDRV', ...
Manually register the 37 non-validated symbols:
inspect_result = bt.Gene.inspect(all_genes, bt.Gene.symbol)
nonval_genes = []
for g in inspect_result.non_validated:
nonval_genes.append(bt.Gene(symbol=g))
ln.save(nonval_genes)
! received 14696 unique terms, 1 empty/duplicated term is ignored
! 38 terms (0.30%) are not validated for symbol: MTRNR2L10, AZF1, LOC102723475, DGS2, DUX3, DUX5, CCL3L1, TAS2R36, TRA, MTRNR2L1, FOXL3-OT1, LOC122319436, MTRNR2L11, MTRNR2L7, LOC102724159, MTRNR2L3, MTRNR2L5, CCL4L1, MDRV, LOC100653049, ...
couldn't validate 38 terms: 'TAS2R36', 'MTRNR2L5', 'LOC344967', 'LOC102724159', 'MTRNR2L3', 'DUX5', 'MTRNR2L6', 'FOXL3-OT1', 'SEPTIN14P20', 'CCL3L1', 'MTRNR2L11', 'LOC122513141', 'DUX3', 'CCL4L1', 'MTRNR2L10', 'MTRNR2L12', 'MDRV', 'TRA', 'LOC100653049', 'IGL', ...
→ if you are sure, create new records via Gene() and save to your registry
Link pathway to genes¶
Now that we are tracking all pathways and genes records, we can link both of them to make the pathways even more queryable.
symbols_gene_records = {record.symbol: record for record in gene_records}
for pathway_record in pathway_records:
pathway_genes = go_bp_parsed.get(pathway_record.ontology_id)[1]
pathway_genes_records = [symbols_gene_records.get(gene) for gene in pathway_genes]
pathway_record.genes.set(pathway_genes_records)
Now genes are linked to pathways:
pathway_record.genes.list("symbol")
['XIAP', 'CST7', 'CARD18', 'CAST', 'CARD8']
pathway_record.genes.list("ensembl_gene_id")
['ENSG00000101966',
'ENSG00000077984',
'ENSG00000255501',
'ENSG00000153113',
'ENSG00000105483']
Move on to the next analysis: Standardize metadata on-the-fly