没听错吧！单细胞中也可以进行免疫浸润分析？SpatialCells带给你新的体验！

pip install -r requirements.txtpip install .

conda env create --name spatialcells --file=conda.yamlpip install .conda.yaml 指定如下：name: spatial-cells-env               channels:  - conda-forge  - defaultsdependencies:  - python>=3.7, <3.11  - ipykernel  - matplotlib>=3.7.0  - pandas>=2.0.3  - seaborn>=0.12.2  - shapely>=2.0  - tqdm  - pip  - pip:      - anndata==0.9.2      - scanpy==1.9.4           加载相关的python包import pandas as pdimport matplotlib.pyplot as pltimport anndata as adimport spatialcells as spc             

adata = ad.read_h5ad("../../data/MEL1_adata.h5ad")spc.prep.setGate(adata, "KERATIN_cellRingMask", 6.4, debug=True)spc.prep.setGate(adata, "SOX10_cellRingMask", 7.9, debug=True)spc.prep.setGate(adata, "CD3D_cellRingMask", 7, debug=True)

marker = ["SOX10_cellRingMask_positive"]communitycolumn = "COI_community"ret = spc.spatial.getCommunities(adata, marker, eps=60, newcolumn=communitycolumn)fig, ax = plt.subplots(figsize=(10, 8))spc.plt.plotCommunities(    adata, ret, communitycolumn, plot_first_n_clusters=10, s=2, fontsize=10, ax=ax)    ax.invert_yaxis()plt.show()

communityIndexList = [6, 3, 14, 51, 29, 47, 39, 44, 22]boundary = spc.spatial.getBoundary(    adata, communitycolumn, communityIndexList, alpha=130)boundary = spc.spa.pruneSmallComponents(boundary, min_edges=50, holes_min_edges=500)roi_boundary = spc.spa.getExtendedBoundary(boundary, offset=2000)markersize = 1    fig, ax = plt.subplots(figsize=(10, 7))## all pointsax.scatter(    *zip(*adata.obs[["X_centroid", "Y_centroid"]].to_numpy()),    s=markersize,    color="grey",    alpha=0.2)# Points in selected commnitiesxy = adata.obs[adata.obs[communitycolumn].isin(communityIndexList)][    ["X_centroid", "Y_centroid"]].to_numpy()ax.scatter(xy[:, 0], xy[:, 1], s=markersize, color="r")# Bounds of points in selected commnitiesspc.plt.plotBoundary(boundary, ax=ax, label="Boundary", color="b")spc.plt.plotBoundary(roi_boundary, ax=ax, label="ROI boundary", color="g")ax.invert_yaxis()plt.show() 

spc.spatial.assignPointsToRegions(    adata,    [boundary, roi_boundary],    ["Tumor", "Tumor_ROI"],    assigncolumn="region",    default="BG",)point_size = 1fig, ax = plt.subplots(figsize=(10, 7))for region in sorted(set(adata.obs["region"])):        tmp = adata.obs[adata.obs.region == region]    ax.scatter(        *zip(*tmp[["X_centroid", "Y_centroid"]].to_numpy()),        s=point_size,        alpha=0.7,        label=region    )# Bounds of points in selected commnitiesspc.plt.plotBoundary(boundary, ax=ax, label="Boundary", color="purple")spc.plt.plotBoundary(roi_boundary, ax=ax, label="ROI boundary", color="r")plt.legend(loc="upper right")ax.invert_yaxis()plt.show()

def merge_pheno(row):    if row["phenotype_large_cohort"] in [        "T cells",        "Cytotoxic T cells",        "Exhausted T cells",    ]:        return "T cells"    elif row["phenotype_large_cohort"] in ["Melanocytes"]:        return "Tumor cells"    else:        return "Other cells"def cell_type(row):    if row["SOX10_cellRingMask_positive"]:        return "SOX10+"    elif row["CD3D_cellRingMask_positive"]:        return "CD3D+"    else:        return "Other cells"# Applying the function to create the new columnsadata.obs["pheno1"] = pd.Categorical(adata.obs.apply(merge_pheno, axis=1))    adata.obs["Cell Types"] = pd.Categorical(adata.obs.apply(cell_type, axis=1))spc.msmt.getRegionComposition(adata, "pheno1")

melano = adata[    (adata.obs.SOX10_cellRingMask_positive) & (adata.obs.region.isin(["Tumor_ROI", "Tumor"]))]tcells = adata[    (adata.obs.CD3D_cellRingMask_positive)    & (adata.obs.region.isin(["Tumor_ROI", "Tumor"]))]           fig, ax = plt.subplots(figsize=(10, 7))ax.invert_yaxis()ax.set_aspect("equal")plt.scatter(    tcells.obs["X_centroid"],    tcells.obs["Y_centroid"],    s=0.5,    label="T cells",        color="green",    alpha=0.5,)spc.plt.plotBoundary(roi_boundary, ax=ax, label="ROI boundary", color="r")           plt.legend(loc="upper right", markerscale=5)plt.show()

tumor = adata[adata.obs.region.isin(["Tumor_ROI", "Tumor"])]communitycolumn = "CD3D_cellRingMask_positive"communityIndexList = [True]           immune_boundary = spc.spatial.getBoundary(    tumor, communitycolumn, communityIndexList, alpha=130    )immune_boundary = spc.spa.pruneSmallComponents(    immune_boundary, min_edges=25, holes_min_edges=30, min_area=30000)markersize = 0.1fig, ax = plt.subplots(figsize=(10, 7))           ## all pointsax.scatter(    *zip(*adata.obs[["X_centroid", "Y_centroid"]].to_numpy()),    s=markersize,    color="grey",    alpha=0.2)           # Points in selected commnitiesxy = tumor.obs[tumor.obs[communitycolumn].isin(communityIndexList)][    ["X_centroid", "Y_centroid"]].to_numpy()ax.scatter(xy[:, 0], xy[:, 1], s=markersize, color="green", alpha=1, label="T cells")               # Bounds of points in selected commnitiesspc.plt.plotBoundary(    immune_boundary, ax=ax, label="Immune Cell Region Boundary", color="k", linewidth=1)           spc.plt.plotBoundary(roi_boundary, ax=ax, label="ROI boundary", color="r")# ax.set_xlim(0, 20000)# ax.set_ylim(0, 13000)ax.invert_yaxis()ax.set_aspect("equal")ax.set_axis_off()# plt.legend(loc="upper right", markerscale=5, fontsize=13.5)# plt.savefig("immune_cell_region1.png", dpi=400)plt.show() 

spc.spatial.assignPointsToRegions(    melano, [immune_boundary], ["T"], assigncolumn="tumor_isolated_region", default="F")

point_size = 0.5fig, ax = plt.subplots(figsize=(10, 7))ax.scatter(    *zip(*adata.obs[["X_centroid", "Y_centroid"]].to_numpy()),    s=markersize,    color="grey",        alpha=0.2)           colors = ["red", "orange"]labels = ["Immune-isolated Tumor Cells", "Immune-rich Tumor Cells"]for i, region in enumerate(sorted(set(melano.obs["tumor_isolated_region"]))):    tmp = melano.obs[melano.obs.tumor_isolated_region == region]    ax.scatter(        *zip(*tmp[["X_centroid", "Y_centroid"]].to_numpy()),        s=point_size,        alpha=0.5,        color=colors[i],        label=labels[i]    )           # Bounds of points in selected commnitiesspc.plt.plotBoundary(    immune_boundary, ax=ax, label="Immune Cell Region Boundary", color="k", linewidth=1)spc.plt.plotBoundary(roi_boundary, ax=ax, label="ROI boundary", color="b")    ax.invert_yaxis()ax.set_axis_off()# plt.savefig("roi_region1.png", dpi=400)plt.show()

print("Percentage of tumor cells in immune-isolated regions: ")melano.obs["tumor_isolated_region"].value_counts() / len(melano.obs)

roi_area = spc.msmt.getRegionArea(roi_boundary)tumor_area = spc.msmt.getRegionArea(boundary)immune_area = spc.msmt.getRegionArea(immune_boundary)           tumor_immune_overlap = boundary.intersection(immune_boundary)    overlap_area = spc.msmt.getRegionArea(tumor_immune_overlap)           print(f"Area of ROI: {roi_area:.2f}")print(f"Area of main tumor cell region: {tumor_area:.2f}")print(f"Area of immune cell region: {immune_area}")print(f"Area of overlap between tumor and immune cell regions: {overlap_area:.2f}")print(    f"Percentage of tumor region that has overlap with "    f"immune cell region: {overlap_area / tumor_area:.3f}")           识别与免疫细胞相邻的肿瘤细胞。dists = spc.msmt.getMinCellTypesDistance(melano, tcells)adata.obs.loc[    (adata.obs.SOX10_cellRingMask_positive) & (adata.obs.region.isin(["Tumor_ROI", "Tumor"])), "dist"] = diststhreshold = 20adata.obs["dist_binned"] = adata.obs["dist"] <= thresholdinfiltrated = adata.obs[    (adata.obs.SOX10_cellRingMask_positive)        & (adata.obs.region == "Tumor")    & (adata.obs.dist_binned == True)]non_infiltrated = adata.obs[    (adata.obs.SOX10_cellRingMask_positive)    & (adata.obs.region == "Tumor")    & (adata.obs.dist_binned == False)]           fig, ax = plt.subplots(figsize=(10, 6))ax.invert_yaxis()region = adata.obs[(adata.obs.region == "Tumor")]plt.scatter(region["X_centroid"], region["Y_centroid"], s=1, alpha=0.2, color="grey")plt.scatter(    infiltrated["X_centroid"],    infiltrated["Y_centroid"],    s=1,    alpha=0.5,    color="green",    label=f"Infiltrated (distance to t-cells <= {threshold}um)",)    plt.legend(markerscale=10)plt.show()