Infer Functional Cell Cell Communication

This notebook walks through infering cell-cell communication events that drive specific cell states

SpaceTravLR can model how a genetic perturbation in one cell modulates the gene expression in neighboring cells

%load_ext autoreload
%autoreload 2
%config InlineBackend.figure_format = 'retina'

import pandas as pd
import scanpy as sc
import sys 
import seaborn as sns
import matplotlib.pyplot as plt
sys.path.append('/Users/koush/Projects/jishnulab/SpaceTravLR/src')

from SpaceTravLR.beta import Betabase

adata = sc.read_h5ad('/Volumes/SSD/training_data/snrna_human_melanoma.h5ad')
adata

bdb = Betabase(
    adata=adata, 
    folder='/Volumes/SSD/lasso_runs/human_melanoma_v2', 
    auto_load=False
)

positive_interactions = bdb.collect_interactions(cell_type='CD8+ T', aggregate='positive')
negative_interactions = bdb.collect_interactions(cell_type='CD8+ T', aggregate='negative')

Unraveling genes in CD8+ T | mode: positive 100%|██████████| 2846/2846 [18:49<00:00, 2.52 parquet/s]

Unraveling genes in CD8+ T | mode: negative 100%|██████████| 2846/2846 [18:55<00:00, 2.51 parquet/s]

activation_genes = [
    'CD8A', 'GZMK', 'IFNG', 'PRF1', 
    'CD69', 'TBX21', 'GZMA', 'GZMB'
]

positive_interactions.sort_values(
    by='beta', 
    ascending=False, 
    key=abs).query('gene in @activation_genes')

	interaction	gene	beta	interaction_type
CD8+ T
12745	beta_THRB	GZMK	6.914885e+00	tf
303406	beta_ANGPT2$ITGB1	GZMA	5.198300e+00	ligand-receptor
281935	beta_IL7$IL2RG	GZMA	3.815622e+00	ligand-receptor
20738	beta_GRHL2	GZMA	3.753025e+00	tf
12678	beta_THRB	CD8A	3.655421e+00	tf
...	...	...	...	...
265644	beta_SEMA4A$PLXNA3	TBX21	2.035810e-08	ligand-receptor
272200	beta_THY1$ITGAV	CD8A	1.479332e-08	ligand-receptor
436033	beta_COL1A2$ITGA9	IFNG	8.826670e-09	ligand-receptor
167988	beta_COL6A2$ITGA10	GZMA	5.215077e-09	ligand-receptor
21963	beta_IFNG#GLIS3	TBX21	6.802419e-10	ligand-tf

2738 rows × 4 columns

df = positive_interactions.sort_values(
    by='beta', 
    ascending=False, 
    key=abs).query(
        'gene in @activation_genes and interaction_type == "ligand-receptor"')

cd8_means = adata.to_df(layer='imputed_count')[adata.obs['cell_type'] == "CD8+ T"].mean()
genes = df.interaction.str.split("$").str.get(1)
receptor_expr = genes.map(cd8_means)

df['weighted_beta'] = df['beta'] * receptor_expr

df_cd8 = df.copy()
df_cd8.interaction = df_cd8.interaction.str.replace('beta_', '')
df_cd8.interaction = df_cd8.interaction.str.replace('$', '-')

top_interactions = (
    df_cd8[['interaction', 'weighted_beta']].groupby(
        'interaction').mean().sort_values(by='weighted_beta', ascending=False)
      .head(15)
)

sns.set_theme(style="white", context="talk")
plt.figure(figsize=(9, 6))

ax = sns.barplot(
    data=top_interactions,
    x="weighted_beta",
    y="interaction",
    hue="interaction",      
    palette="Greens_r",
    linewidth=2,
    edgecolor='black',         
    legend=False
)

ax.set_title("Top Interactions\nPromoting CD8 Activation", fontweight='bold', pad=20)
ax.set_xlabel("Weighted Beta (Interaction Strength)", fontsize=14)
ax.set_ylabel("")

sns.despine(ax=ax, left=True, bottom=True)

plt.tight_layout()
plt.show()