spateo.alignment.methods.paste ============================== .. py:module:: spateo.alignment.methods.paste Functions --------- .. autoapisummary:: spateo.alignment.methods.paste.paste_pairwise_align spateo.alignment.methods.paste.center_NMF spateo.alignment.methods.paste.paste_center_align spateo.alignment.methods.paste.generalized_procrustes_analysis Module Contents --------------- .. py:function:: paste_pairwise_align(sampleA: anndata.AnnData, sampleB: anndata.AnnData, layer: str = 'X', genes: Optional[Union[list, numpy.ndarray]] = None, spatial_key: str = 'spatial', alpha: float = 0.1, dissimilarity: str = 'kl', G_init=None, a_distribution=None, b_distribution=None, norm: bool = False, numItermax: int = 200, numItermaxEmd: int = 100000, dtype: str = 'float32', device: str = 'cpu', verbose: bool = True) -> Tuple[numpy.ndarray, Optional[int]] Calculates and returns optimal alignment of two slices. :param sampleA: Sample A to align. :param sampleB: Sample B to align. :param layer: If `'X'`, uses ``sample.X`` to calculate dissimilarity between spots, otherwise uses the representation given by ``sample.layers[layer]``. :param genes: Genes used for calculation. If None, use all common genes for calculation. :param spatial_key: The key in `.obsm` that corresponds to the raw spatial coordinates. :param alpha: Alignment tuning parameter. Note: 0 <= alpha <= 1. When α = 0 only the gene expression data is taken into account, while when α =1 only the spatial coordinates are taken into account. :param dissimilarity: Expression dissimilarity measure: ``'kl'`` or ``'euclidean'``. :param G_init: Initial mapping to be used in FGW-OT, otherwise default is uniform mapping. :param a_distribution: Distribution of sampleA spots, otherwise default is uniform. :param b_distribution: Distribution of sampleB spots, otherwise default is uniform. :param norm: If ``True``, scales spatial distances such that neighboring spots are at distance 1. Otherwise, spatial distances remain unchanged. :param numItermax: Max number of iterations for cg during FGW-OT. :param numItermaxEmd: Max number of iterations for emd during FGW-OT. :param dtype: The floating-point number type. Only float32 and float64. :param device: Equipment used to run the program. You can also set the specified GPU for running. E.g.: '0'. :param verbose: If ``True``, print progress updates. :returns: Alignment of spots. obj: Objective function output of FGW-OT. :rtype: pi .. py:function:: center_NMF(n_components, random_seed, dissimilarity='kl') .. py:function:: paste_center_align(init_center_sample: anndata.AnnData, samples: List[anndata.AnnData], layer: str = 'X', genes: Optional[Union[list, numpy.ndarray]] = None, spatial_key: str = 'spatial', lmbda: Optional[numpy.ndarray] = None, alpha: float = 0.1, n_components: int = 15, threshold: float = 0.001, max_iter: int = 10, numItermax: int = 200, numItermaxEmd: int = 100000, dissimilarity: str = 'kl', norm: bool = False, random_seed: Optional[int] = None, pis_init: Optional[List[numpy.ndarray]] = None, distributions: Optional[List[numpy.ndarray]] = None, dtype: str = 'float32', device: str = 'cpu', verbose: bool = True) -> Tuple[anndata.AnnData, List[numpy.ndarray]] Computes center alignment of slices. :param init_center_sample: Sample to use as the initialization for center alignment; Make sure to include gene expression and spatial information. :param samples: List of samples to use in the center alignment. :param layer: If `'X'`, uses ``sample.X`` to calculate dissimilarity between spots, otherwise uses the representation given by ``sample.layers[layer]``. :param genes: Genes used for calculation. If None, use all common genes for calculation. :param spatial_key: The key in `.obsm` that corresponds to the raw spatial coordinates. :param lmbda: List of probability weights assigned to each slice; If ``None``, use uniform weights. :param alpha: Alignment tuning parameter. Note: 0 <= alpha <= 1. When α = 0 only the gene expression data is taken into account, while when α =1 only the spatial coordinates are taken into account. :param n_components: Number of components in NMF decomposition. :param threshold: Threshold for convergence of W and H during NMF decomposition. :param max_iter: Maximum number of iterations for our center alignment algorithm. :param numItermax: Max number of iterations for cg during FGW-OT. :param numItermaxEmd: Max number of iterations for emd during FGW-OT. :param dissimilarity: Expression dissimilarity measure: ``'kl'`` or ``'euclidean'``. :param norm: If ``True``, scales spatial distances such that neighboring spots are at distance 1. Otherwise, spatial distances remain unchanged. :param random_seed: Set random seed for reproducibility. :param pis_init: Initial list of mappings between 'A' and 'slices' to solver. Otherwise, default will automatically calculate mappings. :param distributions: Distributions of spots for each slice. Otherwise, default is uniform. :param dtype: The floating-point number type. Only float32 and float64. :param device: Equipment used to run the program. You can also set the specified GPU for running. E.g.: '0'. :param verbose: If ``True``, print progress updates. :returns: - Inferred center sample with full and low dimensional representations (W, H) of the gene expression matrix. - List of pairwise alignment mappings of the center sample (rows) to each input sample (columns). .. py:function:: generalized_procrustes_analysis(X, Y, pi) Finds and applies optimal rotation between spatial coordinates of two layers (may also do a reflection). :param X: np array of spatial coordinates. :param Y: np array of spatial coordinates. :param pi: mapping between the two layers output by PASTE. :returns: Aligned spatial coordinates of X, Y and the mapping relations.