import scomv
import stlearn as st
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches
import cv2
import pandas as pd
import numpy as np
import anndata
import scanpy as sc
import seaborn as sns
from adjustText import adjust_text
from scipy.cluster.hierarchy import linkage, dendrogram, ClusterNode, to_tree
# Download tutorial dataset
#!mkdir tutorial_data
#!wget -P tutorial_data/ https://raw.githubusercontent.com/RyosukeNomural/SpatialCompassV/main/docs/tutorials/tutorial_data/xenium_data/cell_feature_matrix.h5
#!wget -P tutorial_data/ https://raw.githubusercontent.com/RyosukeNomural/SpatialCompassV/main/docs/tutorials/tutorial_data/xenium_data/cells.csv.gz
# Load Xenium data using stlearn
adata = st.ReadXenium(
    feature_cell_matrix_file="tutorial_data/cell_feature_matrix.h5",
    cell_summary_file="tutorial_data/cells.csv.gz",
    library_id="example data",
    image_path=None,
    scale=1,
    spot_diameter_fullres=10
)
Warning: Using default pixel size of 0.2125 microns. Consider providing experiment_xenium_file for accurate pixel size.
# Gridding at 10μm interval using stlearn
N_COL = int((adata.obs.imagecol.max() - adata.obs.imagecol.min()) / 10)
N_ROW = int((adata.obs.imagerow.max() - adata.obs.imagerow.min()) / 10)
grid = st.tl.cci.grid(adata, n_row=N_ROW, n_col=N_COL, n_cpus=10, verbose=False)
from scomv.preparation.skny_calc_distance import calculate_distance
## apply SKNY
grid = calculate_distance(
    grid, pos_marker_ls=['CDH1',"EPCAM"],
)
# Figure 2B
fig, ax = plt.subplots(figsize=(8, 6), dpi=150)
ax.imshow(cv2.cvtColor(grid.uns["marker_median_delineation"], cv2.COLOR_BGR2RGB),
          interpolation="nearest")
ax.axis("off")
plt.show()

#plt.savefig("Figure_2B.png", dpi=300)
../../_images/dbabf6d7b35995f646ca6915588682263042099d722526f6265f40483b2af4ba.png
import plotly.express as px
import plotly.io as pio
import cv2
import numpy as np

# Convert the OpenCV BGR image to RGB for Plotly visualization
img = cv2.cvtColor(grid.uns["marker_median_delineation"], cv2.COLOR_BGR2RGB)

# Minimum coordinates of the original spatial data
x_min = float(adata.obs.imagecol.min())
y_min = float(adata.obs.imagerow.min())

# Grid bin size used when generating the grid image
bin_size = 10

# Create Plotly image figure
fig = px.imshow(img)

h, w = img.shape[:2]

# Convert grid coordinates back to absolute spatial coordinates
xx_abs = (np.arange(w)[None, :] * bin_size + x_min).repeat(h, axis=0)
yy_abs = (np.arange(h)[:, None] * bin_size + y_min).repeat(w, axis=1)

# Round coordinates to the nearest multiple of 10 for display
xx_round10 = np.round(xx_abs / 10) * 10
yy_round10 = np.round(yy_abs / 10) * 10

# Stack the coordinates so they can be accessed in hover data
custom = np.dstack([xx_round10, yy_round10])

# Axis tick configuration
# step = 20 in grid space corresponds to ~200 units in absolute coordinates
step = 20
x_tickvals = list(range(0, w, step))
y_tickvals = list(range(0, h, step))

# Generate axis labels rounded to multiples of 10
x_ticktext = [str(int(round((v * bin_size + x_min) / 10) * 10)) for v in x_tickvals]
y_ticktext = [str(int(round((v * bin_size + y_min) / 10) * 10)) for v in y_tickvals]

fig.update_xaxes(
    tickmode="array",
    tickvals=x_tickvals,
    ticktext=x_ticktext
)

fig.update_yaxes(
    tickmode="array",
    tickvals=y_tickvals,
    ticktext=y_ticktext
)

# Customize hover information to display rounded absolute spatial coordinates
fig.update_traces(
    customdata=custom,
    hovertemplate=
    "imagecol: %{customdata[0]:.0f}<br>"
    "imagerow: %{customdata[1]:.0f}"
    "<extra></extra>"
)

fig.update_layout(
    xaxis_title="imagecol",
    yaxis_title="imagerow",
    height=700
)

# Use the interactive notebook renderer
pio.renderers.default = "notebook_connected"

fig.show()
# annotation each section to obs object
df_shotest = getattr(grid, "shortest")
df_grid = grid.to_df()

# extract grid info
df_grid = pd.merge(
    pd.DataFrame(index=["grid_" + str(i+1) for i in range(N_ROW * N_COL)]),
    df_grid, right_index=True, left_index=True, how="left"
).fillna(np.nan)

# extract section info
df_region = pd.DataFrame(
    np.array(df_shotest["region"]).reshape(N_ROW, N_COL).T.reshape(N_ROW * N_COL),
    index=["grid_" + str(i+1) for i in range(N_ROW * N_COL)], columns=["region"]
)

# marge
df_grid_region = pd.merge(
    df_grid, df_region,
    right_index=True, left_index=True, how="left"
)
df_grid_region = df_grid_region.dropna()

# add to obs
grid.obs = pd.merge(
    grid.obs, df_grid_region[["region"]],
    right_index=True, left_index=True, how="left"
)

# shaping
grid.obs["region_10"] = [str(i*10) for i in grid.obs["region"]]
grid.obs["region_10"] = ["("+str(int(float(i.split(", ")[0][1:])))+", "+str(int(float(i.split(", ")[-1][:-1])))+"]" if i != "nan" else np.nan for i in grid.obs["region_10"]]
# exclude because of small number
grid.obs["region_10"] = grid.obs["region_10"].replace(
    {"(-150, -120]": np.nan}
)
from scomv.preparation.choose_roi import extract_roi, contour_regions
# select ROI
roi = (2400, 3400, 2400, 3800)

subset_grid, filtered_shortest, xy_list = extract_roi(
    grid=grid,
    roi=roi,
    bin_size=10,
    region_col="region_10",
)
g_x_cont, g_y_cont, g_x_inside, g_y_inside = contour_regions(
    filtered_shortest, adata,
    min_x=0, max_x=4000, min_y=0, max_y=4000,
    #out_dir_base="./roi_plots",
    show=True,
)
../../_images/c9d3fa38c331b8ab17aed3bbf368273d5f77bd134bbd8dc07b3a78cc3d1fa1c0.png
from scomv.preparation.scomv_calc_vector import compute_min_vectors_polar
outline_points = list(zip(g_x_cont, g_y_cont))
inside_points  = list(zip(g_x_inside, g_y_inside))

min_vector_df = compute_min_vectors_polar(
    xy_list=xy_list,
    outline_points=outline_points,
    inside_points=inside_points,
    invert_y=True,
    make_inside_negative=True,
)
from scomv.gene_view_tool import plot_gene_polar_hist2d, plot_2d_expression
## Draw polar maps
genes = list(subset_grid.var.index)
A_counts, total_n, idx = plot_gene_polar_hist2d(
    subset_grid=subset_grid,
    min_vector_df=min_vector_df,
    gene="CD3E",
)
CD3E: n=3431
../../_images/a4f36eb2240214f4172d586c3ed79d9d5d95a0ec76e65b4a4c4ab6e910c31aa1.png
## Draw 2D pictures
genes = ["CD3E"] # the genes you want to draw
for g in genes:
    plot_2d_expression(g, subset_grid)
Please specify a valid `library_id` or set it permanently in `adata.uns['spatial']`
../../_images/9f41989b1edeb889f6f56509cb88a8800970747bd16e87936dd00e8a15342011.png
from scomv.gene_pipeline import SCOMVPipeline

gene_pipe = SCOMVPipeline(adata=adata, grid=grid, bin_size=10, region_col="region_10")

gene_out =gene_pipe.run(
    roi=roi,
    subset_grid=subset_grid,
    filtered_shortest=filtered_shortest,
    xy_list=xy_list,
    outline_points=outline_points,
    inside_points=inside_points,
)
/Users/riuyamas/miniconda3/envs/scomv-pypi-install-2/lib/python3.11/site-packages/scanpy/metrics/_morans_i.py:105: UserWarning:

1 variables were constant, will return nan for these.
/Users/riuyamas/miniconda3/envs/scomv-pypi-install-2/lib/python3.11/site-packages/numba/core/decorators.py:282: RuntimeWarning: nopython is set for njit and is ignored
  warnings.warn('nopython is set for njit and is ignored', RuntimeWarning)
import plotly.io as pio
gene_pipe.plot_pcoa_plotly(gene_pipe.coords, gene_pipe.selected_genes)
pcoa_two_dims = gene_out["coords"][["PC1", "PC2"]]
df = pcoa_two_dims

# --- Gene sets ---
t_cells  = {"CD3D", "CD3E", "CD4", "CD247", "CD8A", "GZMA", "TRAC"}
micro_f  = ["CD68", "CD80", "CD86", "CD163"]
b_cells  = ["CD19", "BANK1", "MS4A1", "TCL1A", "CD79A"]
mono     = ["CD14"]
epi      = ["CDH1", "EPCAM", "KRT8"]
smooth   = ["ACTA2", "MYH11", "ACTG2"]
endo     = ["CD93", "PECAM1", "VWF"]
fib      = ["POSTN", "LUM", "MMP2"]
focus    = ["PLVAP", "RHOA", "SMAD6", "TP53", "VEGFA", "MKI67", "ADGRL4", "RGS5",
            "PI3", "CAV1", "EGFR", "AQP3", "C15orf48"]

# --- Others (everything except T-cell set) ---
others = sorted(list(set(df.index) - set(t_cells)))

# --- Styling config (keeps your colors/labels/sizes) ---
groups = [
    ("T cell",          t_cells,  "#1f77b4", 200, True),
    ("B cell",          b_cells,  "cyan",    200, True),
    ("epithelial cell", epi,      "#ffbb78", 200, True),
    ("macrophage",      micro_f,  "#2ca02c", 200, True),
    ("monocyte",        mono,     "#ff9896", 200, True),
    ("smooth muscle",   smooth,   "#aec7e8", 200, True),
    ("endotherial",     endo,     "#8B4513", 200, True), 
    ("fibroblast",      fib,      "red",     200, True),
    (None,              focus,    "black",   200, True),  # no legend label
]

# --- Plot ---
plt.figure(figsize=(20, 15))
sns.set_theme(style="white")
texts = []

# Others: gray points (no labels)
for g in others:
    if g in df.index:
        plt.scatter(df.loc[g, "PC1"], df.loc[g, "PC2"], color="gray", alpha=0.3, s=150)

# Highlighted groups
for label, gene_list, color, s, add_text in groups:
    for g in gene_list:
        if g in df.index:
            if label is None:
                plt.scatter(df.loc[g, "PC1"], df.loc[g, "PC2"], color=color, s=s)
            else:
                plt.scatter(df.loc[g, "PC1"], df.loc[g, "PC2"], color=color, label=label, s=s)
            if add_text:
                texts.append(plt.text(df.loc[g, "PC1"], df.loc[g, "PC2"], g, fontsize=25))

plt.xlabel("PCoA1", fontsize=30)
plt.ylabel("PCoA2", fontsize=30)
plt.xticks(fontsize=25)
plt.yticks(fontsize=25)
plt.title("PCoA Plot", fontsize=40)

# Legend (remove duplicates)
handles, labels = plt.gca().get_legend_handles_labels()
by_label = dict(zip(labels, handles))
plt.legend(by_label.values(), by_label.keys(), loc="best", prop={"size": 23}, ncol=1)

adjust_text(
    texts,
    arrowprops=dict(arrowstyle="->", color="gray", lw=0.5),
    only_move={"points": "xy", "text": "xy"},
    force_text=1.5,
    force_points=1.0,
    expand_text=(1.2, 1.2),
    expand_points=(1.2, 1.2),
    lim=200,
    precision=0.01,
    autoalign="y",
)

plt.show()
Looks like you are using a tranform that doesn't support FancyArrowPatch, using ax.annotate instead. The arrows might strike through texts. Increasing shrinkA in arrowprops might help.
../../_images/c03a057b7b58b5490cdf4317fac7f7121dd6e5d7fc69f569cad1f649dcaf7549.png
from scomv.plot_3d import plot_3d
import plotly.io as pio

fig = plot_3d(
    adata=gene_out["subset_grid"],
    genes=("POSTN", "CD3E"),
    anchor_gene="CDH1",
    bin_size=10,
    x_range=(2400, 3100),
    y_range=(2400, 3400),
    height_scale=3,
    anchor_scale=1,
    threshold=4,
    sigma=1.0,
)
pio.renderers.default = "notebook_connected"
fig.show()
from scomv.dendrogram import dendrogram2newick
z = linkage(gene_out["dist_df"], metric='euclidean', method='ward')
tree= to_tree(z)
newick = dendrogram2newick(tree, tree.dist, gene_out["selected_genes"])
/var/folders/rp/843_v26x35314skd4rnkz6jm0000gp/T/ipykernel_37409/1097620655.py:2: ClusterWarning:

scipy.cluster: The symmetric non-negative hollow observation matrix looks suspiciously like an uncondensed distance matrix
# The Newick output can be used to visualize a dendrogram in iTOL(https://www.google.com/url?sa=t&source=web&rct=j&opi=89978449&url=https://itol.embl.de/&ved=2ahUKEwiY7vnzxoOSAxWRrlYBHUuiHisQFnoECBoQAQ&usg=AOvVaw2ipidUpg8SauabV53lQ6Jh).
#newick