Repackaging RISE

RISE/factorization.py

@@ -6,7 +6,6 @@
 from scipy.stats import gmean
 from sklearn.decomposition import PCA
 from sklearn.linear_model import LinearRegression
-from sklearn.preprocessing import scale
 from tqdm import tqdm
 
 
@@ -95,8 +94,10 @@ def pf2(
         - X.uns["Pf2_B"]: Eigen-state factors (shape: rank, rank)
         - X.varm["Pf2_C"]: Gene factors (shape: n_genes, rank)
         - X.obsm["projections"]: Cell projections (shape: n_cells, rank)
-        - X.obsm["weighted_projections"]: Weighted cell projections (shape: n_cells, rank)
-        - X.obsm["X_pf2_PaCMAP"]: PaCMAP embedding (shape: n_cells, 2) if doEmbedding=True
+        - X.obsm["weighted_projections"]: Weighted cell projections
+          (shape: n_cells, rank)
+        - X.obsm["X_pf2_PaCMAP"]: PaCMAP embedding (shape: n_cells, 2)
+          if doEmbedding=True
     """
     pf_out, _ = parafac2_nd(
         X, rank=rank, random_state=random_state, tol=tolerance, n_iter_max=max_iter
@@ -145,7 +146,8 @@ def rise_pca_r2x(X: anndata.AnnData, ranks):
         r2x_rise[index] = R2X
 
     # Mean center because this is done within RISE
-    XX = scale(XX.todense(), with_mean=True, with_std=False)
+    XX = XX.toarray()
+    XX = XX - np.mean(XX, axis=0)
 
     pca = PCA(n_components=ranks[-1])
     pca.fit(XX)

RISE/plotting/__init__.py

@@ -0,0 +1,53 @@
+"""
+Plotting and visualization functions for RISE.
+"""
+
+from .factors import (
+    plot_condition_factors,
+    plot_eigenstate_factors,
+    plot_gene_factors,
+)
+from .general import (
+    avegene_per_status,
+    cell_count_perc_df,
+    cell_count_perc_lupus_df,
+    gene_plot_cells,
+    plot_avegene_per_category,
+    plot_avegene_per_celltype,
+    plot_cell_gene_corr,
+    plot_r2x,
+    rotate_xaxis,
+    rotate_yaxis,
+)
+from .pacmap import (
+    plot_gene_pacmap,
+    plot_labels_pacmap,
+    plot_wp_pacmap,
+)
+from .stability import (
+    calculateFMS,
+    plot_fms_diff_ranks,
+    plot_fms_percent_drop,
+)
+
+__all__ = [
+    "plot_condition_factors",
+    "plot_eigenstate_factors",
+    "plot_gene_factors",
+    "avegene_per_status",
+    "cell_count_perc_df",
+    "cell_count_perc_lupus_df",
+    "gene_plot_cells",
+    "plot_avegene_per_category",
+    "plot_avegene_per_celltype",
+    "plot_cell_gene_corr",
+    "plot_r2x",
+    "rotate_xaxis",
+    "rotate_yaxis",
+    "plot_gene_pacmap",
+    "plot_labels_pacmap",
+    "plot_wp_pacmap",
+    "calculateFMS",
+    "plot_fms_diff_ranks",
+    "plot_fms_percent_drop",
+]

RISE/figures/commonFuncs/plotFactors.py → RISE/plotting/factors.py

@@ -22,7 +22,8 @@ def plot_condition_factors(
     color_key=None,
     group_cond=False,
 ):
-    """Plot condition factors as a heatmap showing how conditions contribute to components.
+    """Plot condition factors as a heatmap showing how conditions contribute to
+    components.
 
     This visualization shows how each experimental condition (rows) contributes to
     each RISE component (columns). High values indicate strong association between
@@ -46,7 +47,8 @@ def plot_condition_factors(
         Series mapping conditions to group labels for colored row annotations.
         Useful for grouping related conditions (e.g., drug classes, patient cohorts).
     ThomsonNorm : bool, optional (default: False)
-        If True, normalizes factors using only control conditions (those containing 'CTRL').
+        If True, normalizes factors using only control conditions (those
+        containing 'CTRL').
     color_key : list, optional (default: None)
         Custom colors for condition group labels. If None, uses default palette.
     group_cond : bool, optional (default: False)
@@ -162,7 +164,8 @@ def plot_gene_factors(data: anndata.AnnData, ax: Axes, weight=0.08, trim=True):
 
     This visualization reveals coordinated gene modules by showing which genes (rows)
     are highly weighted in each component (columns). The weight parameter filters out
-    genes with low contributions, focusing on the most important genes for interpretation.
+    genes with low contributions, focusing on the most important genes for
+    interpretation.
 
     Parameters
     ----------

RISE/figures/commonFuncs/plotGeneral.py → RISE/plotting/general.py

@@ -4,7 +4,7 @@
 import seaborn as sns
 from matplotlib.axes import Axes
 
-from ...factorization import rise_pca_r2x
+from ..factorization import rise_pca_r2x
 
 
 def plot_r2x(data, rank_vec, ax: Axes):

RISE/figures/commonFuncs/plotPaCMAP.py → RISE/plotting/pacmap.py

@@ -75,7 +75,7 @@ def plot_gene_pacmap(gene: str, X: anndata.AnnData, ax: Axes, clip_outliers=0.99
     points = np.array(X.obsm["X_pf2_PaCMAP"])
 
     canvas = _get_canvas(points)
-    data = pd.DataFrame(points, columns=("x", "y"))
+    data = pd.DataFrame(points, columns=["x", "y"])  # type: ignore
 
     values -= np.min(values)
     values /= np.max(values)
@@ -126,7 +126,7 @@ def plot_wp_pacmap(X: anndata.AnnData, cmp: int, ax: Axes, cbarMax: float = 1.0)
     cmap = sns.diverging_palette(240, 10, as_cmap=True)
 
     canvas = _get_canvas(points)
-    data = pd.DataFrame(points, columns=("x", "y"))
+    data = pd.DataFrame(points, columns=["x", "y"])  # type: ignore
 
     values /= np.max(np.abs(values))
 
@@ -195,7 +195,7 @@ def plot_labels_pacmap(
     labels = labels.iloc[indices]
 
     canvas = _get_canvas(points)
-    data = pd.DataFrame(points, columns=("x", "y"))
+    data = pd.DataFrame(points, columns=["x", "y"])  # type: ignore
 
     data["label"] = pd.Categorical(labels)
     aggregation = canvas.points(data, "x", "y", agg=ds.count_cat("label"))

RISE/figures/commonFuncs/seriate.py → RISE/plotting/seriate.py

@@ -41,20 +41,22 @@ def seriate(
     dists: numpy.ndarray, approximation_multiplier: int = 1000, timeout: float = 2.0
 ) -> list[int]:
     """
-    Order the elements of a set so that the sum of sequential pairwise distances is minimal.
+    Order the elements of a set so that the sum of sequential pairwise
+    distances is minimal.
 
     We solve the Travelling Salesman Problem (TSP) under the hood.
     Reference: http://nicolas.kruchten.com/content/2018/02/seriation/
 
-    :param dists: Either a condensed pdist-like or a symmetric square distance matrix.
-    :param approximation_multiplier: Multiply by this number before converting distances \
-                                     to integers.
-    :param timeout: Maximum amount of time allowed to spend for solving the TSP, in seconds. \
-                    This value cannot be less than 0. If timeout is 0 it will try to solve \
-                    the problem with timeout = 1, and double the timeout every time it fails \
-                    until a valid solution is found.
-    :return: list with ordered element indexes, the same length as the number of elements \
-             involved in calculating `dists`.
+    :param dists: Either a condensed pdist-like or a symmetric square distance
+                  matrix.
+    :param approximation_multiplier: Multiply by this number before converting
+                                     distances to integers.
+    :param timeout: Maximum amount of time allowed to spend for solving the TSP,
+                    in seconds. This value cannot be less than 0. If timeout is 0
+                    it will try to solve the problem with timeout = 1, and double
+                    the timeout every time it fails until a valid solution is found.
+    :return: list with ordered element indexes, the same length as the number of
+             elements involved in calculating `dists`.
     """
     _validate_data(dists)
     if timeout > 0:
@@ -83,17 +85,20 @@ def _seriate(
     dists: numpy.ndarray, approximation_multiplier=1000, timeout=2.0
 ) -> list[int]:
     """
-    Order the elements of a set so that the sum of sequential pairwise distances is minimal.
+    Order the elements of a set so that the sum of sequential pairwise
+    distances is minimal.
 
     We solve the Travelling Salesman Problem (TSP) under the hood.
     Reference: http://nicolas.kruchten.com/content/2018/02/seriation/
 
-    :param dists: Either a condensed pdist-like or a symmetric square distance matrix.
-    :param approximation_multiplier: Multiply by this number before converting distances \
-                                     to integers.
-    :param timeout: Maximum amount of time allowed to spend for solving the TSP, in seconds.
-    :return: list with ordered element indexes, the same length as the number of elements \
-             involved in calculating `dists`.
+    :param dists: Either a condensed pdist-like or a symmetric square distance
+                  matrix.
+    :param approximation_multiplier: Multiply by this number before converting
+                                     distances to integers.
+    :param timeout: Maximum amount of time allowed to spend for solving the TSP,
+                    in seconds.
+    :return: list with ordered element indexes, the same length as the number of
+             elements involved in calculating `dists`.
     """
     assert dists[dists < 0].size == 0, "distances must be non-negative"
     assert timeout > 0

RISE/figures/figureS4.py → RISE/plotting/stability.py

@@ -1,5 +1,5 @@
 """
-Figure S4
+Decomposition stability analysis and plotting functions.
 """
 
 import anndata
@@ -12,26 +12,9 @@
 from tlviz.factor_tools import factor_match_score as fms
 
 from ..factorization import pf2
-from .common import getSetup, subplotLabel
 
-# from ..imports import import_thomson
 
-
-def makeFigure():
-    ax, f = getSetup((6, 3), (1, 2))
-    subplotLabel(ax)
-
-    # X = import_thomson()
-    # percentList = np.arange(0.0, 8.0, 5.0)
-    # plot_fms_percent_drop(X, ax[0], percentList=percentList, runs=3)
-
-    # ranks = list(range(1, 3))
-    # plot_fms_diff_ranks(X, ax[1], ranksList=ranks, runs=3)
-
-    return f
-
-
-def calculateFMS(A: anndata.AnnData, B: anndata.AnnData):
+def calculateFMS(A: anndata.AnnData, B: anndata.AnnData) -> float:
     """Calculate Factor Match Score (FMS) between two RISE decompositions.
 
     Factor Match Score measures the similarity between two tensor decompositions
@@ -80,7 +63,7 @@ def calculateFMS(A: anndata.AnnData, B: anndata.AnnData):
         )
     )
 
-    return fms(A_CP, B_CP, consider_weights=False, skip_mode=1)  # type: ignore
+    return fms(A_CP, B_CP, consider_weights=False, skip_mode=1)
 
 
 def plot_fms_percent_drop(

RISE/tests/test_parafac2.py

@@ -5,8 +5,9 @@
 import numpy as np
 import pandas as pd
 
+from analysis.imports import import_thomson
+
 from ..factorization import pf2, rise_pca_r2x
-from ..imports import import_thomson
 
 
 def test_factor_thomson_reprod():

analysis/__init__.py

@@ -0,0 +1,3 @@
+"""Analysis package containing paper-specific datasets, gating, regression,
+and figure-generation scripts.
+"""

analysis/figures/__init__.py

@@ -0,0 +1 @@
+"""Figures subpackage for paper analysis."""

RISE/figures/common.py → analysis/figures/common.py

@@ -68,7 +68,7 @@ def genFigure():
     start = time.time()
     nameOut = "figure" + sys.argv[1]
 
-    exec(f"from RISE.figures.{nameOut} import makeFigure", globals())
+    exec(f"from analysis.figures.{nameOut} import makeFigure", globals())
     ff = makeFigure()  # type: ignore # noqa: F821
 
     if ff is not None:

analysis/figures/commonFuncs/__init__.py

@@ -0,0 +1 @@
+"""Common functions and plotting utilities for figures."""

RISE/figures/commonFuncs/plotLupus.py → analysis/figures/commonFuncs/plotLupus.py

@@ -5,12 +5,12 @@
 from matplotlib.axes import Axes
 from sklearn.metrics import RocCurveDisplay
 
-from ...logisticReg import (
+from analysis.logisticReg import (
     predaccuracy_lupus,
     predaccuracy_ranks_lupus,
     roc_lupus_fourtbatch,
 )
-from ..commonFuncs.plotGeneral import cell_count_perc_lupus_df
+from RISE.plotting import cell_count_perc_lupus_df
 
 
 def samples_only_lupus(X: anndata.AnnData):

RISE/figures/figure2b_e.py → analysis/figures/figure2b_e.py

@@ -7,8 +7,9 @@
 import pacmap
 from sklearn.decomposition import PCA
 
+from RISE.plotting import plot_labels_pacmap
+
 from .common import getSetup, subplotLabel
-from .commonFuncs.plotPaCMAP import plot_labels_pacmap
 
 
 def makeFigure():

RISE/figures/figure2f_h.py → analysis/figures/figure2f_h.py

@@ -6,13 +6,14 @@
 import numpy as np
 import pandas as pd
 
-from .common import getSetup, subplotLabel
-from .commonFuncs.plotFactors import (
+from RISE.plotting import (
     plot_condition_factors,
     plot_eigenstate_factors,
     plot_gene_factors,
 )
 
+from .common import getSetup, subplotLabel
+
 
 def makeFigure():
     """Get a list of the axis objects and create a figure."""

RISE/figures/figure3.py → analysis/figures/figure3.py

@@ -8,22 +8,22 @@
 import seaborn as sns
 from matplotlib.axes import Axes
 
-from .common import (
-    getSetup,
-    subplotLabel,
-)
-from .commonFuncs.plotFactors import (
-    plot_condition_factors,
-    plot_gene_factors,
-)
-from .commonFuncs.plotGeneral import (
+from RISE.plotting import (
     cell_count_perc_df,
     gene_plot_cells,
     plot_avegene_per_category,
     plot_avegene_per_celltype,
     plot_cell_gene_corr,
+    plot_condition_factors,
+    plot_gene_factors,
+    plot_labels_pacmap,
+    plot_wp_pacmap,
+)
+
+from .common import (
+    getSetup,
+    subplotLabel,
 )
-from .commonFuncs.plotPaCMAP import plot_labels_pacmap, plot_wp_pacmap
 from .figure2f_h import groupDrugs
 
 

RISE/figures/figure4d.py → analysis/figures/figure4d.py

@@ -12,10 +12,11 @@
 from sklearn import preprocessing
 from sklearn.linear_model import LogisticRegression
 
-from ..factorization import correct_conditions
+from analysis.figures.commonFuncs.plotLupus import samples_only_lupus
+from RISE.factorization import correct_conditions
+from RISE.plotting import rotate_xaxis, rotate_yaxis
+
 from .common import getSetup, subplotLabel
-from .commonFuncs.plotGeneral import rotate_xaxis, rotate_yaxis
-from .commonFuncs.plotLupus import samples_only_lupus
 
 
 def makeFigure():