Integrating Multi-Modal Spatially-Resolved Technologies with LIANA+¶
Here, we apply of LIANA+ on a spatially-resolved metabolite-transcriptome dataset from a recent murine Parkinson’s disease model Vicari et al., 2023. We demonstrate LIANA+’s utility in harmonizing spatially-resolved transcriptomics and MALDI-MSI data to unravel metabolite-mediated interactions and the molecular mechanisms of dopamine regulation in the striatum.
Two particular challenges with this data are:
The unaligned spatial locations of the two omics technologies
The untargeted nature of the MALDI-MSI data, which results in a large number of features with unknown identities. Only few of which were previously identified as specific metabolites.
Here, we show untargeted modelling of known and unknown metabolite peaks and their spatial relationships with transcriptomics data. Specifically, we use a multi-view modelling strategy MISTy to decipher global spatial relationships of metabolite peaks with cell types and brain-specific receptors. Then, we use LIANA+’s local metrics to pinpoint the subregions of interaction.
We also show strategies to enable spatial multi-omics analysis from diverse omics technologies with unaligned locations and observations.
[ ]:
import numpy as np
import liana as li
import mudata as mu
import scanpy as sc
from matplotlib import pyplot as plt
from adjustText import adjust_text
[2]:
# set global figure parameters
kwargs = {'frameon':False, 'size':1.5, 'img_key':'lowres'}
Obtain and Examine the Data¶
First, we obtain the a single slide from the dataset. This slide has already been preprocessed, including filtering and normalisation. We have log1p transformed the RNA-seq data and total-ion count normalised the MALDI-MSI data.
We have also pre-aligned the images from the two technologies, though the observations are not aligned - an issue that we will address in this notebook. For image and coordinate transformations, we refer the users to SpatialData or stAlign.
We have additionally deconvoluted the cell types in the RNA-seq data using Tangram.
So, in total we have three modalities: the MALDI-MSI data, the RNA-seq data, and the cell type data:
[3]:
# let's download the data
rna = sc.read("sma_rna.h5ad", backup_url="https://figshare.com/ndownloader/files/44624974?private_link=4744950f8768d5c8f68c")
msi = sc.read("sma_msi.h5ad", backup_url="https://figshare.com/ndownloader/files/44624971?private_link=4744950f8768d5c8f68c")
ct = sc.read("sma_ct.h5ad", backup_url="https://figshare.com/ndownloader/files/44624968?private_link=4744950f8768d5c8f68c")
[ ]:
# and create a MuData object
mdata = mu.MuData({'rna':rna, 'msi':msi, 'ct':ct})
mdata
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/mudata/_core/mudata.py:1531: FutureWarning: From 0.4 .update() will not pull obs/var columns from individual modalities by default anymore. Set mudata.set_options(pull_on_update=False) to adopt the new behaviour, which will become the default. Use new pull_obs/pull_var and push_obs/push_var methods for more flexibility.
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/mudata/_core/mudata.py:1429: FutureWarning: From 0.4 .update() will not pull obs/var columns from individual modalities by default anymore. Set mudata.set_options(pull_on_update=False) to adopt the new behaviour, which will become the default. Use new pull_obs/pull_var and push_obs/push_var methods for more flexibility.
MuData object with n_obs × n_vars = 6041 × 17782
3 modalities
rna: 3036 x 16486
obs: 'in_tissue', 'array_row', 'array_col', 'x', 'y', 'lesion', 'region', 'n_genes_by_counts', 'log1p_n_genes_by_counts', 'total_counts', 'log1p_total_counts', 'pct_counts_in_top_50_genes', 'pct_counts_in_top_100_genes', 'pct_counts_in_top_200_genes', 'pct_counts_in_top_500_genes', 'total_counts_mt', 'log1p_total_counts_mt', 'pct_counts_mt', 'n_genes', 'n_counts'
var: 'gene_ids', 'feature_types', 'genome', 'mt', 'n_cells_by_counts', 'mean_counts', 'log1p_mean_counts', 'pct_dropout_by_counts', 'total_counts', 'log1p_total_counts', 'n_cells'
uns: 'lesion_colors', 'log1p', 'region_colors', 'spatial'
obsm: 'spatial'
layers: 'counts'
msi: 3005 x 1248
obs: 'x', 'y', 'array_row', 'array_col', 'leiden', 'n_counts', 'index_right', 'region', 'lesion'
var: 'mean', 'std', 'mz', 'max_intensity', 'mz_raw', 'annotated'
uns: 'leiden', 'leiden_colors', 'log1p', 'neighbors', 'pca', 'spatial'
obsm: 'X_pca', 'spatial'
varm: 'PCs'
layers: 'raw'
obsp: 'connectivities', 'distances'
ct: 3036 x 48
obs: 'in_tissue', 'array_row', 'array_col', 'x', 'y', 'lesion', 'region', 'n_genes_by_counts', 'log1p_n_genes_by_counts', 'total_counts', 'log1p_total_counts', 'pct_counts_in_top_50_genes', 'pct_counts_in_top_100_genes', 'pct_counts_in_top_200_genes', 'pct_counts_in_top_500_genes', 'total_counts_mt', 'log1p_total_counts_mt', 'pct_counts_mt', 'n_genes', 'n_counts', 'uniform_density', 'rna_count_based_density'
uns: 'lesion_colors', 'log1p', 'overlap_genes', 'region_colors', 'spatial', 'training_genes'
obsm: 'spatial', 'tangram_ct_pred'[5]:
fig, axes = plt.subplots(1, 3, figsize=(12, 3))
sc.pl.spatial(rna, color='log1p_total_counts', ax=axes[0], **kwargs, show=False)
sc.pl.spatial(msi, color='Dopamine', cmap='magma', ax=axes[1], **kwargs, show=False)
sc.pl.spatial(ct, color='MSN1', cmap='viridis', ax=axes[2], **kwargs, show=False)
fig.subplots_adjust(wspace=0, hspace=0)
fig.tight_layout()
If you look closely here, you will notice that the metabolite locations are not aligned in a grid-like manner as the 10X Visium Data. We will address this in the next steps.
Experimental Design¶
Of note, mice in this data were subjected to unilateral 6-hydroxydopamine-induced lesions in one hemisphere while the other remained intact. These 6-Hydroxydopamine-induced lesions selectively destroy substantia nigra-originated dopaminergic neurons, thereby impairing dopamine-mediated regulatory mechanisms of the striatum - an area of the brain crucial for movement coordination.
Along with annotations of the lesioned and intact hemispheres, we also have annotations for the striatum:
[6]:
sc.pl.spatial(rna, color=['lesion', 'region'], **kwargs, wspace=0.25)
Multi-view Modelling Metabolite Intensities¶
Next, we will model the metabolite intensities jointly using cell types and brain specific receptors as predictors.
Remove features with little-to-no variation¶
Done to avoid fitting noisy readouts in our model. Ideally, one could instead use e.g. Moran’s I to identify spatially-variable features.
[8]:
sc.pp.highly_variable_genes(rna, flavor='cell_ranger', n_top_genes=5000)
sc.pp.highly_variable_genes(msi, flavor='cell_ranger', n_top_genes=150)
ct.var['cv'] = ct.X.toarray().var(axis=0) / ct.X.toarray().mean(axis=0)
ct.var['highly_variable'] = ct.var['cv'] > np.percentile(ct.var['cv'], 20)
[9]:
msi = msi[:, msi.var['highly_variable']]
rna = rna[:, rna.var['highly_variable']]
ct = ct[:, ct.var['highly_variable']]
Note that we use simple coefficient of variation and highly-variable gene functions, but one may easily replace this step with e.g. the spatially-informed Moran’s I from Squidpy.
Additional processing steps¶
Scale the intensities and cap the max value
[10]:
sc.pp.scale(msi, max_value=5)
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/scanpy/preprocessing/_scale.py:318: UserWarning: Received a view of an AnnData. Making a copy.
Obtain brain-specific metabolite-receptor interactions from MetalinksDB.
[11]:
metalinks = li.rs.get_metalinks(tissue_location='Brain',
biospecimen_location='Cerebrospinal Fluid (CSF)',
source=['CellPhoneDB', 'NeuronChat']
)
metalinks.head()
Downloading database...
Database downloaded and saved to /data/ddimitrov/repos/liana-py/docs/source/notebooks/metalinksdb.db.
[11]:
| hmdb | uniprot | gene_symbol | metabolite | mor | transport_direction | type | source | |
|---|---|---|---|---|---|---|---|---|
| 0 | HMDB0000870 | P25021 | HRH2 | Histamine | 0 | None | lr | CellPhoneDB |
| 1 | HMDB0000870 | P25021 | HRH2 | Histamine | 1 | None | lr | CellPhoneDB |
| 2 | HMDB0000870 | P35367 | HRH1 | Histamine | 0 | None | lr | CellPhoneDB |
| 3 | HMDB0000870 | P35367 | HRH1 | Histamine | 1 | None | lr | CellPhoneDB |
| 4 | HMDB0000870 | Q9H3N8 | HRH4 | Histamine | 0 | None | lr | CellPhoneDB |
Convert to murine symbols using orthology knowledge from HCOP. If you use this function, please reference the original database.
[12]:
map_df = li.rs.get_hcop_orthologs(columns=['human_symbol', 'mouse_symbol'],
min_evidence=3
).rename(columns={'human_symbol':'source',
'mouse_symbol':'target'})
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/resource/_orthology.py:199: DtypeWarning: Columns (0) have mixed types. Specify dtype option on import or set low_memory=False.
Translate the Receptors to Murine symbols
[13]:
metalinks = li.rs.translate_column(resource=metalinks,
map_df=map_df,
column='gene_symbol',
one_to_many=1)
metalinks.head()
[13]:
| hmdb | uniprot | gene_symbol | metabolite | mor | transport_direction | type | source | |
|---|---|---|---|---|---|---|---|---|
| 0 | HMDB0000870 | P25021 | Hrh2 | Histamine | 0 | None | lr | CellPhoneDB |
| 1 | HMDB0000870 | P25021 | Hrh2 | Histamine | 1 | None | lr | CellPhoneDB |
| 2 | HMDB0000870 | P35367 | Hrh1 | Histamine | 0 | None | lr | CellPhoneDB |
| 3 | HMDB0000870 | P35367 | Hrh1 | Histamine | 1 | None | lr | CellPhoneDB |
| 4 | HMDB0000870 | Q9H3N8 | Hrh4 | Histamine | 0 | None | lr | CellPhoneDB |
[14]:
metalinks
[14]:
| hmdb | uniprot | gene_symbol | metabolite | mor | transport_direction | type | source | |
|---|---|---|---|---|---|---|---|---|
| 0 | HMDB0000870 | P25021 | Hrh2 | Histamine | 0 | None | lr | CellPhoneDB |
| 1 | HMDB0000870 | P25021 | Hrh2 | Histamine | 1 | None | lr | CellPhoneDB |
| 2 | HMDB0000870 | P35367 | Hrh1 | Histamine | 0 | None | lr | CellPhoneDB |
| 3 | HMDB0000870 | P35367 | Hrh1 | Histamine | 1 | None | lr | CellPhoneDB |
| 4 | HMDB0000870 | Q9H3N8 | Hrh4 | Histamine | 0 | None | lr | CellPhoneDB |
| ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 259 | HMDB0000077 | O75469 | Nr1i2 | Dehydroepiandrosterone | 1 | None | lr | CellPhoneDB |
| 260 | HMDB0000077 | P03372 | Esr1 | Dehydroepiandrosterone | 0 | None | lr | CellPhoneDB |
| 261 | HMDB0000077 | P03372 | Esr1 | Dehydroepiandrosterone | 1 | None | lr | CellPhoneDB |
| 262 | HMDB0000077 | Q07869 | Ppara | Dehydroepiandrosterone | 1 | None | lr | CellPhoneDB |
| 263 | HMDB0000181 | P51810 | Gpr143 | L-Dopa | 1 | None | lr | CellPhoneDB |
254 rows × 8 columns
Intersect the RNA modality with the receptors
[15]:
receptors = np.intersect1d(metalinks['gene_symbol'].unique(), rna.var_names)
rec = rna[:, receptors].copy()
Compute Spatial Proximies for the Multi-view Model¶
We use the metabolite modality as reference to which we align the other modalities. We use the spatial_neighbors function to compute the spatial proximity from cell types and brain-specific receptors to the metabolite intensities.
[16]:
plot, _ = li.ut.query_bandwidth(coordinates=rna.obsm['spatial'], start=0, end=1500)
plot
We then pick 1,000 as it’s roughly equivalent to the 6 nearest neighbors from a hexagonal grid.
[17]:
bandwidth = 500
cutoff = 0.1
# distances of metabolties to RNA
reference = mdata.mod["msi"].obsm["spatial"]
Compute distances for each extra modality
[18]:
li.ut.spatial_neighbors(ct, bandwidth=bandwidth, cutoff=cutoff, spatial_key="spatial", reference=reference, set_diag=False, standardize=False)
li.ut.spatial_neighbors(rec, bandwidth=bandwidth, cutoff=cutoff, spatial_key="spatial", reference=reference, set_diag=False, standardize=False)
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/utils/spatial_neighbors.py:141: ImplicitModificationWarning: Setting element `.obsm['spatial_connectivities']` of view, initializing view as actual.
Construct and Run the Multi-view model¶
Specifically, we use the metabolite intensities as the intraview - i.e. targets of prediction. While the cell types and brain-specific receptors are the extraview - i.e. spatially-weighted predictors.
[19]:
# MISTy
mdata.update_obs()
misty = li.mt.MistyData({"intra": msi, "receptor": rec, "ct": ct}, enforce_obs=False, obs=mdata.obs)
misty
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/mudata/_core/mudata.py:1429: FutureWarning: From 0.4 .update() will not pull obs/var columns from individual modalities by default anymore. Set mudata.set_options(pull_on_update=False) to adopt the new behaviour, which will become the default. Use new pull_obs/pull_var and push_obs/push_var methods for more flexibility.
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/mudata/_core/mudata.py:1531: FutureWarning: From 0.4 .update() will not pull obs/var columns from individual modalities by default anymore. Set mudata.set_options(pull_on_update=False) to adopt the new behaviour, which will become the default. Use new pull_obs/pull_var and push_obs/push_var methods for more flexibility.
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/mudata/_core/mudata.py:1429: FutureWarning: From 0.4 .update() will not pull obs/var columns from individual modalities by default anymore. Set mudata.set_options(pull_on_update=False) to adopt the new behaviour, which will become the default. Use new pull_obs/pull_var and push_obs/push_var methods for more flexibility.
[19]:
MuData object with n_obs × n_vars = 6041 × 207
obs: 'rna:in_tissue', 'rna:array_row', 'rna:array_col', 'rna:x', 'rna:y', 'rna:lesion', 'rna:region', 'rna:n_genes_by_counts', 'rna:log1p_n_genes_by_counts', 'rna:total_counts', 'rna:log1p_total_counts', 'rna:pct_counts_in_top_50_genes', 'rna:pct_counts_in_top_100_genes', 'rna:pct_counts_in_top_200_genes', 'rna:pct_counts_in_top_500_genes', 'rna:total_counts_mt', 'rna:log1p_total_counts_mt', 'rna:pct_counts_mt', 'rna:n_genes', 'rna:n_counts', 'msi:x', 'msi:y', 'msi:array_row', 'msi:array_col', 'msi:leiden', 'msi:n_counts', 'msi:index_right', 'msi:region', 'msi:lesion'
var: 'highly_variable'
3 modalities
intra: 3005 x 150
obs: 'x', 'y', 'array_row', 'array_col', 'leiden', 'n_counts', 'index_right', 'region', 'lesion'
var: 'mean', 'std', 'mz', 'max_intensity', 'mz_raw', 'annotated', 'highly_variable', 'means', 'dispersions', 'dispersions_norm'
uns: 'leiden', 'leiden_colors', 'log1p', 'neighbors', 'pca', 'spatial', 'hvg'
obsm: 'X_pca', 'spatial'
varm: 'PCs'
layers: 'raw'
obsp: 'connectivities', 'distances'
receptor: 3036 x 19
obs: 'in_tissue', 'array_row', 'array_col', 'x', 'y', 'lesion', 'region', 'n_genes_by_counts', 'log1p_n_genes_by_counts', 'total_counts', 'log1p_total_counts', 'pct_counts_in_top_50_genes', 'pct_counts_in_top_100_genes', 'pct_counts_in_top_200_genes', 'pct_counts_in_top_500_genes', 'total_counts_mt', 'log1p_total_counts_mt', 'pct_counts_mt', 'n_genes', 'n_counts'
var: 'gene_ids', 'feature_types', 'genome', 'mt', 'n_cells_by_counts', 'mean_counts', 'log1p_mean_counts', 'pct_dropout_by_counts', 'total_counts', 'log1p_total_counts', 'n_cells', 'highly_variable', 'means', 'dispersions', 'dispersions_norm'
uns: 'lesion_colors', 'log1p', 'region_colors', 'spatial', 'hvg'
obsm: 'spatial', 'spatial_connectivities'
varm: 'weighted'
layers: 'counts'
ct: 3036 x 38
obs: 'in_tissue', 'array_row', 'array_col', 'x', 'y', 'lesion', 'region', 'n_genes_by_counts', 'log1p_n_genes_by_counts', 'total_counts', 'log1p_total_counts', 'pct_counts_in_top_50_genes', 'pct_counts_in_top_100_genes', 'pct_counts_in_top_200_genes', 'pct_counts_in_top_500_genes', 'total_counts_mt', 'log1p_total_counts_mt', 'pct_counts_mt', 'n_genes', 'n_counts', 'uniform_density', 'rna_count_based_density'
var: 'cv', 'highly_variable'
uns: 'lesion_colors', 'log1p', 'overlap_genes', 'region_colors', 'spatial', 'training_genes'
obsm: 'spatial', 'tangram_ct_pred', 'spatial_connectivities'
varm: 'weighted'By using the metabolite modality as a reference, we can now predict the metabolite intensities from the spatially-weighted cell types and brain-specific receptors.
Specifically, we do so, by modelling the lesioned and intact hemispheres separately - enabled via the maskby parameter. This masking procedure, simply masks the observations according to the lesion annotations, following spatially-weighting the extra views.
[20]:
misty(model=li.mt.sp.LinearModel, verbose=True, bypass_intra=True, maskby='lesion')
Now learning: 313.17 masked by intact: 3%|▎ | 5/150 [00:05<01:27, 1.65it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 313.17 masked by lesioned: 3%|▎ | 5/150 [00:05<01:27, 1.65it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 381.17 masked by intact: 3%|▎ | 5/150 [00:05<01:27, 1.65it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 3-MT masked by lesioned: 11%|█ | 16/150 [00:07<00:30, 4.47it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 554.26 masked by intact: 15%|█▍ | 22/150 [00:08<00:26, 4.81it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 556.05 masked by intact: 17%|█▋ | 25/150 [00:09<00:24, 5.17it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 582.2 masked by lesioned: 21%|██ | 31/150 [00:10<00:25, 4.69it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 673.25 masked by lesioned: 23%|██▎ | 35/150 [00:11<00:25, 4.45it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 681.38 masked by intact: 25%|██▌ | 38/150 [00:11<00:27, 4.13it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 700.5 masked by intact: 30%|███ | 45/150 [00:13<00:22, 4.64it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 700.5 masked by lesioned: 30%|███ | 45/150 [00:13<00:22, 4.64it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 755.16 masked by intact: 35%|███▍ | 52/150 [00:15<00:21, 4.57it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 757.33 masked by intact: 35%|███▌ | 53/150 [00:15<00:21, 4.53it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 757.33 masked by lesioned: 35%|███▌ | 53/150 [00:15<00:21, 4.53it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 773.53 masked by intact: 36%|███▌ | 54/150 [00:15<00:21, 4.37it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 785.45 masked by lesioned: 39%|███▉ | 59/150 [00:16<00:17, 5.18it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 789.31 masked by lesioned: 40%|████ | 60/150 [00:16<00:16, 5.44it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 791.99 masked by lesioned: 41%|████▏ | 62/150 [00:17<00:16, 5.22it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 807.32 masked by intact: 45%|████▍ | 67/150 [00:18<00:17, 4.69it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 808.31 masked by lesioned: 47%|████▋ | 70/150 [00:18<00:15, 5.07it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 811.33 masked by lesioned: 47%|████▋ | 71/150 [00:18<00:15, 5.08it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 819.36 masked by intact: 49%|████▉ | 74/150 [00:19<00:16, 4.54it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 820.33 masked by intact: 50%|█████ | 75/150 [00:19<00:15, 4.79it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 824.29 masked by intact: 53%|█████▎ | 79/150 [00:20<00:14, 5.02it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 825.35 masked by intact: 55%|█████▌ | 83/150 [00:21<00:13, 5.03it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 826.28 masked by lesioned: 56%|█████▌ | 84/150 [00:21<00:12, 5.33it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 836.33 masked by lesioned: 58%|█████▊ | 87/150 [00:22<00:13, 4.77it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 838.3 masked by intact: 59%|█████▊ | 88/150 [00:22<00:12, 4.91it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 838.3 masked by lesioned: 59%|█████▊ | 88/150 [00:22<00:12, 4.91it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 839.16 masked by lesioned: 59%|█████▉ | 89/150 [00:22<00:11, 5.22it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 841.15 masked by lesioned: 61%|██████ | 91/150 [00:22<00:11, 4.93it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 844.52 masked by intact: 66%|██████▌ | 99/150 [00:24<00:09, 5.41it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 846.54 masked by intact: 67%|██████▋ | 100/150 [00:24<00:09, 5.34it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 846.54 masked by lesioned: 67%|██████▋ | 100/150 [00:24<00:09, 5.34it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 932.25 masked by intact: 74%|███████▍ | 111/150 [00:26<00:08, 4.75it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 932.25 masked by lesioned: 74%|███████▍ | 111/150 [00:26<00:08, 4.75it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 938.24 masked by lesioned: 75%|███████▌ | 113/150 [00:27<00:07, 4.79it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 957.23 masked by intact: 80%|████████ | 120/150 [00:28<00:07, 4.02it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 957.23 masked by lesioned: 80%|████████ | 120/150 [00:28<00:07, 4.02it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 958.24 masked by lesioned: 81%|████████ | 121/150 [00:29<00:06, 4.24it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 966.27 masked by lesioned: 81%|████████▏ | 122/150 [00:29<00:06, 4.13it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 970.24 masked by lesioned: 82%|████████▏ | 123/150 [00:29<00:06, 4.33it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 984.25 masked by intact: 84%|████████▍ | 126/150 [00:30<00:05, 4.72it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 984.25 masked by lesioned: 84%|████████▍ | 126/150 [00:30<00:05, 4.72it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 985.34 masked by lesioned: 85%|████████▍ | 127/150 [00:30<00:05, 4.47it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 1001.31 masked by intact: 87%|████████▋ | 131/150 [00:31<00:04, 4.26it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 1024.21 masked by intact: 91%|█████████▏| 137/150 [00:32<00:02, 4.34it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 1025.22 masked by intact: 92%|█████████▏| 138/150 [00:32<00:02, 4.60it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 1025.22 masked by lesioned: 92%|█████████▏| 138/150 [00:33<00:02, 4.60it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 1026.22 masked by lesioned: 93%|█████████▎| 139/150 [00:33<00:02, 4.83it/s]/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 1042.29 masked by intact: 95%|█████████▌| 143/150 [00:33<00:01, 4.52it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 1046.34 masked by intact: 97%|█████████▋| 145/150 [00:34<00:01, 4.60it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 1048.34 masked by intact: 99%|█████████▊| 148/150 [00:35<00:00, 4.84it/s] /data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/liana/method/sp/_misty/_Misty.py:318: RuntimeWarning: invalid value encountered in divide
Now learning: 1048.38 masked by lesioned: 100%|██████████| 150/150 [00:35<00:00, 4.23it/s]
We can see that in the intact hemisphere there are several metabolites peaks which are relatively well predicted (R2 > 0.5), and among them are Dopamine and 3-Methoxytyramine (3-MT). Of note, while we don’t focus on the unannotated metabolite peaks, we can see that some of them are also well predicted. This is a good indication of which metabolite peaks might be of interest for further investigation, and subsequently, identification.
[21]:
li.pl.target_metrics(misty, stat='multi_R2', return_fig=True, top_n=20, filter_fun=lambda x: x['intra_group']=='intact')
On the other hand, we don’t see those peaks in the lesioned hemisphere, which is consistent with the 6-Hydroxydopamine-induced lesions, and the absence of Dopamine in this hemisphere.
[22]:
li.pl.target_metrics(misty, stat='multi_R2', return_fig=True, top_n=20, filter_fun=lambda x: x['intra_group']=='lesioned')
Within the intact hemisphere we can see that in the top predictors of Dopamine are MSN1/2 as well as Drd1/2 receptors. MSN1/2 are Medium Spiny Neurons, which are the main cell type in the striatum, and Drd1/2 are Dopamine receptors. This is consistent with the known biology of the striatum, where Dopamine is a key neurotransmitter, and the Drd1/2 receptors are the main receptors for Dopamine in the striatum.
[24]:
interactions = misty.uns['interactions']
[25]:
interactions = interactions[(interactions['intra_group'] == 'intact') & (interactions['target'] == 'Dopamine')]
# Create scatter plot
plt.figure(figsize=(5, 4))
# rank rank by abs importances
interactions['rank'] = interactions['importances'].rank(ascending=False)
plt.scatter(interactions['rank'], interactions['importances'], s=11,
c=interactions['view'].map({'ct': '#008B8B', 'receptor': '#a11838'}))
# add for top 10
top_n = interactions[interactions['rank'] <= 10]
texts = []
for i, row in top_n.iterrows():
texts.append(plt.text(row['rank'], row['importances'], row['predictor'], fontsize=10))
adjust_text(texts, arrowprops=dict(arrowstyle="->", color='grey', lw=1.5))
plt.tight_layout()
/tmp/ipykernel_2555126/360477328.py:5: SettingWithCopyWarning:
A value is trying to be set on a copy of a slice from a DataFrame.
Try using .loc[row_indexer,col_indexer] = value instead
See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy
Identifying Local Interactions¶
Focusing on Dopamine, we can next use LIANA+’s local metrics to identify the subregions of interactions with MSN1/2 cells.
While the transformation of spatial locations is done internally by MISTy, we need to transform the metabolite intensities to a grid-like manner, so that we can use the local metrics. To do so, we use the interpolate_adata function, which interpolates one modality to another. Here, we will interpolate the metabolite intensities to the RNA-seq data, so that we can use the spatial locations of the RNA-seq data to identify the local interactions.
[26]:
metabs = li.ut.interpolate_adata(target=msi, reference=rna, use_raw=False, spatial_key='spatial')
Notice that the Metabolite observations now resemble the grid-like structure of Visium data:
[27]:
sc.set_figure_params(dpi=80, dpi_save=300, format='png', frameon=False, transparent=True, figsize=[5,5])
[28]:
sc.pl.spatial(metabs, color='Dopamine', cmap='magma', **kwargs)
Let’s rebuild a MuData object with these updated metabolite intensities
[29]:
mdata = mu.MuData({'msi': metabs, 'rna':rna, 'deconv':ct}, obsm=rna.obsm, obs=rna.obs, uns=rna.uns)
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/mudata/_core/mudata.py:1531: FutureWarning: From 0.4 .update() will not pull obs/var columns from individual modalities by default anymore. Set mudata.set_options(pull_on_update=False) to adopt the new behaviour, which will become the default. Use new pull_obs/pull_var and push_obs/push_var methods for more flexibility.
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/mudata/_core/mudata.py:1429: FutureWarning: From 0.4 .update() will not pull obs/var columns from individual modalities by default anymore. Set mudata.set_options(pull_on_update=False) to adopt the new behaviour, which will become the default. Use new pull_obs/pull_var and push_obs/push_var methods for more flexibility.
and now we can again calculate the spatial proximities, but this time without a reference as all observations within the MuData have the same spatial locations.
[30]:
li.ut.spatial_neighbors(mdata, bandwidth=bandwidth, cutoff=cutoff, set_diag=True)
Define interactions of interest:
[31]:
interactions = metalinks[['metabolite', 'gene_symbol']].apply(tuple, axis=1).tolist()
Let’s calculate the local metrics for the Dopamine intensities with and Drd1/2 receptors.
[32]:
lrdata = li.mt.bivariate(mdata,
local_name='cosine',
x_mod='msi',
y_mod='rna',
x_use_raw=False,
y_use_raw=False,
verbose=True,
mask_negatives=True,
n_perms=1000,
interactions=interactions,
x_transform=sc.pp.scale,
y_transform=sc.pp.scale,
)
Using `.X`!
Using `.X`!
Using `.X`!
Converting to sparse csr matrix!
Transforming msi using scale
Transforming rna using scale
/data/ddimitrov/software/miniforge3/envs/kiara2/lib/python3.10/site-packages/anndata/_core/anndata.py:430: FutureWarning: The dtype argument is deprecated and will be removed in late 2024.
Using provided `interactions`.
100%|██████████| 1000/1000 [00:01<00:00, 508.10it/s]
Here, we will plot the permutation-based local P-values
[33]:
sc.pl.spatial(lrdata,
color=['Dopamine^Drd1', 'Dopamine^Drd2'],
cmap='cividis_r', vmax=1, layer='pvals',
**kwargs)
... storing 'x' as categorical
We see that interactions with Dopamine as largely anticipated are predominantly located within the Striatum of the intact hemisphere, and are typically absent in the lesioned hemisphere.