mechanochat.SEM3#

class mechanochat.SEM3(ne_per_cell: int, re: float, rd_ratio: float, adata: AnnData | None = None, cluster_key: str = 'leiden', spatial_key: str = 'spatial', xc: ndarray[Any, dtype[_ScalarType_co]] | None = None, ctype: ndarray[Any, dtype[_ScalarType_co]] | None = None, sim_name: str = 'untitled', seed: int = 1)[source]#

Bases: CellBase

Subcellular Element Method 3D

Create a SEM object

Parameters:

ne_per_cell (int) – Number of elements per cell
re (float) – Element radius
rd_ratio (float) –
Cell radius-distance ratio

rd_ratio>2: cell radius < cell distance/2, tissue with gaps

rd_ratio=2: cell radius = cell distance/2, no gaps (confluent tissue)

rd_ratio<2: cell radius > cell distance/2, overcrowded
adata (Anndata) –
Anndata with .obsm[spatial_key] for cell coordinates, .obs[cluster_key] for cell types

If not provided, xc and ctype are required
xc (Optional[np.ndarray]) – cell coordinates. Ignored, if adata.obsm[spatial_key] is provided
ctype (Optional[np.ndarray]) – cell types. Ignored, if adata.obs[cluster_key] is provided
cluster_key (str, default: 'leiden') – Key for cell type in .obs
spatial_key (str, default: 'spatial') – Key for spatial coordinates in .obsm
sim_name (str, default: 'untitled') – Simulation name

f_cc: csr_matrix#: cell-cell force matrix (symmetric, abs)

pairs: ndarray[Any, dtype[_ScalarType_co]]#

interaction element pairs [i,j]

pairs[k] = [i,j]

f_pairs: ndarray[Any, dtype[_ScalarType_co]]#

f_pairs[k] stores force j -> i

pairs[k] = [i,j]

scale: float#: normalization scaling factor

deltax: NDArray#: normalization scaling offset. shape=(1,dim)

xc: NDArray#: cell coordinates. shape=(nc,dim)

ne: int#: element number

ecid: NDArray#

(ne)

Type:: id of cell to which each element belongs. shape

ceidn: NDArray#

first elements id of each cell (indptr). shape=(nc+1)

Elements of cell i can be retrieve by ceidn[i]:ceidn[i+1].

Number of elements of cell i is ceidn[i+1]-ceidn[i]

xe: NDArray#: element coordinates. shape=(ne,dim)

sim_name: str#: simulation name

t: int#: time step

param: dict#: simulation parameters (Deprecated)

ns_default: int#: number of augment points for alpha-shape

get_param(df_adh: DataFrame, df_stiffness: DataFrame, embedding_key: str = 'X_pca', adh_range: Tuple[float, float] = (0.25, 4.), stiffness_range: Tuple[float, float] = (0.25, 0.5), w: float = 0.75, return_data: bool = False)[source]#

Compute cell mechanical parameters alpha (adhesion) and beta (stiffness) from gene expression data and store them in self.alpha_csr and self.beta.

Computation pipeline

Neighbor detection — build a KD-tree on cell centroids self.xc and find all cell pairs within a Chebyshev distance of 60 (simulation units). Only neighboring pairs contribute to adhesion and similarity terms.
Adhesion score (adhesion_data) — for each ligand-receptor (LR) pair in df_adh, the per-cell ligand activity l and receptor activity r are computed as the product of min-shifted subunit expression levels (supports multi-subunit ligands/receptors encoded as 'gene1_gene2' in the table). The symmetric interaction score for neighbor pair (i, j) is:
```
LR_ij = r_i * l_j + l_i * r_j
```
Scores from all LR pairs are summed, capped at the 95th percentile to reduce outlier influence, then scaled to adh_range.
Gene expression similarity (sim_data) — cosine similarity between cell embedding vectors (default: PCA, adata.obsm[embedding_key]) for each neighbor pair, computed as:
```
cos_sim_ij = (X_i · X_j) / (||X_i|| * ||X_j||)
```
Values are scaled to adh_range.
Alpha matrix — the cell-cell adhesion parameter for pair (i, j) is a weighted combination of expression similarity and LR-based adhesion:
```
alpha_ij = sim_data_s * (1 - a) + adhesion_data * a
```
Stored as a symmetric sparse matrix self.alpha_csr of shape (nc, nc).
Beta vector — mean expression of cytoskeleton-related genes in df_stiffness per cell, capped at the 95th percentile and scaled to stiffness_range. Stored as self.beta of shape (nc,).

Parameters:

df_adh (DataFrame) – Ligand-receptor pair table. Column 0 = ligand gene(s), column 1 = receptor gene(s). Multi-subunit complexes should be joined with '_' (e.g. 'ITGA1_ITGB1').
df_stiffness (DataFrame) – Cytoskeleton/stiffness gene table. Column 0 contains gene names whose mean expression determines per-cell stiffness.
embedding_key (str, optional) – Key in adata.obsm for the cell embedding used to compute gene expression similarity. Default 'X_pca'.
adh_range (tuple of float, optional) – (min, max) range that both adhesion_data and sim_data are scaled to before combining. Default (0.25, 4.0).
stiffness_range (tuple of float, optional) – (min, max) range that stiffness values are scaled to. Default (0.25, 0.5).
w (float, optional) – Mixing weight for adhesion vs. similarity in the alpha parameter. w=1 → pure LR adhesion; w=0 → pure gene expression similarity. Default 0.75.
return_data (bool, optional) – If True, return intermediate arrays for inspection. Default False.

Returns:

None – Sets self.alpha_csr (sparse (nc, nc) float32) and self.beta ((nc,) float32) in place.
row_indices, col_indices, lr, sim_data (only when ``return_data=True``) – row_indices, col_indices — COO indices of neighbor pairs (symmetric); lr — raw (unscaled) symmetric LR adhesion scores; sim_data — raw (unscaled) symmetric cosine similarity values.

sim_gpu(rm_intra, rm_inter, sigma, dt, T, enable_tqdm=True) → None[source]#

Implement the SEM simulation on the GPU.

Parameters:

rm_intra (float) – Equilibrium distance of the intracellular element-element potential. Sets the packing of the elements inside a cell.
rm_inter (float) – Equilibrium distance of the intercellular element-element potential. Also the default contact threshold used by compute_contact and the alpha-shape radius (_get_e_radius).
sigma (float) – Amplitude of the thermal noise.
dt (float) – Time step.
T (int) – Number of steps.
enable_tqdm (bool, default: True) – Show the progress bar.

Returns:

None
Sets the following fields in `sem` (update `.xe (element coordinates)` and .xc (cell centroids),) – advances .t by T, stores rm_intra / rm_inter in .param, and marks the alpha-shapes as outdated.

Notes

Requires get_param() to have been run, which provides .alpha_csr and .beta.

Intercellular interactions are only evaluated for the cell pairs stored in .alpha_csr and, inside such a pair, only for the elements closer than 30 simulation units along each axis.

NUMBA_CUDA_USE_NVIDIA_BINDING=1 must be set before numba is imported,

compute_force_chunked(alpha=None, beta=None, gamma=None, chunk_size=300000000) → None[source]#

Compute the element-element forces of the current configuration.

All element pairs closer than 30 simulation distance (Chebyshev) are enumerated, and the force of each pair is evaluated on the GPU with the same parameters as sim_gpu. Pairs are processed in chunks to bound the GPU memory.

Parameters:

chunk_size (int, default: 300000000) – Number of element pairs sent to the GPU per launch. Lower it if the GPU runs out of memory.

Returns:

None
Sets the following fields in `sem` –
- .pairs : (n_pairs, 2) element pair indices, pairs[k] = [i, j]
- .f_pairs : (n_pairs, dim) the force from j exerted on i for pair k
- .f_cc and adata.obsp[‘f_cc’]: (n_cell, n_cell) cell-cell mechanical interaction matrix.

Notes

Requires .param.

sim_gpu() must have been run first.

compute_cell_forces() → None[source]#

Aggregate the element-element forces into a cell-cell mechanical interaction matrix.

The magnitudes of the forces of all element pairs that belong to two different cells are summed, which gives a symmetric measure of how strongly two cells push against or pull on each other. Entries below 0.1 are dropped.

Returns:

None
Sets the following fields in `sem` – Adds .f_cc, a symmetric (nc, nc) sparse matrix, and copies it to adata.obsp[‘f_cc’] when the SEM object is linked to an AnnData. Called automatically at the end of compute_force().

f_pairs_tocsr() → Tuple[csr_array, csr_array, csr_array][source]#

Convert force pairs to CSR sparse matrices. Run .compute_force() first

fx_csr[i, j] = fx_{j->i}

row index i : under force column index j : exert force

Returns:

fx_csr (csr_array) – Sparse matrix of x-component forces
fy_csr (csr_array) – Sparse matrix of y-component forces

compute_shape(alpha: ndarray[Any, dtype[_ScalarType_co]] | float | None = None, xe_data: ndarray[Any, dtype[_ScalarType_co]] | None = None, r: float | None = None, recompute: bool = False)[source]#

Reconstruct the 3D shape of each cell with an alpha-shape.

Every element is replaced by the 12 vertices of an icosahedron of radius r * scale, giving the point cloud its finite element thickness, and an AlphaShape3 is fitted per cell.

Parameters:

alpha (float or numpy.ndarray, optional) – Alpha radius, a scalar shared by all cells or one value per cell. If None, each cell uses max(r * scale * 1.01, 1.1 * alpha_best), where alpha_best is the smallest radius that keeps the shape connected.
xe_data (numpy.ndarray, optional) – (ne, …) per-element values (e.g. a signal or a pressure) carried over to the vertices of the mesh so that they can be used to color the surface.
r (float, optional) – Element radius used to build the icosahedra. Default: _get_e_radius() / 2, i.e. half of rm_inter.
recompute (bool, default: False) – Rebuild the shapes even if they are already cached. When the shapes exist and only alpha changed, the existing meshes are updated instead of rebuilt.

Returns:

None – Adds .shp_list (one AlphaShape3 per cell), fills .alpha_radius and updates .alphashape_info.