课前准备--CODEX（IMC）基础分析梳理

import sys
import subprocess
import numpy as np
import pandas as pd
import scanpy as sc
import anndata
from anndata import AnnData
from matplotlib import pyplot as plt
import matplotlib as mpl
import seaborn as sns
import os

import monkeybread as mb

sc.pp.normalize_total(adata, inplace=True)
sc.pp.log1p(adata)
sc.pp.pca(adata)

# Display high-level cell types
fig, ax = plt.subplots(1, 1, figsize=(5,5))
sc.pl.embedding(
    adata,
    "spatial",
    color = 'Cell Type',
    s=1,
    ax=ax,
    palette=sc.pl.palettes.vega_20_scanpy,
    show=False
)
plt.show()


# Display more granular cell subtypes
fig, ax = plt.subplots(1, 1, figsize=(5,5))
sc.pl.embedding(
    adata,
    "spatial",
    color = 'Cell Subtype',
    s=1,
    ax=ax,
    palette=sc.pl.palettes.godsnot_102,
    show=False
)
plt.show()

# Display expression dotplot
sc.tl.dendrogram(adata, groupby='Cell Subtype')
dp = sc.pl.dotplot(
    adata, 
    var_names={
        'T cell': [
            'CD3D'
        ],
        'CD8 T': [
            'GZMB',
            'CD8A'
        ],
        'CD4 T': [
            'CD4'
        ],
        'Treg': [
            'FOXP3', 
            'TIGIT'
        ],
        'NK': [
            'NCAM1', 
            'KLRB1'
        ],
        'B cell': [
            'MS4A1',
            'CD79A',
            'CD79B',
            'CD19'
        ],
        'mono/mac': [
            'LYZ', 
            'CD14', 
            'FCGR3A', 
            'CD163',
            'CD68'
        ],
        'dendritic': [
            'ITGAX',
            'FCER1A',
            'CD1C'
        ],
        'mregDC': [
            'CCR7',
            'LAMP3',
            'CD83'
        ],
        'pDC': [
            'IL3RA',
            'CLEC4C'
        ],
        'mast': [
            'CTSG',
            'KIT'
        ],
        'epithelial': [
            'EPCAM',
            'MUC1'
        ],
        'endothelial': [
            'VWF', 
            'PECAM1'
        ],
        'stromal': [
            'MMP2',
            'COL1A1',
            'COL5A1'
        ],
        'misc': [
            'CXCL9',
            'CXCL10',
            'SIGLEC1',
            'PDCD1',
            'PRF1'
        ]
    },
    groupby=['Cell Subtype'],
    standard_scale='var',
    cmap='bwr',
    dendrogram=True,
    return_fig=True
)
dp.add_totals().show()

density_key = mb.calc.cell_density(
    adata_zoom,
    'Cell Type',
    'B lineage',
    bandwidth=25,
    approx=False,
    radius_threshold=250
)
mb.plot.cell_density(
    adata_zoom,
    density_key,
    cmap='Blues',
    title=f'B cell density'
)

fig, _ = mb.plot.location_and_density(
    adata, 
    'Cell Subtype',                   
    [
        ['B lineage'],  # The first group of cell types to plot
        ['CD4 T (Tfh)'] # The second group of cell types to plot
    ],      
    [
        'B cell', # Name of the first group of cell types
        'CD4 Tfh' # Name of the second group of cell types
    ],
    dot_size=[1, 1], # Dot sizes for each group
    title='B cell vs. CD4 Tfh', 
    grid=True, 
    n_grids=5,
    show=False
)
plt.tight_layout()
plt.show()

niche_subtypes = []
other_cell_types = ['other', 'malignant']
for cell, ct in zip(adata.obs.index, adata.obs['Cell Subtype']):
    if ct in other_cell_types:
        niche_subtypes.append('malignant/other')
    else:
        niche_subtypes.append(ct)
adata.obs['niche_subtypes'] = niche_subtypes
adata.obs['niche_subtypes'] = adata.obs['niche_subtypes'].astype('category')

print("Cell subtypes considered in niche analysis:")
print(set(adata.obs['niche_subtypes']))
# Only perform niche analysis on immune cells
immune_mask = ~adata.obs['niche_subtypes'].isin([
    'malignant/other', 'endothelial', 'stromal'
])

# Compute niches
adata_neighbors = mb.calc.cellular_niches(
    adata,
    cell_type_key='niche_subtypes',
    radius=75,
    normalize_counts=True,
    standard_scale=True,
    clip_min=-5,
    clip_max=5,
    mask=immune_mask,
    n_neighbors=100,
    resolution=0.25,
    min_niche_size=300,
    key_added='niche',
    non_niche_value='malignant/other'
)

# Subset cells
adata_neighbors_sub = sc.pp.subsample(
    adata_neighbors, fraction=1/2, copy=True
)

# Generate UMAP plots
sc.tl.umap(adata_neighbors_sub)
sc.pl.umap(
    adata_neighbors_sub, 
    color='niche_subtypes', 
    palette=sc.pl.palettes.vega_20_scanpy
)
sc.pl.umap(
    adata_neighbors_sub, 
    color='niche', 
    palette=sc.pl.palettes.vega_20_scanpy
)

# Map each niche to a color so that plots are consistent
niche_to_color = {
    val: mb.plot.monkey_palette[i]
    for i, val in enumerate(sorted(set(adata_neighbors.obs['niche'])))
}
niche_to_color['malignant/other'] = 'lightgrey'

fig, ax = plt.subplots(1,1,figsize=(6,6))
sc.pl.embedding(
    adata,
    "spatial",
    color = 'niche',
    palette=niche_to_color,
    s=1,
    ax=ax,
    show=False
)

mb.plot.neighbors_profile_matrixplot(
    adata_neighbors,
    'niche',
    include_niches=[ # Exclude the miscellaneous niches from the plot
        niche 
        for niche in set(adata_neighbors.obs['niche'])
        if niche != 'malignant/other'
    ],
    clustermap_kwargs={
        'linewidths': 0.5,
        'linecolor': 'black',
        'cmap': 'bwr',
        'clip_on': False,
        'vmin': -3.5,
        'vmax': 3.5,
        'figsize': (8,8)
    }
)