Disentangling Interpretable Factors with Supervised Independent Subspace Principal Component Analysis
Authors: Jiayu Su
David A. Knowles
Raul Rabadan
Venue: NeurIPS 2024 poster (2024)
Keywords: Principal Component Analysis (PCA)
Disentanglement Representation Learning
Computational Biology
Primary area: machine_learning_for_healthcare
TL;DR: A multi-subspace extension of PCA to disentangle interpretable subspaces of variations with supervision.

Abstract: The success of machine learning models relies heavily on effectively representing high-dimensional data. However, ensuring data representations capture human-understandable concepts remains difficult, often requiring the incorporation of prior knowledge and decomposition of data into multiple subspaces. Traditional linear methods fall short in modeling more than one space, while more expressive deep learning approaches lack interpretability. Here, we introduce Supervised Independent Subspace Principal Component Analysis ($\texttt{sisPCA}$), a PCA extension designed for multi-subspace learning. Leveraging the Hilbert-Schmidt Independence Criterion (HSIC), $\texttt{sisPCA}$ incorporates supervision and simultaneously ensures subspace disentanglement. We demonstrate $\texttt{sisPCA}$'s connections with autoencoders and regularized linear regression and showcase its ability to identify and separate hidden data structures through extensive applications, including breast cancer diagnosis from image features, learning aging-associated DNA methylation changes, and single-cell analysis of malaria infection. Our results reveal distinct functional pathways associated with malaria colonization, underscoring the essentiality of explainable representation in high-dimensional data analysis.

Disentangling Interpretable Factors with Supervised Independent Subspace Principal Component Analysis
Jiayu Su 12.5t David A. Knowles2:4.5t 2.4.tt Raul Rabadan 1,2,31
Program for Mathematical Genomics; *Departeent of Systems Biology. Columbia University DDepartment of Biomedical Informatics, Columbia University
*Department of Computer Science, Columbia University New York Genome Center i Correspondence to [js5756, rr2579) columbia. dak2173@colume edu
Abstract
The success machine learning models relies heavily on effectively representing high-dimensional data. However, ensuring data representations capture human- understandable concepts remains difficult, often requiring the incorporation prior knowledge and decomposition of data into multiple subspaces. Traditional linear methods fall short in modeling more than one space, while more expressive deep learning approaches lack interpretability. Here, we introduce Supervised Independent Subspace Principal Component Analysis (sisPCA), a PCA extension designed for multi-subspace learning. Leveraging the filbert-Schmidt Indepen- dence Criterion (HSIC), sisPCA incorporates supervision and simultaneously ensures subspace disentanglement. We demonstrate sisPCA's connections with autoencoders and regularized linear regression and showcase its ability to identify and separate hidden data structures through extensive applications, including breast cancer diagnosis from image features, learning aging aging-associated DNA methyla- tion changes, and single-cell analysis of malaria infection. Our results reveal distinct functional pathways associated with malaria colonization, underscoring the essentiality explainable representation high-dimensional data analysis.
Introduction
High-dimensional data generated by complex biological mechanisms encapsulate an ensemble of patterns. A prime example is single-cell RNA sequencing (scRNA-seq) data (Fig. [ These datasets measure the expression of tens of thousands of genes across potentially millions of cells, creating rich tapestries woven from interacting cellular pathways, gene dynamics, cell states, and inherent measurement noise. To unravel these patterns and reveal hidden relationships, it is necessary to decompose the data into meaningful, lower-dimensional subspaces.
Linear representation learning methods, such Principal Component Analysis (PCA) (Hotelling 1933) and Independent Component Analysis (ICA) [Comon] 1994), extract latent spaces from data using explainable linear transformations. These widely employed unsupervised tools learn a single latent space or a union of one-dimensional subspaces. Independent Subspace Analysis (ISA) extends ICA by extracting multidimensional components as independent subspaces [Cardoso 1998 Yet. the unsupervised nature of these methods precludes knowledge integration, restricting their utility and sometimes even identifiability [Theis] 2006). the context of Fig. [ the representation learned without supervision fails to separate temporal variability from technical batch effects.
Conversely, recent advancements in deep generative models, especially semi-supervised approaches Kingma et al. 20144), have shown promise in disentangling diverse latent spaces and retaining relevant information under supervision. However, challenges remain in ensuring subspace independence within variational autoencoder (VAE). While models like 3-VAE [Higgins et al. 2016), HCV
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

Temporaldynamiic
Sequencing batch
Cell intrinsic quality
M tochandial counts percentage
Pact-iifacion time point
Experiment
Eentrll hi
Aug20
Feb20
40
24
Jan22
hgl
NeviB
20
Figure 1: Example scRNA-seq dataset from Afriat al. [022]] Each dot represents the gene expression vector F E 8R23003 of a cell, visualized in 2D and colored by cell properties (Ym). Variability in the dataset X arises from multiple sources: (left to right) temporal dynamics of infection, technical batch effects, and cell quality. Incorporating supervisory information Y, such as time points, allows for the extraction of patterns in distinct subspaces {Zm, that correspond to different sources of variability. Moreover, the linear mapping {Um X + Zm} directly quantifies the relationship between gene expression and the property of interest, enabling discoveries such as the identification of genes underlying the persistent defense against infection. The disentanglement is particularly important ensure minimal confounding effects. See Section 4.4 for details.
Lopez et al., 2018] and biolord [Piran et al. 2024 attempt to address this, inference for deep generative models remains challenging and the learned representations are not interpretable.
To bridge the gap, we propose Supervised Independent Subspace Principal Component Analysis (sisPCA | an innovative method extending PCA multiple subspaces. By incorporating the Hilbert- Schmidt Independence Criterion (HSIC), sisPCA effectively decomposes data into explainable independent subspaces that align with target supervision. It thus reconciles the simplicity and clarity of linear methods with the nuanced spute modeling capabilities advanced generative models. In summary, contributions include:
A multi- multi-subspace extension of PCA for disentangling linear latent subspaces in high- dimensional data. We additionally show that supervising subspaces with a linear target kernel can be conceptualized linear regression regularized akin Zellner's g-prior.
alternating optimization algorithm to compute SisPCA subspaces. The learning process with linear kernels resembles matrix factorization, which may potentially benefit from desirable local geometric properties (Conjecture 3.1
An efficient eigendecomppsition-basse
Demonstrated effectiveness and interpretability of sisPCA in various applications. This includes identifying diagnostic image features for breast cancer, dissecting aging signatures in human DNA methylation data, and unraveling inne-independent transcriptomic changes in mouse liver cells upon malaria infection.
2 Background
2.1 Hilbert-Schmidt Independence Criterion (HSIC)
The Hilbert-Schmidt Independence Criterion (HSIC) serves as a methodology for testing the indepen- dence of two random variables X and Y [Gretton et al. 2005a]. It operates by embedding probability distributions into reproducing kernel Hilbert spaces (RKHS) and quantifying independence by dis- tance between the joint distribution and the product of its marginals. Specifically, HSIC builds on the cross-covariance operator CXY G - F, which is generalization of the cross-cor covericov matrix to RKHS (i.e., the feature spaces of X and Y after transformation).
CXY == Exy (ø(X) @ (Y)] - Ex [x(X] @ [y(v(Y)]-
Here :=
= span((o( and G == span ( (Y))) are the RKHSS with feature maps ø and t for X and respectively. The tensor product operator f so g maps G F, such that (f 8 g)h fg, hig for all h E G. The HSIC is then defined as,
HSIC(X,Y:
= ||xx Ilus ((Cry u2)
A Python implementation of SisPCA is available on GitHub at https://github..c http:/ghwww.ooo/iaauuuuuuu/ssaie (DOI 10.5281/zemodo. wtw.ood........... 13932660). The repository also includes notebooks to reproduce results in this paper.

---

Goals
HSIC Interpretation of PCA and its extensions
Samples (0)
PCA
Supervised PCA (sPCA)
Y
Features (s)
Supeniisions
Linear
projections (U,, Us...
Superviedd Independeent Subsspaces
Supervised Independent Subspace PCA (sisPCA)
Samplee (1)
Zy
Za
HSIC maximization
|-. Subepacce expansion
Embeddings (d., do...)
HSIC mimizatioon
(intrinsc)
Figure 2: Overview of SisPCA and its relationship with other PCA models.
where o and My are the orthogonal bases of 9 and F. and Y are independent and only HSIC(X, FF 0. In practice, given finite sample ((1s, b)) from the joint distribution Px,y, the empirical HSIC can be computed as,
HSIC,(X,Y) = (n- 1)2 tr (x
where Kx and Ly are matrices of kernel evaluations over the samples [m in F and {w:) in g, and H the centering matrix defined as = In 1,1,, with I, being the identity matrix.
2.2 Related work
The HSIC, as non-parametric criterion for independence, has been effectively integrated into numerous representation learning models, especially for disentangling complex data structures. One of the earliest applications of HSIC is ICA [Comon 1994), where the goal is to recover unmixed and statistically independent sources. Gretton et al. [005a] showed that minimizing HSIC via gradient descent outperforms specialized linear ICA algorithms. Alternatively, HSIC can also be maximized to encode specific information learning tasks. supervised PCA (sPCA), HSIC is deployed to guide the identification of principal subspaces with maximum dependence on target variables (Barshan al. 201T]. More recently, Ma al. (2020) used HSTC as a supervised learning objective for deep neural networks to bypass back proppaappagg,onn and Li al. 2021 subsequently extend it to the context of self-supervised learning.
Our primary interest lies in extending these models to identify multiple subspaces, each representing independent, meaningful signatures. this direction, Cao et al. (2015) suggested the inclusion of HSIC as a diversity term in multi-view subspace clustering, encouraging representations from different views to capture complementary information. Building on similar idea, Lopez al. [2018] incorporated HSIC-based regularization into VAE architectures to promote subspace separation. However, there has yet to linear multi-space model for more interpretable data decomposition.
The presented work is also related to contrastive representation learning, introduced in Abid et al. (2018) and|Abid and Zou [2019] Along this line of research, recent studies have also applied HSIC to regularize contrastive VAE subspaces (Tu et al. 2024 Qiu et al. 2023]. See Appendix A for more discussions connections and differences.
3 The sisPCA model
We introduce Supervised Independent Subspace Principal Component Analysis (sisPCA), linear model for disentangling independent data variation (Fig. 2 The model's linearity ensures explicit interpretability and enables regression-batte based extensions such as sparse feature selection, We formally discuss the connection between sisPCA and regularized linear regression in Section 3.2
3.1 Problem formulation
As motivated in Fig- associated target variables {Yr ym mmnn JeE, we aim to find separate subspace representations of
given a dataset (X, E XP), with n observations and p features and m

---

the data {[2] C CR 1/ dependence with one target variable while minimizing dependence with other subspaces. For simplicity, we assume Fuclidean space, {X,} == X E xR= and (Y:) == E Ruxm Other types of data can be easily handled using appropriate kernels, which will be discussed later. The linear projection to the j-th subspace U, R + R& is represented by the matrix Un E RRxx, with 2j : XU; E R being the data representation in this j-th subspace. Similar to the concept of PCA loading, the projection Uj depicts linear combinations of original features in and thus can be directly interpreted as feature importance scores.
with latent dimensions {d, [d, Each subspace should maximize
Our overall objective find the set of subspace projections {U1, Um that solves
argmax Σ II(U,,Y)) - A ΣE X(XU,. Us...,, U
j=1
j=
under some constraints. Here I is measure of dependence and À penalizes overlapping subspaces. Using HSIC as the dependence measure, sisPCA solves the constrained optimization
argmax U Um
Σ tr(K2, HKy, H) X ΣE j=1
tr tr(K2, HK, H)
(I)
subject to UTUU=1, Vj E m)
where K2, and Ky; kernels defined the j-th subspace 2j and the j-th target variable } respect. This formulation differs from kernel PCA, where in the first term of (I) the kernel is defined over Z rather than X (the kernel extension of sisPCA will be discussed later).
In the special case where K2, = =2,4 is linear for subspaces, the optimization becomes
argmax U Un j=1
Σ UXUXXU HKy, H)
to (XXUXX
(2)
subject to U3U,- = 1, Vj E {I,..., m).
The first term is the supervised PCA objective [Barshan et al. 2011]. We can thus view the above formulation, termed SisPCA- -innear as an extension of supervised PCA to multiple subspaces with additional regularization for subspace independence (Fig. 2 It comes with an appealing property:
Remark 3.1. Maximizing the sisPCA- 1 linear objective 2 is equivalent to minimizing the recon- struction error of a linear autoencoder plus regularization. See Appendix B
We now examine the second HSIC term in (2). Consider two subspaces Z, XU E Ruxd, and Z. = XV E Rnnx with centered X. The HSIC regularization is tr(XUU X" XVVT x tr tr( 7,
Z Z. 11 2 0.
This term equals zero if and only if 2u and Zu are orthogonal. While not convex in U and V jointly, the coupling of U and V solely through the matrix product 2 Z. indicates a well-behaved local geometry, as explored by [Sun and Luo (2016). Indeed, it allows SisPCA- theoretical insights function ff over the matrix UTV in our context f = |2 2. globally optimal under mild conditions achievable through proper regularization.
to benefit from matrix factorization. For example, Ge et al. (017)) showed that for quadratic all local minima are also
Returning to (2), we see that spurious local optima may only emerge from subspace imbalance the first symmetry-breaking supervision term, where some subspaces may contribute more to the overall objective. This leads to the following conjecture on the optimization landscape of (2):
Conjecture 3.1
(informal): The sisPCA- - linear objective (2) has no spurious local optima under balanced supervision, which is achievable through proper regularization; With unbalanced supervision, the global optima can still be recovered using local search algorithms following simple initialization based on the relative supervision strength of each subspace.
Intuitively, the balanced supervision condition is to ensure that the local optima induced by subspace symmetry and interchangeability are also global optima. For more discussions the optimization landscape, the balance condition, and Conjecture 3.1 see Appendix C
We now consider an iterative optimization approach to solve 2

---

Remark 3.2. The optimization problem in 2 can be solved via alternating optimization, where each iteration has an analytical update for each subspace.
Using basic matrix algebra, the objective of simplifies to,
m
m
m
m
Lr Ir(Kz, HK = Σ 2 (((U 1=1
2
==1
j=1 i
where K3 = H(K,, E K2 and Kz, = Z, XU.U? X Given the set
{U U, can be updated by maximizing tr(UT X XUJU), leading to the update (/((+)) t Qd,- where Qa, are the columns of Q corresponding the d, largest eigenvalues from the eigendecompo- sition X² == QAQ".
The full optimization process is outlined in Algorithm 1 Appendix B convergence to an optimum as the objective bounded and non-decreasing in every iteration. We further implement initialization step find the path (subspace update order) towards the global optimum as proposed in Conjecture 3.1 Briefly, we compare subspace contributions to the supervision loss and prioritize updates for subspaces under stronger supervision.
This procedure guarantees
While convenient, a zero HSIC regularization loss with a linear kernel independent subspaces. For strict independence, the subspace kernel Kz in needs to be universal (e.g., Gaussian). We refer to this as SisPCA- sisPCA-general and solve using gradient descent (Algorithm 2 Appendix D). The naive implementation with Gaussian kernel has complexity O(n²) in contrast to O(?²) with a linear kermel. See Appendix D for performance difference discussions.
does not guarantee
Both SisPCA- linear and -eeeralleerral are linear methods where {Um] } measures direct con- tribution of original features to each subspace. Nevertheless, the sisPCA framework can easily extended to incorporate nonlinear feature interactions, analogous to kernel PCA. See Appendix E
3.2 Kernel selection for different target variables
The of kernel independence measure in I allows SisPCA to accommodate various data types through flexible kernel choices for Kx (data), Kz (latent subspace), and Ky (target).
For categorical variables Y ji=1 (e.g., cancer types), we use the Dirac delta kernel, Ky (i, = _y,Y,,
It is also possible to other general graph kernels for categorical variables with an intrinsic hierarchical structure, c.g., subtypes and stages [Smola and Kondor 20033.
For continuous variables Y E RR", we use the linear kernel
Ky = YYT.
When Y e xn is multivariate, the kernel is the sum of der-dimension kernels Ky Ei Y.Y."" Remark 3.3. Maximizing the sisPCA CeaPCA-linnar linear objective (2) with linear kernels on the target space is equivalent to performing regularized regression against the target. See Appendix B
In a nutshell, we show that SisPCA- -inear can be viewed as approximating the target Y with Z = Xu. The particular regularization on u corresponds to zero-mean multivariate Gaussian prior, which is related to Zellner's g-prior in Bayesian regression.
3.3 Learning an unknown subspace without supervision
In practical applications, common goal would be to recover both subspaces linked to known attributes (supervised) and to unknown attributes (unsupervised) simultancously. sisPCA, this is achieved by setting the target kernel for the unknown subspace to the identity matrix (Ky = I),
corresponding to unsupervised PCA (Fig. 2 Section 4.1 provides an example this process.
However, the absence of external supervision introduces a potential identifiability issue. As indicated Appendix eq. eq. eeach supervised subspace is driven by two forces to (1) align with the target and (2) capture major variations in the data. This dual objective can lead to scenarios where unknown attributes, ideally retained the residual unsupervised subspace, being inadvertently presented in supervised subspaces. Fig. in Appendix D gives an example where sisPCA -arcCA--eeeaall fails.
5

---

Calagorical
51: Calegonical
Calagoniaa
51: Categoricl
51: Categonicl
X
PN
Grout
Group
PCLIABTCA)
susrccl
52: Continuuus
52: Comtnnus
52: Continous
$2: Continuuus
52: Continuuus
Srapp
Unsupervised
53: Unsuparvised
53: Unsuparvisec
$3. Unsupeniised
S2: Unsuparvised
RCCMMPCAL
Group
F2: GEPC4)
PCLISSPCA
(OPPC)
(a) Ground truth
sPCA
sisPCA-linear (A
Figure 3: Example application of recovering latent space with three subspaces (rows in panel a) embedded in a high-dimensional space. The first subspaces (rows) of sPCA (panel b) and sisPCA (panel c) are supervised by the corresponding target variables.
Applications
Unless otherwise specified, in the following sections sisPCA refers to sisPCA-linear, where the objective (2) is solved using Algorithm| For baseline comparisons, consider linear models including PCA and sPCA. While non-linear VAE counterparts such as HCV (Lopez et al.) (2018 are included for quantitative performance benchmarking in Section 4.4 they are not considered for interpretability analyses due to their inherent complexity. The rationale and details for baseline selection are provided in Appendix F
4.1 Recovering supervised and unsupervised subspaces in simulated data
We first consider learning latent subspaces associated with known and unknown attributes using simulated data. The dataset reflects ground truth 6-dimensional latent space, comprising three distinct 2D subspaces (Fig. 3a) S1 with two Gaussian distributions, S2 with noisy 2D grid and S3 with a ring structure. The defining manifold characteristics á S3 remain unknown to the model, representing the unsupervised component. These subspaces were concatenated and linearly projected to dimensinss using a x 20 matrix with entries uniformly distributed on [0, 1].
Both unsupervised PCA (targeting S3) and supervised PCA (S1 and S2) subspaces capture structures heavily influenced by S2 (Fig. 3h and Fig. 11 in Appendix [H]. due S2's pronounced variations. In contrast, sisPCA markedly improves the disentanglement of these subspaces, especially the unsuper- vised S3 (Fig. (..) Despite S2's dominant influence, sisPCA isolates the effects of each subspace, resulting in clearer separation of the three independent signals. The two supervised subspaces S1 and S2 are distinctly characterized by patterns exclusively associated with the supervision attributes. In the unsupervised subspace S3, although Principal Component 2 (PC2) picks up some categorical information from SI, SisPCA successfully uncovers the underlying circular structure.
4.2 Learning diagnostic subspaces from breast cancer image features
We apply sisPCA to the Kaggle Breast Cancer Wisconsin Date demonstrate its utility in data compression and feature extraction. The dataset contains 569 samples with 30 summary features from breast mass imaging. Our goals ("Malignant' or 'Benign', agnostic during training), and (2) understand the relationship between original features and how they contribute the learned representation and diagnosis potential.
(1) learn compressed subspaces for predicting disease status
CC BY-NC-SA 4.0 license.

---

Symmety SSitouste 0408)
PCA IS ihouette 3.234)
PC1 loading
(Sinouste 0473)
(Siiroust 0.518
Bouest -8034)
800
5.0
Deggoosh
25
Berign
-10
-25
PCS
PC1
PC1
(a) PCA
(b) sPCA
(c) sisPCA (A = 10)
Figure Feature extraction on the breast cancer dataset. The two top PC1 contributors in PCA (panel a) are used as supervisions to construct the 'radius' and 'symmetry' subspaces (panel b and c).
The diagnosis label, invisible to all models, is used to measure subspace quality via the mean silhouette score n. E, sii = :/2 (bi) - a(i)) max b(i),, where a(i) and b(i) are mean intra-cluster and nearest-cluster distances for sample i. Higher scores indicate larger diagnostic
potential. In addition, we use the Geodesic Grassmann distance d( (i, 2j) = (E 0% to
1/2
measure subspace separateness, where =
Zi, dim Zj} and {6,.} the principal angles [Mio and Isn-ls 1992]. Higher scores indicate better disentanglement. Inputs for all model are zero-centered and geriance-standardin
In the PCA space, samples are well- well-separated based on diagnosis along PCI (Fig. try_mean" and 'radius madius_meann are the top two features negatively contributing to PCI, motivating us to construct separate subspaces to reflect nuclei size (using 'radius mean' and 'radius sd' as targets) and shape (using 'symmetry mean' and 'symmetry_ss __sd' as targets). The remaining 26 features are projected onto these subspaces using sPCA (Fig. 4b) and sisPCA (Fig. 4b). In sPCA, both subspaces better explain diagnosis status than PCA but remain highly entangled (Grassmann distance: 1.493, Pearson correlation of PC2 loading: 0.850). However, with SisPCA's explicit disentanglement, the symmetry subspace loses its predictive power as the two spaces separate further (Grassmann dis- tance: 2.710). SisPCA subspaces are constructed from distinet feature sets (PC2 loading correlation -0.203), with 'area' and 'perimeter' contributing more to the radius subspace and 'compactness' and 'smoothness' the symmetry one. The radius subspace also gains additional predictive power by repulsing further from the symmetry space (Silhouette scores: 0.516 in SisPCA, 0.470 in sPCA).
"symme-
Our SisPCA results suggest that cell nuclear size is more informative for breast cancer diagnosis than nuclear shape. We confirm this by measuring directly the predictive potential of target variables (Silhouette scores: 0.457 for 'radius meaasmmman and "radius_ssi sd', 0.092 for 'symmetry mean' and "sym- metry_sd'). Our conclusion also aligns with previous clinical observations (Kashyap et al. 2018). In contrast, PCA and sPCA, while capable extract new diagnostic features (PC1), cannot faithfully capture feature relationships without disentanglement and potentially overestimate symmetry-related features' relevance in malignancy.
4.3 Separating aging-dep DNA methylation changes from tumorigenic signatures
Tumorigenesis and aging are two intricately linked biological processes, resulting in cancer omics data that often display patterns of both. DNA methylation (DNAm) exemplifies this complexity, undergoing genome- wide alterations during aging while also exhibiting cancer-specte changes in particular regions, presumably silencing tumor suppressor genes activating oncogenes.
4.3.1 Problem and dataset description
The Cancer Genome Atlas (TCGA offers a comprehensive collection of DNAm datasets from patients with various cancer types. In TCGA DNAm data, methylation status is probed across genomic locations using the Illumina Infinium 450K array and quantified as beta values ranging from 0 (unmethylated) to (methylated). The resulting data matrix X presents challenges due to the use of three different probe types and highly correlated features, typically requiring careful preprocessing. For illustration purpose, we use downsampled dataset comprising the first 5,000 non-constant and non-NA CpG sites from 9,725 TCGA tumor samples across 33 cancer types.
Data access is controlled through the NCI Genomic Data Commons (ttps::/portal.ddc.cmmeeee
7

---

Our goal is to disentangle tumorigenic signatures from age-dependent methylation dynamics in this pan-cancer DNAm data. Traditional methylation analyses often employ regression- based methods to learn site-specific statistics, later aggregated by gene high-level genomic annotations [Bock] 2012]. However, these approaches may suffer from high dimensionality and multi-collineariyy due to potential redundancy in CpG methylation activity. We propose using SisPCA to address these limitations, which allows us to (1) learn compressed, low-dimensional representations that retain biological information, and (2) minimize confounding factors not controlled in simple regression models by enforcing disentanglement.
Specifically, we aim to learn two subspaces: one aligning with chronological age (CA, aging subspace, supervised with linear kernel) and another with TCGA cancer categories (cancer subspace, supervised with a delta kernel). We evaluate the quality of learned representations using information density measured by the Silhouette score and subspace separateness by the Grassmann distance, For the rank-one aging subspace (rank(Ky) rank(YY = 1), we also measure information density using the maximum absolute Spearman correlation, maxae[1.4) ( (l((Z))) each axis of the subspace Zj E Rxx,,
CA)]]; between CA and
4.3.2 Quantitative performance of subspace quality
We first validate the disentanglement effect of formulation 2 Unlike models such as Lopez et al. 20188] that use Gaussian kemel, SisPCA minimizes the HSIC regularization with linear kernel. This approach trades strict statistical guarantees for improved computational efficiency (detailed in Appendix D). Our experiments demonstrate that as À increases, the two subspaces show increasing divergence (Table [T. Notably, while HSIC-Gaussian is not explicitly optimized, decreases in tandem with HSIC-Linear. The generally larger values of HSIC-linear potentially offer advantages in optimization and help mitigate numerical rounding errors.
Table 1: Separateness of the aging and cancer subspaces inferred by sisPCA.
sisPCA
PCA
À = 0 (sPCA)
A A=
X 10
HSIC-Linear (in the objective of sisPCA) USIC-Gaussian Grassmann distance
481.6
189.7
2.0c-4
2.4e-05
1.le-2
6.5e-3
7.0e-4
7.0e-4
0
3,09
4.97
4.97
Furthermore, the separation of aging and cancer subspaces leads to moderate increase in target information density and a decrease in confounding information (Table 2}. However, stronger regular- ization does not always equate to better representations. This is partly due to the inherent coupling between aging and tumorigenesis. Efforts to remove aging signals inevitably result in information loss on cancer type. We discuss the tuning of À more generally in Appendix G
Table 2: Information density in each SisPCA subspace.
sisPCA
Subspace
PCA
A
10
Maximum Spearman correlation with age
age cancer
0.213
0.278
0.286
0.294
0.213
0.233
0.103
0.115
Silhouette score with cancer type
age 0.074
-0.183
-0.221
-0.230
cancer
0.074
0.106
0.107
0.097
4.4 Disentangling infection-induced changes in the mouse single-cell atlas of the Plasmodium liver stage
Malaria, transmitted by mosquitoes carrying the Plasmodium parasite, involves critical liver stage where the parasite colonizes and replicates within host hepatocytes. This section examines the intricate host-parasite interactions at single-cell resolution during this stage.

---

(5ihoucti (144))
PCA
Infection (Sihosete 00515
Time 15 Incueto 0710
Infecto d5 Woouete -..01
di illouute MM99)
par FALSE
Infection (5thouete -- 097)
Time iS iosuett a 3151)
Infectio iSitouette 0055)
Time houre 0 582)
[Sinnueee 0,111]
Fost-inlecion time point
Control
24 Tpl
35
(a) PCA
(b) sisPCA (A = 10)
(c) hsVAE (A = 100)
Figure 5: UMAP visualizations of scRNA-seq data. Each column shows different learned subspace: (a) PCA, (b) ALsPCA-infection and SisPCA-time, and (c) hs VsV---infection and VAE--timme See Fig. 12 for other models. Cells are colored by either infection status (top row) or post-infection time (bottom row). In an optimal pair of subspaces, each property (infection status or time) should be more distinguishable in its corresponding subspace while showing less separation in the other.
4.4.1 Problem and dataset description
We analyze ScRNA-seq data of mouse hepatocytes from Afriat et al. 2022
Our goal is to distinguish genes directly involved in parasite harboring from those associated with broader temporal changes post-infection. The processed dataset comprises gene expression profiles of 19,053 cells collected at five post-infection time points (2, 12, 24, 30, and 36h) and from control mice. Infection status was determined based on GFP expression linked to malaria. We keep the top 2,000 highly variable genes and use normalized expression as model inputs. Time points are treated discrete categories to account for potential nonlinear dynamics. See Appendix FF for full experiment details.
4.4.2 Learning the infection subspace associated with parasite encapsulation
UMAP visualizations of PCA, sPCA, and sisPCA subspaces show that while PCA primarily captures temporal variations (Fig. 5h. both sPCA and SisPCA successfully differentiate between infected and uninfected hepatocytes their infection subspaces. However, sPCA's infection space still exhibits significant temporal effects, suggesting uncontrolled confounding effects (Fig. [2p in Appendix H... sisPCA effectively eliminates this intermingling, yielding cleaner representations where relevant biological information is further enriched in the corresponding spaces (Fig. 5p).
Comparisons with non-linear VAE counterparts (Fig. 5 and Fig. 12 full model description in Appendix F) and quantitative evaluations (Table 3 demonstrate that our HSIC-based supervision formulation
(T) achieves performance comparable to neural network predictors. Notably, sisPCA outperforms HSIC-constraing and uninfected cells. Indeed, hs VAE's performance in the infection subspace is so poor that even under supervision, the representation contains near-random information (Fig. 5b). To address this gap, we developed hs VAE-sc by incorporating additional domain-specific knowledge (See Appendix [F] and Fig. 12 in Appendix [H). This model learns a much improved infection space (Silhouette score: 0.233) while maintaining high distinguishability in the temporal space (Silhouette score: 0.634). The improvement is likely due to the fact that parasite encapsulation
supervised VAE (hs (AE) et al.) 20088) in separating infected
induce slight increase in total RNA counts; by using unnormalized count-level data and explicitly modeling the library size, the model captures this additional information and thus produces enhanced results.
While visualizations and quantitative metrics provide useful estimates, they may not fully capture the biological relevance of the representations. For instance, a subspace with two point masses perfectly separating cells by infection status would have the highest information density, yet offer little new biological insight. Advantageously, linear models like sisPCA are inherently interpretable. Using the learned sisPCA projection U as feature importance scores, we rank genes based on their PC contributions. Chemokine genes, including Cxcl10, emerge as top positive contributors to infection- PCI,, which is elevated in infected cells. This aligns with their established role as acute-
Preppocessed data available http:://bb..ogg/10.00////// www.:iiasaaeeeeeeie.ee. (CC BY 4.0).

---

Table Quantitative evaluation of subspace representation quality.
Linear
Non-linear
Subspace
PCA
sPCA
sisPCA
VAE
SupVAE
hsVAE
Grassmann distance Silhouette infection
3.771
4.824
3.571
3.598
infection time
-0.014
0.207
0.235
-0.036
0.041
-0.015
-0.014
-0.075
-0.097
-0.036
-0.087
-0.069
Silhouette time point
infection
0.311
-0.029
-0.028
0.296
0.052
0,015
time
0.311
0.348
0.355
0.296
0.479
0.582
phase response markers. Gene Ontology (GO) enrichment analysis of genes with significant PC1 loading scores reveals that infection leads to reduced fatty acid metabolism and enhanced stress and defense responses, consistent with known Plasmodium harboring effects. These results remain highly consistent across wide range of hyperparameters, demonstrating the robustness of our approach. See Appendix G for full examination on the effect of A.
Top 100 positive contributors Inf
Top 100 negative contributors to Infection-PC1 flower activity infected cells)
[higher activity in infected colls)
Cytdine-meddtte Response vine
Ferty acid metebolic process
pathway
Response to melezu bacaenil oro on
P.adjust
serodd mebsbole Organic acis brosytthete process
P.adjust
Calluar response melecure aacteaaa
cel ular response to blods samu
Defence resgonse symbiont
ipopelyaac harktee
25
Carboryki and hessytheee peocess
Carbooficc ac cabablic piocess Kenobicti metabolic piocess Response metabolic
Lit 25
220-28
Ce Regulation resporic lo hrlt samutus Process susophaglc Autophagy
S6s
Celitar response to senotiote stimulas
0.00 0:0 0.10 0.15
0.00
0.05
0.30 015 0.20
Gene rado
Gene 1200
Figure 6: GO biological process enrichment results of top genes contributing to the sisPCA-infection subspace. Genes are ranked by their PCI loading and are grouped by effect direction.
Discussion
This study presents SisPCA, novel extension of PCA for disentangling multiple subspaces. We showcase its capability in interpretable analyses for learning complex biological patterns, such as aging dynamics methylation and transcriptomic changes during Plasmodium infection. To enhance usability, we have implemented an automatic hyperparameter tuning pipeline for À using grid search and spectral clustering, similar to contrastive PCA [Abdd et al. 2018 (Appendix G
The linear nature of SisPCA may miss non-linear feature interactions. potentially underperforming on more complicated datasets. While nonlinear extensions are possible (Appendix [D. they come at the cost of reduced computational efficiency and interpretability. Linear kernel HSIC: The HSIC-linear regularization, while compu- tationally convenient, does not guarantee complete subspace independence. However, minimizing HSIC-linear tends reduce HSIC with Gaussian kernel (T) which suggests that the issue is less of a concern in practice. Subspace identifiability: Our formulation (I) relies on external supervision to differentiate subspaces, which could lead to identifiability issues the supervisions are too similar, or when one subspace is unsupervised.
Still, sisPCA has several limitations: Linearity constraints:
Despite these limitations, sisPCA's 's ability provide interpretable, disentangled representations of complex biological data makes it a valuable addition to the toolkit of biologists working with high- dimensional datasets. We envision future applications on larger-scale omics datasets and potential new biomedical discoveries.
Acknowledgements
thank Bianca Dumitrascu for early discussions and anonymous reviewers for their valuable comments. This work was funded by the National Institutes of Health, National Cancer Institute (R35CA253126, U0ICA261822, and U0ICA243073 R. R.) and the Edward Evans Center for MDS Columbia University (to J. S. and R. R.).
10

---

References
Amichay Afriat, Vanessa Luzsrte-LLufs, Keren Bahar Halpern, Lisa Buchauer, Sofia Marques, An- gelo Ferreira Chora, Aparajita Lahree, Ido Amit, Maria M Mota, and Shalev Itzkovitz. A spa- tiotemporally resolved single-cell atlas of the plasmodium liver stage. Nature, 61 611(7936):663-565 2022.
Harold Hotelling. Analysis of a complex of statistical variables into principal components. Journal of educational psychology, 24(6):417, 1933.
Pierre Comon. Independent component analysis, new concept? Signal processing, 36(3):287-314, 1994.
J-F Cardoso. Multidimensional independent component analysis. In Proceedings of the 1998 IEEE International Conference on Acoustics, Speech and Signal Processing, ICASSP'98 (Cat. No. 98CH36181), volume 4, pages 1941-1944. IEEE, 1998.
Fabian Theis. Towards a general independent subspace analysis. Advances in Neural Information Processing Systems, 19, 2006.
Durk P Kingma, Shakir Mohamed, Danilo Jimenez Rezende, and Max Welling. Semi-supervised learning with deep generative models. Advances in neural information processing systems, 27, 2014.
Irina Higgins, Loic Matthey, Arka Pal, Christopher Burgess, Xavier Glorot, Matthew Botvinick, Shakir Mohamed, and Alexander Lerchner. beta-vae: Learning basic visual concepts with a constrained variational framework. In International conference learning representations, 2016.
Romain Lopez, Jeffrey Regier, Michael I Jordan, and Nir Yosef. Information constraints on auto- encoding variational bayes. Advances in neural information processing systems, 31, 2018.
Zoe Piran, Niv Cohen, Yedid Hoshen, and Mor Nitzan. Disentanglement of single-cell data with biolord. Nature Biotechnology, pages 1-6, 2024.
Arthur Gretton, Olivier Bousquet, Alex Smola, and Bernhard Schölkopf. Measuring statistical dependence with hilbert-schmidt norms. In Proceedings of the 16th international conference on Algorithmic Learning Theory, pages 63-77, 2005a.
Elnaz Barshan, Ali Ghodsi, Zohreh Azimifar, and Mansoor Zolghadri Jahromi. Supervised principal component analysis: Visualization, classification and regression on subspaces and submanifolds. Pattern Recognition, 44(7):1357-1377 2011.
Wan-Duo Kurt Ma, JP Lewis, and W Bastiaan Kleijn. The hsic bottleneck: Deep learning without back-propagation. In Proceedings of the AAAI conference artificial intelligence, volume 34-04,
pages 5085-5092, 2020.
Yazhe Li, Roman Pogodin, Danica Sutherland, and Arthur Gretton. Self-supervised learning with kernel dependence maximization. Advances in Neural Information Processing Systems, 34: 15543-15556, 2021.
Xiaochun Cao, Changqing Zhang, Huazhu Fu, Si Liu, and Hua Zhang. Diversity-induced multi-view subspace clustering. In Proceedings of the IEEE conference on computer vision and pattern recognition, pages 586-594, 2015.
Abubakar Abid, Martin J Zhang, Vivek K Bagaria, and James Zou. Exploring patterns enriched in dataset with contrastive principal component analysis. Nature communications, 91):2113,, 2018.
Abubakar Abid and James Zou. Contrastive variational autoencoder enhances salient features. arXiv preprint 19t.:/002200000. 2019.
Xinming Tu, -hh-hhiistian Hütter, Zitong Jerry Wang, Takamasa Kudo, Aviv Regev, and Romain Lopez. A supervised contrastive framework for learning disentangled representations of cellular perturbation data. In David A. Knowles and Sara Mostafavi, editors, Proceedings of the 18th Machine Learning in Computational Biology meeting. volume 240 of Proceedings of Machine Learning Research, pages 90-100. PMLR, 30 Nov-01 Dec 2024. UtR.tttss::/// hhttts://prrocedinng mlr, wwwpppesss///ttuiiuuu wtt.0//t24a...a.... html

---

Wei Qiu, Ethan Weinberger, and Su-In Lee. Isolating structured salient variations in single-cell transcriptomic data with strastivevi. bioRxiv, pages 2023-10, 2023.
Ruoyu Sun and Zhi-Quan Luo. Guaranteed matrix completion via non-convex factorization. JEEE Transactions Information Theory, (2(11))65356657) 2016.
Rong Ge, Chi Jin, and Yi Zheng. No spurious local minima in nonconvex low rank problems: A unified geometric analysis. In International Conference on Machine Learning, pages 1233-1242. PMLR, 2017.
Alexander Smola and Risi Kondor. Kernels and regularization on graphs. In Learning Theory and Kernel Machines: 16th Annual Conference on Learning Theory and 7th Kernel Workshop, COLT/Kernel 2003, Washington, DC, USA, August 24-27, 2003, Proceedings, pages 144-158. Springer, 2003.
Jianming Miao and Adi Ben-Israel. On principal angles between subspaces in m. Linear algebra and its applications, 171:81-98, 1992.
Anamika Kashyap. Manjula Jain, Shailaja Shukla, and Manoj Andley. Role of nuclear morphometry in breast cancer and its correlation with cytomorphological grading of breast cancer. A study of 64 cases. Journal of cytology, 35(1):41-45, 2018.
Christoph Bock. Analysing and interpreting dna methylation data. Nature Reviews Genetics, 13(10): 705-719, 2012.
Ethan Weinberger, Chris Lin, and Su-In Lee. Isolating salient variations of interest in single-cell data with contrastivevi. Nature Methods, 20(9):1336-151 1345, 2023.
Arthur Gretton, Alexander Smola, Olivier Bousquet, Ralf Herbrich, Andrei Belitski, Mark Augath, Yusuke Murayama, Jon Pauls, Bernhard Schölkopf. and Nikos Logothetis. Kernel constrained covariance for dependence measurement. In International Workshop Artificial Intelligence and Statistics, pages 112-119. PMLR, 2005b.
Bernhard Schölkopf, Alexander Smola, and Klaus-Robert Müller. Kernel principal component analysis. In International conference artificial neural networks, pages 583-588. Springer, 1997.
12

---

A Connections with contrastive representation learning
Contrastive representation learning, as introduced Abid et al. 2018 aims to find representations that capture the difference between explicit positive and negative samples. For example, contrastive PCA (cPCA) takes the pair of a target dataset of interest and a background dataset as inputs and finds the principal components that account for large variances in the target but small variances in the background (Abid et al., 2018). cVAE [Abid and Zou 2019] and contrastive [Weinberger et al. 2023] apply the idea to VAE architectures, where the target and background data are assumed to have different generative processes and the corresponding representations lie in separate subspaces. The learning of the target and background subspaces achieved by explicitly forcing target and background samples to use different spaces. More recently. Tu et al. (2024) and QQu et al. [203] adopt an HSIC-based regularization approach to encourage the disentanglement of the two subspaces.
Despite the similarity, our proposed sisPCA model conceptually different from these contrastive models. SisPCA (and the non- oin-iinear HCV (Lopez et 2018)) focuses on the decomposition of
a single dataset, as opposed to a pair of target and background datasets. The sense of contrast in sisPCA comes from explicit supervision, of which each subspace either aligns with or is repulsed from. we concatenate the pair of target-background into one dataset and consider the case-control information supervision, the associated sisPCA subspace can also reflect the difference similar to cPCA (the objective is still not the same).
B Properties of the sisPCA-lineat optimization problem
B.1 Proof of Remark 3.1 on the equivalence of sisPCA-lineat and autoencoder
Here we show that thesisPCA SrrissCAA-linen -ineer optimization 2 is equivalent to solving a regularized linear autoencoder. First note the well-known fact of PCA that maximizing the variance of the projection ||UU| is equivalent to minimizing the reconstruction error of linear orthogonal autoencoder, since
|| XX-XXUTII 11 lr - 2x" XUUUTUU X) = -rXX))
Now consider supervised PCA with positive semidefinite kernel K2 and its Cholesky decomposition Kg L. Note L : I in vanilla unsupervised PCA. That is, PCA can be viewed aligning to supervision kernel where every sample forms a group of itself. Without loss of generality, we assume X and Ky to be centered, such that the first term in eq per subspace becomes
tr(XUU (IXUU"X"K; K",)) = tr(U" x" XyX = ||XUU||
Maximizing the above quantity is thus equivalent to minimizing the weighted reconstruction error || - -XUUT)) where L from supervision controls the aspects to focus on.
B.2
The alternating optimization algorithm for sisPCA-linear
Below we describe the Algorithm to solve 2
B.3
Proof of Remark 3.3 on the equivalence of SisPCA-linews and regression
Here we show that the sisPCA-linear problem is equivalent to regularized linear regres- sion. For simplicity, we assume both X and Y
centered and Since X"Ky XTYYTX X is rank-one, we need only to consider the largest eigenvalue and corresponding eigenvector of X Ky Therefore, maximizing the first cross-covariance term tr(UTX XT XU) = ||YT XU
Rxx
univariate.
I in is equivalent to finding unit vector E Rpx1 by
argmax] |||||=1
Y Xu 2 = argmin([|y -Xull² [ ||
(3)
-argmin(]]y ||]|==
Xull2 + ||X Xun xuu""l13).
(4)
13

---

Algorithm Solving sisPCA-linear using alternating eigendecomposition
input: Data X E R"x,, kernels for target variables (Y 1....,m, m), dimensions of subspaces [dj.] = .... m), independence regularization A.
Randomly initialize and orthonormalize U, and, optionally, determine the optimal order of subspace updates through one complete update cycle. Assuming with update order {I,..., m... repeat
K2, f- XU,U? X for i E {I...,,) for j E {1,.,, m] do A Lm K2 35+HH1 H(y,, 2 Li=l,zy H:
Update Uj + the first d, eigenvectors of X KXK Update Kz, XU,U/X*;
end for until converge
return: Linear projections {U,} and subspace representations of data {2 == XU;) for (.....m..
That is, maximizing the linear HSIC(Xu, Y) solves the regression problem of approximating the target Y with serving as the regression coefficients, while also maximizing variations of the fitted value y == Xu (or, equivalently, minimizing an additional self eeffssuperviee reconstruction loss). The particular regularization term ( corresponds
zero-mean multivariate Gaussian X x into p SSP where P is orthogonal and = (ai, J contains the singular values of in decreasing order. Denote == Pu E Rxx11 Now we have
(MVN) prior on u. To see this, we first diagonalize X
P Σ
-11Xu|1 = (PPu)
S(Pu) -(Pu) oR
(o - -oi-o?)k
-22--111A41 - of11111 ||l|E=
||X w| of
where == (o,I-S)) Subsequently, becomes
argmin(I|Y - Xull + ||xu|19),
||
and solving (5) is equivalent to finding the maximum a posteriori estimator of the regression co- efficients u under a MVN prior with zero mean and inverse covariance X X. Directions corre- sponding smaller singular values will receive stronger shrinkage effects. The formulation of (5) is also closely related to Zellner's s g g-prior in Bayesian regression where the coefficient 3 follows (V,(00 goo (X"X) for some positive and noise variance a². Finally, when fixing other subspaces }. HSTC(Xu, u, Z.) == (Xu)² Kz, (Xu), now a function of u, can also be integrated into (5) The new quadratic regularization term is Q(u) == w" (X" (X"'K, X + 8 X yu, which again corresponds to a zero-centered MVN prior.
C Optimization landscape of sisPCA-linear and Conjecture 3.1
We first provide motivating example to see how the first supervision term in affects the optimization landscape. Assume no supervision Y (i.e. KY I) and all the latent subspaces {Zm) share the same dimension d. In the absence of disentanglement regularization, every subspace will simply be the same vanilla PCA space, leading to complete overlap. The introduction of the disentanglement penalty in 22) maintains the objective's symmetry with respect Um- In other words, subspaces are interchangeable the pure unsupervised setting, and consequently the symmetry gives rise to multiple local optima. Nevertheless, the presented challenge is trivial, in that all local optima are also global, and the number of local optima, considered in terms invariant groups up to column bign-flipping depends only on the number of subspaces
Under supervision, some subspaces may be selectively favored. For example, the scaling of subspace supervision kernels {Ky... can break the symmetry of the overall objective, allowing information to preferentially preserved in one subspace while being removed from others. In such scenarios (referred to as unbalanced supervision), local optima that are closer the supervised PCA results of
14

---

subspaces with larger Ky,. also yield higher objective values. The balanced supervision condition can be evaluated by calculating the relative ratio of kernel scales, or equivalently, the ratio of objective gradient with respect to each subspace kernel.
Another important observation is that, because the disentanglement penalty is symmetric, tuning the hyperparameter À will not alter the ranking of local optima. This suggests straightforward initialization strategy to navigate towards the global solution by evaluating the relative supervision strength, HSIC(ZZ, ( Ym), for each subspace. Here, Zm represents the initial, pre-disentanglement solution for subspace m, which can be efficiently computed using supervised PCA.
We can directly visualize the optimization landscape in a simplified scenario (Fig. 77 Consider
X -1 0) which consists of two 2D vectors separated along the first axis. With supervisions
Y = -1.-1)" and Y2 = cV1 differing only by a scalar c, our goal is to determine the projections U (u, u2) and V = (2, )² that transform X into subspace representations Z1 XU = (u, -u)
11"
and 7/ = XV = (u, - Adopting a linear kernel for Ky, the sisPCA- CMPPC--linear objective ( is
now equivalent to the following simplified problem
argmax f(u, v) = u + cv Auv. u,
Since X only has one effective dimension, both U and align with that axis (i.e., u and v 1, the unregularized solution) when the disentanglement penalty À is low (A < 1), leading to a single optimal solution (u. v) = (1, 1). Under strong regularization, two local maxima (1, 0) and (0. emerge, and their relative order depends solely on the scalar and is unaffected by X. Specifically, each maximum represents a preference for capturing the significant axis in within one subspace, allowing for easy decision initialization based on comparison between and co. In the context of alternating optimization methods like Algorithm I it can be simply done by fixing the subspace V (v 1) and updating U first to minimize the disentanglement penalty (from u the unregularized solution, to u = 0). The optimization path stays near the global optima regardless of X.
Weak disentanglement
Balanced superversion
Unbalanced superversion
(c= Vambda 0.5)
(c 1, Vambda
33=3, Vambda 5)
1.00
1.00
1.00
0.75
0.75
Nu. v)
1tu. v)
Nu. v)
1.5
1.0
0.0
0.25
0.00
0.00
0.00
0.25
0.35
Figure 7: A simplified example of the optimization landscape of sisPCA-linear. Under balanced supervision (c = 1), the symmetry-induee two local maxima are both global regardless of X. When # 1, one solution is favored than the other, and the relative order is also independent of A.
This phenomenon also resembles the situation described in Ge et al. (217)) on S(UVT), where an additional regularizer || VTV is necessary maintain balance when dealing with asymmetric matrices. The first half of Conjecture 3.1 is indeed motivated by the above result. For sisPCA, a similar regularizer in the form of ||, Z,Z7" I| may be enough ensure balanced supervision. Moreover, given that the multi-solution challenge arises from subspace interchangeability. we hypothesize that more discriminative Y could alleviate the issue. Since there are only pairwise interactions between subspaces in (2) Ge et al.] [2017)ss framework can be adapted to study the behavior of more than two subspaces. Intuitively, the proof of Conjecture 3.1 would be to show that the local gradient direction of (2) aligns with the global descent direction in region near global optima, and proper regularization or initialization can make sure the optimizer stays in that region. Following the procedures described in section 5 of Ge et al. [017]] for asymmetric matrices, we may concatenate all subspaces into W := [Z,,.., 7m and reformulate the objective as a quadratic function of N := WWT plus some regularization.
15

---

D Design and applications of sisPCA-general with Gaussian kernel
In this section, we discuss the sisPCA-gemeral optimization problem in detail. Recall that it is motivated by the limitation of sisPCA-linear that zero HSIC disentanglement regularization does ensure statistical independence. In sisPCA- -generall we want solve ( by employing HSIC with with universal kernels, e.g. the Gaussian kernel.
Definition D.1. Given compact metric space X' and a continuous kernel K on X, K called universal if for every continuous function 9 E C(A) and all c >> 0, there exists an f E H. such that Ilf - <.
nice result on universal kernel from the Constrained Covariance (COCO) framework [Gretton al., 2005b) is that can approximate any continuous function on compact metric space. In other words, the HSIC with a universal kernel will be zero only if, among all possible continuous transformation, the covariance of two random variables are always zero. In sisPCA-general. the latent representation is still encoded by a linear projection Zj = XU,. This indicates that Uj can be interpreted in the same way as PCA and SisPCA- -ineer to extract feature importance.
We solve the SisPCA- sisPCA-general problem with Gaussian kernel.
k(zz, x²) exp(-wOl|z - x211²), -
using gradient descent (Algorithm 2). Since the problem is constrained, the algorithm unfortunately is not guaranteed converge. The HSIC with Gaussian kernel also affects the optimization landscape. Empirically, we observed that tive tuning of training parameters such as the learning rate and stopping criteria. The Gaussian kernel scale is also important. Ma et al. [2020) reported in their experiments that the performance of HSIC-based optimization moderately sensitive to w and thereby proposed a multiple-scale solution that essentially aggregates results from different w. In our experiments, we followed Gretton et al. 2005a by setting the kernel scale to the median distance between data points in the input space.
tends to have many local optima, requiring exhaus-
Another notable disadvantage of SisPCA- geceerlll is its cubic-scaling computational complexity. Here we focus on scenarios where the number of features p is fixed and the sample size n is large (n > p). which is common when data preprocessing and filtering is in place. Specifically, the computational complexity for calculating HSIC(X, Y) with linear kernel for E R and given target kernel Ky for Y is O(n²p). If both Kx and Ky are linear, the complexity can be further reduced to (mp). Additionally, the eigendecomposition of the (cross) covariance incurs a cost of Ogp"). This means that each update in Algorithm [ has overall complexity of O(n²) wrt the sample size n, aligning with that of supervised PCA. comparison, most dimensionality reduction methods requiring pair-wise distance computation have complexity of at least O(n²). On the other hand, for sisPCA- gennerale computing a Gaussian-kerns HSIC Y) naively is of (n²)) complexity. From engincering perspective, further optimizations for both the MSPPA-linneal (Algorithm M) and sisPCA- -eeeCr-leeerall (Algorithm 2 algorithms are conceivable, for instance, by exploiting data sparsity and by mini-batch or incremental update.
Algorithm 2 Solving
using gradient descent
input: Data X E _ kernels for m target variables [Ky,,j I,..., m). dimensions of subspaces {d,. = ..., m), universal kernel (z((;-) for latent representations, independence regularization X.
Randomly initialize and orthonormalize U, and, optionally, determine the optimal order of subspace updates through one complete update cycle. Assuming with update order {1, m). repeat
K2, (k, r) t- Pz (X,,, X,U for i (1,..,mm:
L - b
(r(ntt((ZZ, 2 H1:; En, tr (Kz, HK2,HKHH H);
L ++
Jointly update {U,} by minimizing L using gradient descent; Project U, the orthonormal space;
until converge
return: Linear projections {Un} and subspace representations of data (2, XU3) for
(1,.., m).
16

---

Recovering supervised and unsupervised subspaces in simulated data We demonstrate the performance and limitation of sisPCA-general using simulated data. The task and simulation details are described in Section 4.1 As illustrated in Fig. sisPCA-general partially recovers the two supervised subspaces S1 and S2. Challenges arise, however, fully eliminating continuous S2 signatures from $1. This difficulty may stem from the smaller scale of HSIC values with Gaussian kernels as compared to those with linear kernels, and a more complicated optimization landscape due to the non-linearity introduced in the regularization. Additionally, the unsupervised subspace S3 appears to collapse to one dimension, which partially illustrate the identifiability issue raised in Section 3.3 The issue may be more pronounced in SisPCA- -.PPA-geneeal -eneralerllll because of the diminished HSIC scale (i.e. weaker connection between subspaces and their target variables).
81: Categorical
81: Categorical
S1: Categorical
X+
PC
P1PC2
Group
51
-10
-10
-22
-10
2
81:
81: PC
S1: PC
52: Continuous
S2: Continuous
52: Continuous
Group
-15
-12
S3: Unsupervised
Unsupervised
Unsupervised
X+Y
PC2
S3:PC2
Group
$3:P0
S3:
S3:
Figure 8: Performance of SisPCA-generes (A = 30) on the simulated dataset. Related to Fig.
E
Extending sisPCA to nonlinear settings using kernel transformation
Following kernel PCA's methodology (Schölkopf et al. 1997]. we extend si sisPCA for nonlinear scenarios where features X interact through a feature map / : R" - F. The space F may have infinite dimension, and we can only explicitly write out the projection U in the finite case F R for example, ((x) = x². In that context, the latent representation Z := & 4(X) can be naively computed, and the sisPCA objective (T) will stay the same. By simply replacing X with ((X), methods such as gradient descent and Algorithm still hold. However, the beauty of the kernel trick that bypass the computation of (xX) obtain Z. Instead of learning the projection matrix U, which is no longer possible when F is infinite-dimensionn we can learn a set of coefficients a from
Z Z= OX}(UU Σ d
with (X,) representing the data points and k( y) == (xx,, (yy)) the corresponding kernel of &. Note that k(.,) .)) is specified wrt to the data X, which is different from the kernel K2 for latent space and the kernel K, for supervision. We leave the practical implementation of this extension future exploration.
17

---

Table 4: General comparison of baseline models.
Linear
Non-linear
PCA
SPCA
sisPCA
VAE
SupVAE
hs
Supervision
HSIC
HSIC
Prediction
Prediction
Disentanglement
HSIC
HSIC
Interpretability
U as feature importance (1,2.3)
Black-box
Hyperparameters
(1)
(1,2)
[1.2]
[1,2.3]
[1,2,3,4]
(I) Subspace dimension
[1] Subspace dimension
Subspace kernel
[2] Autoencoder design
HSIC penalty strength X
[3] Predictor design
[4] HSIC penalty strength
Optimization
Closed form
Simple
Subject general training limita- tions of deep generative models*
*VAEs could not run on 6-dimensional simulated data in Figure ddue to NaNs generated during variational training. This is common issue in training SCVI models likely resulting from parameter explosion.
F Baseline section and experiment details
F.1
Overall evaluation criteria
We evaluate model performance based on two key aspects:
Representation Quality:
This criterion assesses the ability to learn low-dimensional representations that reflect specific and unmixed data properties. A high-quality subspace representation should contain minimal confounding information from other subspaces. We employ both quantitative metries (e.g., information density and subspace separateness) and qualitative 2D visualizations of learned representations.
Interpretability: This criterion evaluates the ease of understanding feature contributions to different data properties. For linear methods (PCA, sPCA, and SisPCA), consider the learned projection U as feature importance and extract top features from the loading (see examples in Sections 4.2 and 4.4). Non-linear models generally lack straightforward feature-to-subspar them inherently less interpretable. While approaches like gradient-based saliency maps are becoming standard for interpreting black-box models such VAEs, they are less straightforward and face challenges in aggregating sample-level gradients into global feature importance scores. Therefore, we only compare sisPCA to other linear models for interpretability analysis (e.g. ranking and selecting top genes according to the loading).
mapping, making
F.2 Baseline model design and general comparison
Table 4 provides a comprehensive comparison of linear and non-linear models used in this study, including the linear PCA, sPCA and sisPCA (proposed work) described in the Methods section, as well as sevral non-linear VAE counterparts inspired by the work of HCV [Lopez et al. 20188:
VA Vanilla VAE with Gaussian likelihood.
VAE: Gaussian VAE with additional predictors for target variables. hsVAE: sup VAE with additional HSIC disentanglement penalty (Gaussian kernel). is essentially a general- purpose HCV [Lopez et al. 2018 with Gaussian likelihood.
In the SisPCA package PCA and sPCA are implemented special cases of sisPCA and solved analytically using cigendecomposition. VAE models and variational inference are implemented using the SCVI f The VAE model is special case of SCVI's VAE module (latent distribuutin "normal". dispersion = "gene"" but with a Gaussian generative model. The VAP is a multi-subspace extension of VAE with additional predictors to predict target supervisions from the corresponding
"https//githu..com//iii hatp://hh..oom/aaaaaauu/ss/ssss
ssevitoolls v12.0 (httss//seeiitooll.orgg)
18

---

latent representations. The last layer of the predictor is either linear transformation for continuous supervision or linear followed by softmax for categorical supervision. The hsVAE model is an extension of supVAE, where we add additional subspace disentanglement penalty (i.e., HSIC with Gaussian kernel) to the objective, as described in HCV [Lopez et al. 2018)-
For single-cell RNA-seq data, we further implement domain-specific hs VAE model, namely hs VAE- sc. Specifically, Vs uses negative binomial as the data generative distribution with an extra autoencoder to learn the library size (the sum of total counts per cell) and is correspondingly trained on louell--eveel data. It is very close to a implemeentation of the se ScViGenQCMode developed in the HCV paper under the latest SCVI framework.
F,3 Hyperparameter selection
As outlined in T SisPCA has three main hyperparameters:
1. Subspace latent dimension d.
2. Subspace supervision kernel Ky.
3. Disentanglement penalty strength
The tuning of (3) will be discussed in Appendix G Here we illustrate the selection of (1) and (2), which are largely pre-determined based on the target variable supplied for supervision.
In our analysis, we consistently use the delta kernel for categorical targets and the linear kernel for continuous targets (Section 3.2) Given supervision kernel Ky of rank k. the effective dimension of the corresponding sisPCA subspace, defined as the number of dimensions with non-zero variance, is determined by the eigendecomposstion step in Algorithm T Specifically, we have the following upper bounds:
For sPCA (A = 0): The effective dimension equals rank Xy S min(k;p) s k, where p is the number of features in X.
For SisPCA: The effective dimension equals rank(X Ry S min(k: + d. R + where Ky <- Ky A is the supervision kernel after the rank-d disentanglement update and d is the number of latent dimensions of Z.
In practice, we observe that the effective dimensions of both sPCA and sisPCA closely approximate the target kernel rank k. For instance:
For 1-D continuous variable with linear kernel (e.g., age in Section 4.3) the effective dimension is always 1.
For a categorical variable of l groups with delta kernel (e.g., time point in Section 4.4). the effective dimension is approximately k.
Consequently, we find that the learned sisPCA subspaces remain consistent regardless of the specified dimension. Any additional axes beyond the effective dimension tend to collapse to zero and are primarily affected by numeric errors in the SVD solver. Since sisPCA orders subspace dimensions by the variance explained (i.e., eigenvalues). it straightforward to remove these extra dimensions by examining the variance explained curve.
In contrast, VAE models are more sensitive to dimension changes and generally have more hyperpa- rameters to tune (Table Designing optimal architectures for autoencoders and target predictors is a topic in itself and beyond the scope of our current work. For benchmarking purposes, we largely follow the default VAE architecture from SCVI:
Autoencoders for latent mean and variance: One hidden layer with 128 hidden units, ReLU activation, batch normalization, and dropout.
Predictor design (adapted from the ScVIGenQCModel in the HCV paper): One hidden layer with 25 neurons, ReLU activation, and dropout. Extra softmax output head for classifier. Training objective: Equal weighting prediction error and reconstruction loss, minimized using variational inference.
we set the subspace dimension to 10 for all VAE models. We use UMAP to project the learned subspace onto 2D and subsequently calculate the Silhouette score in Table 3 on this 2D subspace, ensuring dimensional consistency with the PCA-based spaces (which are also projected onto 2D using UMAP)
In Section 4.4
19

---

F.4 Additional experiment details
Recovering Supervised and Unsupervised Subspaces in Simulated Data The code used to simulate the donut data in Section 4.1 is available in the sisPCA GitHub repository. We have verified that altering the noise distribution (e.g., from Gaussian to uniform noise, or from uniform projection matrix to Gaussian matrix) yields similar results.
Learning Diagnostic Subspaces from Breast Cancer Image Features We obtained the dataset from Kaggle and scaled each real-valued feature to have mean and unit variance. The input for PCA includes all 30 quantitative features. The 'radius' subspace is supervised with 'radius mman" and "radius_ sd' (effective dim 2), while the 'symmetry' space is supervised with 'symmetry_mme meen"" and 'symmetry_sd _sd (effective dim = 2). For sPCA and sisPCA, we use the remaining 26 features as input and project the data onto each subspace. Using the elbow method, dimension to 6, and correspondingly, the dimension for sPCA and sisPCA subspaces to (although their effective dimension approximately 2). The Silhouette scores are computed on the first three principal components for all subspaces.
set the PCA subspace
Separating aging-dependent DNA methylation changes from tumorigenic signatures We downloaded the TCGA methylation data from GDC (controlled access) and retained only the first 5,000 non-constant and non-NA CpG sites for analysis. The data was zero-centered. All subspaces specified with 10 latent dimensions. However, noted in the previous section, the sPCA and sisPCA aging subspaces have an effective dimension of only one.
Disentangling infection-induced changes in the mouse single-cell atlas of the Plasmodium liver stage We used the single-cell data preprocessed by Piran et al. (2024] 7 comprises 19,053 cells and 8,203 mouse genes (all malaria genes have been filtered out), and we keep the top 2,000 variable genes using scanpy's pp..iphlyy........ variable geneseeees'ee function. For all models except hs VAE-sc, we use the normalized and log I p-transformee expression data as inputs. For VAE-Ss, we supply the raw data before normalization. It's important to note that the data underwent background correction step before normalization, where the mean expression in empty wells was subtracted from the observed data. This results in the raw counts being technically non- non-integerss Consequently, in Vs,-s,, a continuous extension of the negative binomial likelihood is used to compute the reconstruction loss. All subspaces are given 10 latent dimensions. Although based on the kernel rank, we know that the sPCA and sisPCA infection subspaces have an effective dimension of 2, while the time subspaces have an effective dimension of 6.
The processed data
vigiabl-variable
G Tuning the disentanglement strength in sisPCA
We adapt the approach from Abid et al. 2018 to develop computational pipeline for automated selection of the sisPCA hyperparameter X. The process involves:
1. Fitting SisPCA models with multiple values.
2. Computing pairwise similarity (affinity) between learned subspaces.
3. Identifying representative subspace clusters using spectral clustering.
We compute subspace affinity based on the principal angle [Miao and Israel 1992]:
daubspace Z1 Z2) d((e).... 1x d(O)
where Z Ro is the orthonormal basis of a dim-k subspace with n observations, and the principal angles O are computed from the SVD of Z Ze = U cos(e)
Through empirical testing, found the Fubini- Study Grassmann distance provides optimal perfor- mance. Thus, we define subspace affinity as:
Affinity subspace (Z1, 2z)
II cos cos(O.)
II t=1
(z Zz).
'See the notebook cemporal-infection
20

---

For stability, we truncate subspace latent dimension to k' := min mii(rank(Kyy).e where Ky is the corresponding supervision kernel (see discussion on subspace effective dimension Appendix [F3). For sisPCA models with multiple subspaces, we compute model-wise affinity as the average of subspace affinities
Affinity model ((AI, AmJ (Bi, Bm))
Affinity subepace (Ai, Bt).
Learning Diagnostic Subspaces from Breast Cancer Image Features We evaluated 20 À values ranging from 0 to 100 and analyzed the pairwise similarity between learned sisPCA subspaces (Fig. 99. The sPCA solution (A 0) shows clear separation from other SisPCA models. As À increases, the symmetry subspace becomes progressively less predictive diagnostic status (Fig. [))), and subspaces stabilize and converge robust solution after À 1. This convergence pattern is also
reflected the elbow of the reconstruction loss curve (Fig. b
0.147
0.95
0.215
0.316
0.464
0.881
1.0
0.85
1.488
2.154
dagnzsa
4.622
6.813
0.75
10.0
34.678
21.544
0,70
31,623
46.416
0.65
0.60
0
0213
03T
0.81
1.
1,468
2.154
3162
6442
6,813
10.0
14.678
21.544
28.32
100.0
105
165
492 0.
203
(a) Pairwise similarity of models with different
(b) Evolution of symmetry subspace
recon_ Jocon_loss
reg_loss
train loss
-88
18
1.0
loss
-20-
0.5
-21
0.0
-21
0
25
50
75
100
25
50
76
100
0
25
50
75
100
lambda__cnne
(c) Training loss as functions of À (total loss reconstruction loss disentanglement disettangdemennreequuuzzaai
Figure 9: Effect of À on the learned subspace structure the breast cancer dataset. Related to Fig.
Disentangling infection-induced changes in the mouse single-cell atlas of the Plasmodium
liver stage We examined 10 values ranging from 0 to 100 assess their effect on subspace characteristics (Fig. 10). Our analysis again reveals the clear separation of sPCA solution (A = 0) from other sisPCA models. As increases, infected cells form a tighter cluster in the infection subspace (Fig. 100a,, while temporal dynamics become less pronounced, with increased mixing of cells from different time points (Fig. [10p). Moreover, sisPCA results are quite robust across À in preserving subspace structure (Fig. 10P). Since the disentanglement effect primarily manifests in PC2 of the infection subspace, all models, even sPCA, extract similar sets of top contributor genes in PC1 of the infection subspace, thus yielding GO enrichment results comparable to those shown in Fig. 6
21

---

The infection subspace
The infection subspace
0227
0.22
01
0227
Fmm
Tites
724
17.29
Fileded
PALSE
0 10
40
10
0
10
10
6:
(a) The infection subspace, colored by infection status
(b) The infection subspace, colored by time
10
100
0.0
0.0
100
71
71
74
70
70
0.1
0.9
0.237
100
80
80
B0
80
80
77
0.562
0.562
100
80
81
79
76
1.34
0.7
3.162
100
85
80
80
0.6
7.499
7,499
200
80
84
17.783
42.17
0.4
42.17
100
80
100.0
100.0
100
00
0,237
0,562
3182
7,499
17.78)
42.17
100.0
0.0
0.237
0.562
1.334
3302
7.499
17.783
42.17
100.0
(c) Model similarity across A
(d) Percentage of shared top genes contributing to PC1
of the infection subspace learned across À
Figure 10: Effect of À on the learned subspace structure in the single-cell malaria infection data.
Related to Fig. 5 and Fig. 6
H Supplementary visualizations of baseline performance in applications.
For completeness, here provide subspace visualization for baseline models not included in the main figures. Specifically, Fig. [1 shows the full sPCA performance on the simulated data, related to Fig. 12 visualizee the single-cell subspaces learned by sPCA, VAE, sup VAE, and hs VAE-sc,
Fig.
related to Fig
22

---

51: Categorical
S1: Categorical
51: Categorical
PC2 (sPCA)
PC2 (cPCA)
PC2 (sPCA)
Group
0.0
-0.5
10
20
-10
PC1 (sP CA)
PC1 (sPCA)
PC1 (sPCA)
S2: Continuous
S2: Continuous
S2: Continuous
10-
+Y
Phi
PC2 (sPCA)
PC2 (sPCA)
(SPCA)
Group
0.0
-5
-0.5
-10
-10
-10
0
PC1 (sP CA)
PC1 (sP
PC1 (sPCA)
S3: Unsupervised
S3: Unsupervised
S3: Unsupervised
X
Phi
PC2 (sPCA)
PC2 (sPCA)
G Group
0.0
-5
-0.5
-10
-10
-10
-30
10
PC1 (sP CA)
PC1 (sP CA)
PC1 (sPCA)
Figure 11: sPCA results on the simulated dataset. Related to Fig. 3b. Rows from top to bottom: sPCA results with SI categorical supervision; with $2 continuous supervision; without supervision (vanilla PCA). In this case, all sPCA subspaces are dominated by S2 since it carries the strongest variability.
23

---

VAE iilhouette -0.036)
Infection (S Sillouuttee 002077)
Time (S (Silhouette = -0.029)
Harboring parasite FALSE TRUE
VAE (S ISilhouette 0.296)
Infection (Sithouelte -0.075)
Time ISithouelle 0.3488)
Post-infection time point
Control hpl
12 hpi
24 hpi
30 hpl
36 hpi
(a) VAE
(b) sPCA
Infection (Silhouette 0.041)
Time (S (Sihooutee -0.087)
Infection (Silhouette 0.233)
Time (S (Sihouutte -0.100)
Infaction (Silhouette 0.052)
Time (Silhouette 0.479)
Infection (S Silhouetee 0.015)
Time (S ISilhouette 0.634)
(c) VAE
(d) hs (A 10)
Figure 12: UMAP visualizations of the ScRNA-seq data of mouse liver upon Plasmodium infection. Subspace representations are learned using unsupervised VAE (a) and supervised sPCA (b), sup (c) and hs VAE-sc (d). Note that the infection subspaces of VAE and sup VAE fail to distinguish infected versus uninfected cells. Moreover, all infection subspaces presented here still exhibit significant temporal patterns (lower left plot in each panel) where cells collected at different time points are not fully mixed.
24

---

NeurIPS Paper Checklist
1. Claims
Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Question:
Answer: [Yes]
Justification: In the paper, we proposed nulsi-subspacce extension of PCA to disentangle interpretable subspaces of variations using self-supervision. We believe our work is highly novel and significant to practitioners using PCA as well as general audience computational biology.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
It is fine to include aspirational goals as motivation as long as is clear that these goals attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations the work performed by the authors?
Answer: [Yes]
Justification: See Section 5
Guidelines:
The answer NA means that the paper has limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect the scope of the claims made, e.g., if the approach was only tested on a few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory Assumptions and Proofs
Question: For each theoretical result, does the paper provide the full of assumptions and complete (and correct) proof?
25

---

Answer: [Yes]
Justification:
The proof of Remark 3.2 is stated following the remark in the main text. The
proof's of Remark 3.1 and 3.3 can be found in Appendix B Detailed discussion 3.1 can be found in the Appendix
Guidelines:
The answer NA means that the paper does not include theoretical results, All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but they appear in the supplemental material, the authors are encouraged to provide short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix supplemental material.
Theorems and Lemmas that the proof relies upon should properly referenced.
4. Experimental Result Reproducibility
Does the paper fully disclose all the information needed to reproduce the main perimental results the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided not)?
Question:
Answer: [Yes]
Justification: Details on data preprocessing and model configurations are provided in the corresponding sections and in Appendix A Python implementation of sisPCA is available on GitHub wwt.giiub.oo/jaauuuuuuuu UU/ssppaa. The repository also contains notebooks reproduce most results in the paper main text and in appendices.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution be necessary to either make it possible for others to replicate the model with the same dataset, provide access the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how replicate the results, access a hosted model (e.g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
specific model and empirical evaluation, may
While NeurlPS does not require releasing code, the conference does require submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how reproduce that algorithm.
(b) If the contribution is primarily
model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (c.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the partieular way they provide for reproducibilitty In the case of closed-source models, it may be that access to the model is limited in
26

---

some way (e.g., to registered users), but it should be possible for other researchers have some path reproducing or verifying the results.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [Yes]
Justification: We do not generate new data in this study and use public datasets for all analyses. See detailed description on data access and license in the cor- responding sections. A Python implementation of sisPCA available on GitHub http://thww..oooaaaauauauuuuaA
The repository also contains notebooks to re- produce most results in the paper main text and in appendices.
Guidelines:
The answer NA means that paper does not include experiments requiring code. Please see the NeurIPS code and data submission guidelines (https://nips... cc/
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https
/ itipss.ccppublic///
for more details.
The authors should provide instructions data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: Details on data preprocessing, model configuration and hyperparameter selec- tion are described in the corresponding sections and in Appendix F
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, or as supplemental material.
7. Experiment Statistical Significance
Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer: [No]
Justification: Statistical significance tests are not included in this paper due to the nature our experiments, which primarily focus on qualitative performance assessments and biological interpretations lacking quantitative null hypotheses. The quantitative performance
27

---

evaluations presented in the tables are self-evident and serve to demonstrate design principles. Moreover, our baseline models (PCA and sPCA) are deterministic, eliminating the need for error bars or confidence intervals that typically account for randomness.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, statistical significance tests, least for the experiments that support the main claims of the paper.
The factors of variability that the error bars capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
It should be clear whether the error bar is the standard deviation the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not to show in tables figures symmetric error bars that would yield results that are out range (e.g. negative error rates).
If error bars are reported in tables or plots. The authors should explain the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information on com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: We discuss the computational complexity of sisPCA in Appendix D The proposed method is highly efficient because of the linearity. We ran all provided notebooks SatayyuSSuPKSSsspaaeeee//miiiii
using personal MI Macbook Air with 16GB RAM and completed most analysis steps in minutes including model training.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary failed experiments that didn't make into the paper).
9. Code Of Ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics
Answer: [Yes]
Justification:
We confirm the work conforms with the code of ethics.
Guidelines:
The answer NA means that the authors have reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special cireumstances that require deviation from the Code of Ethics.
28

---

The authors should make sure to preserve anonymity (e.g., if there is special consid- eration due to laws or regulations their jurisdiction).
10. Broader Impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [NA]
Justification: This paper presents work whose goal is advance the field of Machine Learning. There are many potential societal consequences of our work, none which we feel must be specifically highlighted here.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g.. disinformation, generating fake profiles, surveillance), fairesss considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: This paper poses such risks.
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
29

---

Answer: [Yes]
Justification: The Kaggle Breast Cancer Wisconsin Dataset is publicly available under CC BY-NC-SA 4.0. The TCGA methylation dataset has controlled access due to privacy regulation. The single-cell malaria liver infection data is publicly available under CC-BY 4.0.
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset is used and, if possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g.. website), the copyright and terms of service of that source should be provided.
If assets released, the license, copyright information, and terms of in the package should be provided. For popular datasets, paperswi thcode. com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
this information not available online, the authors are encouraged to reach out to the asset's s creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification: The paper does not release new assets.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent obtained from people whose asset is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then much detail as possible should be included the main paper.
According to the NeurIPS Code of Ethics, workers involved in data collection, curation, other labor should be paid at least the minimum wage in the country of the data collector.
Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
30

---

Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsoureing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with
human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent)
may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such the institution conducting the review.
31


Reviews:
[{"review_id": "UPUTrIr186", "reviewer": "Reviewer_PWKM", "summary": "This paper presents a method for distentangling multiple independent linear latent subspaces that align with a set of response variables. This uses the Hilbert-Schmidt Independence Criterion (HSIC) to measure the dependence each subspace and the targe variable providing supervision for that subspace, a concept which was previously applied to infer a single subspace that was connected to a set of response variables. HSIC is additionally used to enforce independence between the subspaces, to encourage the method to find independent subspaces that correlate with each of the provided dependent variables. This expands on Supervised PCA by enabling inference of an indendent subspace for each response variable, rather than a single subspace which aligns with the vector of response variables, reflecting the fact that different responses may be interacting with different subspaces in the data, e.g., methylation changes reflecting natural human aging and changes related to a disease like cancer.\n\nThe method is evaluated on synthetic data and shown to do a good job at extracting the different subspaces. In particular, Supervised PCA is dominated by one of the two supervised subspaces. \n\nThe method is then evaluated on real methylation data and single-cell gene expression data. In the cancer data, the new method is shown to provide subspaces with greater separatio", "questions": "No questions", "limitations": "The authors describe the well-known limitations of unsupervised PCA, i.e., identifiability, and other limitations implied by their choice of kernel.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The paper tackles an important problem in the analysis of high-dimensional omics data, where there is a need to identify interpretable subspaces that can be mapped onto higher-level biological processes but where different phenotypes (e.g., aging vs cancer) may manifest differently, requiring independent subspaces rather than supervision of a single subspace as per sPCA. \n\nThe paper is well-written and combines simulated data experiments, which help understand the method is behaving as expected, with two applications on real data.", "weaknesses": "One could argue that the work is somewhat incremental. The use of the HSIC is not new and nor is the idea of providing supervision to a latent variable model to ensure the latent subspaces reflected some additional variables. However, I feel the contribution is novel enough to warrant presentation at NeurIPS", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "## **Response to Reviewer PWKM**\nThank you for your summary. We appreciate your evaluation that our technical contribution is novel enough to warrant presentation at NeurIPS.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "8C8UY1OoMU", "reviewer": "Reviewer_UeQu", "summary": "The work proposes a new method that determines the factors of variation bonded to different labels in a supervised way, obtaining a method akin to supervised PCA but in several independent subspaces thorough the Hilbert-Schmidt independent criterion (that it’s employed to maximize both the independence between different subspaces and the dependence between the data subspaces and the labels). The paper focus on a linear kernel application of HSIC, but other kernels are also considered. Then it’s applied to some simulated data and to DNA/RNA sequencing data.", "questions": "Could you please provide a comparison of your method with other SOTA methods? \n\nWould you prove that the method is useful in other kinds of data types? \n\nWhich are the hyperparameters of the method? How do you set them?", "limitations": "I’m not satisfied with the way the limitations are addressed. One should explain them clearly and provide examples in which the method is not working.", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "The theoretical connections are interesting. \n\nThe paper is technically sound, with a well-derived mathematical background. \n\nThe results of genetic sequence analysis are meaningful and interesting.", "weaknesses": "The novelty is limited, it appears to be an extension of previously described methods. \n\nThe comparison is provided only with PCA or sPCA, which are not by any means SOTA methods.  \n\nThe explanation of the simulated data is not clear, what should it be the expected result? Are you learning the dimensions of the subspaces?  \n\nA better description of the experiments is needed. How do you choose the hyperparameters (for instance, the dimension of the subspaces is 10 in both sections 4.2 and 4.3, why?). \n\nIn section 4.2, the use of HSIC as a quality metric is not fair, since it’s already included in your loss.  \n\nThe significance of the paper would be increased if the authors prove that the method is useful for the analysis of other types of real-world datasets. \n\nThe discussion is too short and not too clear. I acknowledge the limitations of space, but the discussion is one of the most important sections and should be properly done (even by moving some mathematical details to SI).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "## **Responses to Reviewer UeQu**\nThank you for your constructive feedback on the strengths and weaknesses of our work. We address your main concerns below:\n\n## **Weaknesses**\n### **W1: Limited novelty**\nWhile our method extends PCA, it offers significant contributions by uniquely combining supervision and disentanglement in a linear, interpretable framework valuable for feature selection. This addresses an important gap in current dimensionality reduction techniques, especially in an era where linearity and interpretability are often overlooked (**See M1**).\n### **W2: Comparison with SOTA methods**\nWe've added comparisons with non-linear VAE-based models (**See M2**), focusing on:\n\n1. **Representation quality:** Low-dimensional representations capturing specific data aspects.\n2. **Interpretability:** Understanding feature contributions to different data aspects.\n\nWhile the interpretability criterion clearly flavors linear models, our results show that linear models like sPCA are indeed strong baselines too in representation quality (**See M2.3, Fig S1, Table S2**). We also confirmed similar insights on supervision and disentanglement regularization effects across both linear and non-linear models.\n\n### **W3: Inadequate discussion of limitations**\nWe apologize for the poor presentation of limitations, which are somewhat scattered throughout the paper. We will concentrate them and improve the discussion section in the final version. Key limitations and their potential impacts include:\n\n1. **The limitation of linear-kernel HSIC on theoretical independence guarantee** (Section 3.1, line 159 - 165). Though in practice this is less concerning, since minimizing HSIC-linear will also reduce HSIC-Gaussian (Section 4.2.3, Table 1).\n2. **Subspace interpretability issue** (Section 3.3, line 189 – 193, and Appendix D on the resulting optimization landscape). This would make it especially challenging when the model needs to learn both supervised and unsupervised subspaces. Figure 7 is an example where unsupervised subspace failed to capture the true dimension.\n3. **Challenges in balancing the disentanglement regularization strength \\lambda** (Section 5). This is a general challenge to all supervised disentangled models. \n4. **Potential performance loss due to the lack of non-linear feature interactions** (Section 3.1, line 166-169; Section 5). This could lead to underperformance on datasets with complicated patterns.\n\n## **General questions**\n### **Q1: How to set model hyperparameters, in particular the dimension of each subspace.**\nWe again apologize for the lack of clarity. The main hyperparameters are subspace dimensions and disentanglement penalty scale λ (**See Table S1, M2.1**). The latter was discussed as a limitation in the Discussion and above in **W3**. For sPCA and sisPCA, subspace dimension generally doesn't affect performance beyond numerical issues from SVD solvers. \n- sPCA's subspace dimensions are determined by supervision kernel rank (Appendix B.3, line 406-408). That is, the effective dimension is always 1 for a 1-D continuous variable with linear kernel (age in Section 4.2) or K for a categorical variable of K groups with delta kernel (cancer type in Section 4.2) regardless of the specified dimension. The extra axes beyond the effective dimension (sPCA-age PC2-10) will collapse to zero since the eigenvalue is zero. \n- The effective dim of sisPCA subspaces is determined by both the supervision kernel rank as well as other subspaces (Remark 3.2, line 144-147). In practice, sisPCA shows similar behavior as sPCA, with extra dimensions typically collapsing to zero. For example, the sisPCA-infection subspace in Section 4.3 has approximately two effective dimensions (Figure 4, 9 and Figure S1).\n- In contrast, VAE models are more sensitive to dimension changes and have more hyperparameters to tune (**M2.1**), which is generally beyond the scope of our work. Due to time constraints, we mostly follow the SCVI’s default in designing VAE models for benchmark.\n\n### **Q2: Better description of experiments**\n#### **Simulated data Figure 3**\nAs mentioned in **Q1**, the dimensions of supervised subspaces are mostly predetermined by their target kernel ranks. Here dim(S1) and dim(S2) are 2 regardless of hyperparameter choice. The dimension of the unsupervised space S3 is unknown, and we aim to recover the ground truth donut structure (dim = 2) as did in Fig 3c. Figure 7 represents a failed example where only one dimension is recovered.\n#### **Subspace dimensions in Section 4.2 and 4.3**\nAddressed in **Q1**. Dimensions in PCA-based models are generally ranked by importance and can be selected using methods like \"the elbow curve\" (Discussion, line 311-313). In **M2** subspace dimensions of VAEs models are set to 10 to ensure fair comparison.\n#### **Applications to other types of real-world datasets**\nWe note that similar to PCA and its extensions [1], our model is also general-purpose and can assist with data exploration tasks (*as pointed out by Reviewer debZ*). The only domain-specific model is the hsVAE-sc tailored for single-cell data with a count-based likelihood. \n\n[1] Abid et al., \"Contrastive principal component analysis\"\n\n### **Q3: The use of HSIC as a quality metric in Section 4.2**\nThe HSIC metrics in Table 1 are to support the claim that minimizing the HSIC-linear can also reduce the computationally more expensive HSIC-Gaussian (not explicitly optimized for in the loss) and thus encourage independence in the representations (line 259-261).", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "cMJrJ7Wv3L", "reviewer": "Reviewer_QSK9", "summary": "This paper proposed a linear dimensionality reduction and subspace extraction method based on Hilbert-Schmidt Independence Criterion (HSIC) and Supervised PCA. In specific, several interpretable subspaces, which are independent from each other, are disentangled from the data observations and the leaned subspaces are enforced to be maximally dependent on the supervision variables. The proposed approach is shown to be effective on both synthetic data and two real datasets.", "questions": "1.\tIn Sec. 3.2, the linear kernel is used to target variables. Does it make a difference if a nonlinear kernel function is used?\n2.\tWhat are A and B in Fig. 2? Are they the target variables? In the same figure, why does it maximize dependence (HSIC) with the components of a subspace? I thought the only two objectives are 1) maximizing dependence between the subspace and its corresponding target variable; 2) minimizing the dependence between the subspaces. How is this maximization reflected in Eq. (2)?\n3.\tWhat does balanced/unbalanced supervision mean mathematically in Conjecture 3.1?", "limitations": "The authors are encouraged to discuss limitation of the work in the main paper.", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "1.\tA new method for subspace learning is proposed using supervised PCA and dependence maximization and minimization with HSIC for extracting independent latent low-dimensional subspaces. The framework is corroborated effective in extracting linearly mixed simulation data, and the learned representations are more interpretable than that of supervised PCA.\n2.\tThe problem formulation and optimization are clearly derived and delivered. The overall manuscript is written in an organized way and the notations are well-defined and relatively easy to follow.", "weaknesses": "1.\tIt is not clear if the proposed method is guaranteed to extract the latent subspaces and/or under what conditions/assumptions the method would work or fail. For example, in the synthetic experiment, the mixing matrix is uniformly drawn from [0, 1]. Does it still work if the matrix is changed to Gaussian? Overall, the theoretical analysis is somewhat lacking, e.g., identifiability analysis for linearly mixed subspaces.\n2.\tBaselines seem missing in the real data applications. Only PCA and supervised PCA are considered here. The authors are encouraged to include more state-of-the-art methods to showcase the superiority of the proposed method, e.g., how does PCA work under such a setting? Are other HSIC-based nonlinear methods, e.g., [1], capable of extracting the independent subspace with interpretability?\n\n[1] Lopez, Romain, Jeffrey Regier, Michael I. Jordan, and Nir Yosef. \"Information constraints on auto-encoding variational bayes.\" Advances in neural information processing systems 31 (2018).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "## **Response to Reviewer QSK9**\n\nThank you for your thoughtful review of our work. We address your main points below:\n\n## **Weaknesses**\n### **W1: Theoretical analysis and subspace recovery**\nWe acknowledge the limitations in our theoretical analysis. However, sisPCA makes minimal assumptions about the data, and its objective function provides straightforward insights into the model's behavior and optimization landscape (Conjecture 3.1 and Appendix D).\n\n#### **Regarding the simulation example**\n\n1.\tWe confirmed that PCA-based models yield similar results with Gaussian-drawn mixing matrices compared to the uniform case. This is because the mixing matrix was normalized to ensure equal contribution from each subspace.\n2.\tAs discussed in Appendix D (lines 449-453), subspace scale and supervision strength are the major factors influencing the optimization landscape. For example, sisPCA will prioritize S2’s structure and ignore other subspaces if we scale up the S2 supervision (X, Y) large enough.\n\n#### **General subspace recoverability depends on supervision quality**\n\n1.\tIn the most intertwined case of identical supervision for two subspaces, only one can be recovered and the other will collapse to rank zero (Figure 6, Appendix D). Here, linearity ensures that the same information is not split between subspaces - it's an all-or-nothing scenario.\n2.\tWhen both supervised and unsupervised subspaces are presented, their dynamics are complicated by the dual role of HSIC supervision indicated in Eq.4, Appendix B.3. In particular, the supervised subspace will also try to expand upon the direction that maximize the (unsupervised) variance, potentially leading to identifiability or multi-optima issues (Section 3.3, lines 190-193, Figure 7). Though we reason that the same information will again concentrate in one subspace due to linearity.\n### **W2: Comparison with state-of-the-art baselines**\nWe appreciate your suggestion to include more state-of-the-art methods in our comparison. In response:\n\n1.\tWe've added comparisons with non-linear VAE-based models (**M2**), focusing on the real single-cell data where the practical need is to learn genes responsible for malaria infection defense (a feature selection task).\n2.\tOur analysis shows that while VAE-based models can extract independent subspaces, they lack a straightforward method for identifying genes upon which each subspace was constructed (**M2.2**).\n3.\tRegarding representation quality, we've confirmed that linear models like sPCA serve as strong baselines, both qualitatively and quantitatively (**see Fig S1, Table S2 in M2.3**). Specifically, quality of the linear infection subspace often surpasses its non-linear counterpart.\n\n## **General questions**\n### **Q1: Effect of non-linear kernels for target variables**\nYes, the choice of target kernel does influence sisPCA's results, primarily affecting the rank of the learned subspace. This is because sPCA's subspace dimensions are determined by supervision kernel rank (Appendix B.3, line 406-408), and sisPCA shows similar behavior. Specifically, for a linear kernel on *K* continuous features, the effective dimension of the corresponding sisPCA subspace is the rank of the sum of a rank-K kernel and a rank-D disentanglement kernel, where *D* depends on other subspace dimensions (line 146-147). A non-linear target kernel can effectively expand the learned subspace while maintaining interpretability, albeit losing the connection with regularized linear regression (Appendix B.3).\n\n### **Q2: Clarification on Fig. 2 and maximization objectives**\nWe apologize for the confusion. To clarify:\n\n- A and B represent different target variables used for supervision. \n- The maximization arrow within each subspace represents a PCA-like objective of variance maximization, which is equivalent to maximizing the HSIC towards an identity kernel (Appendix B.1, line 397-398). \n- In supervised PCAs, the variance maximization objective comes from the HSIC to the target (first component of Eq.2, see Eq.4 in Appendix B.3). Upon decomposition, the arrows between the subspace and A/B indeed represents minimizing the prediction error (first component of Eq.4), and the arrows within each subspace represents maximizing the variance (second component of Eq.4). \n\nWe will revise the figure to make notation consistent.\n\n### **Q3: Definition of balanced/unbalanced supervision**\nBalanced/unbalanced supervision refers to the relative ratio of scale/norm between target kernels (discussed in Appendix D, line 449-473). It can alternatively be defined as the ratio of loss gradient with respect to each supervision kernel. Intuitively speaking, unbalanced scenarios occur when one kernel is scaled up, favoring the corresponding subspace and creating a single global optimum (Figure 6 in Appendix D).", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "Sp7bWsqmZg", "reviewer": "Reviewer_debZ", "summary": "Here the authors propose an extension of PCA that decomposes datasets into multiple independent subspaces that are encouraged to reflect provided covariates of interest. To enforce these dependencies, the authors leverage the Hilbert-Schmidt Independence Criterion (HSIC).", "questions": "Did the authors try applying their method with different values of $\\lambda$ to the real-world datasets?", "limitations": "The greatest limitations I see are (1) the method relying on linear transformations (as opposed to e.g. deep neural networks) and (2) potential difficulties in setting $\\lambda$. For (1) the authors specifically note that they are intentionally trading off expressivity for increased interpretability", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "In short, my views of the manuscript are generally positive. Some specific thoughts below:\n\n* **Clarity**: I found the manuscript well written and easy to follow. Well done!\n* **Motivation**: A recent line of work [1,2] has observed that disentangling separate sources of variation in case-control datasets (e.g. in scRNA-seq studies [3-4]) can assist with data exploration tasks. This work extends similar ideas to settings where no explicit case versus control datasets exist, but samples can still be associated with different covariates of interest. I can imagine this work will similarly prove useful in a number of data exploration tasks.\n* **Novelty**: To my knowledge, the authors proposed method is indeed novel, and seems like a sensible extension of the previous supervised PCA method and recent works on supervised disentanglement (e.g. [1-2]).\n\n[1] Abid et al., \"Contrastive principal component analysis\"\n[2] Abid et al., \"Contrastive variational autoencoder enhances salient features\"\n[3] Weinberger et al., \"Isolating salient variations of interest in single-cell data with contrastiveVI\"\n[4] Jones et al., \"Contrastive latent variable modeling with application to case-control sequencing experiments\"", "weaknesses": "While my views of the manuscript are generally positive, I do have a couple of minor concerns. If the authors are able to address my concerns I would be happy to raise my score.\n\n* **Impact of $\\lambda$**: The authors' method attempts to find a balance between encouraging dependence between subspaces and their corresponding covariates while also minimizing the dependence between individual subspaces, and this balance is tuned via $\\lambda$ in Equation (1). I imagine this parameter could have a significant impact on the embeddings/loadings produced by sisPCA, though this isn't explored in the manuscript. It would be great if the authors could perform additional experiments exploring these changes. For example, do the GO results presented in Figure 5 vary for different values of $\\lambda$? How different do the UMAPs look for the scRNA-seq data for different values of $\\lambda$?\n\n* **Automatic selection of $\\lambda**: Related to my previous point, the authors briefly mention that $\\lambda$ could be selected using a procedure similar to [1]. It would be great to see an example of this in practice: forcing the user to select $\\lambda$ manually seems like it could present a big obstacle in terms of usability, so I think it would be useful to see how the automated method performs.\n\n* **Re: connections with self-supervised learning**: In the abstract the authors note \"We elucidate the mathematical connections of\nsisPCA with self-supervised learning and regularized linear regression\". While I was able to locate the regularized linear regression results in Appendix B.1, I couldn't find any additional mathematical details on the self-supervision connection (the only other mention of self-supervision in the main text was in line 83 in the related works section). Could the authors expand upon this connection? (Please let me know if I missed anything here).\n\nMinor points (these did not affect my score):\n\n* It appears that the authors may have used \\citet accidentally instead of \\citep throughout the manuscript.\n* Figure quality was relatively low-resolution when I printed out the manuscript- perhaps the authors could try increasing the DPI?\n* In line 91 the authors mention \"contrastive learning\"- it may be worth putting explicit citations to e.g. Abid et al. [1] to clarify for the reader that \"contrastive learning\" here isn't referring to contrastive in the sense of e.g. SimCLR.\n* I couldn't find any details on how the real-world biological datasets were preprocessed (e.g. I assume library size normalization was applied to the scRNA-seq data?)\n\n[1] Abid et al., \"Contrastive principal component analysis\"", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "## **Response to Reviewer debZ**\n\nThank you for your positive feedback and insightful comments! We appreciate your thorough review and would very much like to address any remaining concern.\n\n### **W1. Impact of $\\lambda$ selection and an auto-selection pipeline**\nThank you for the suggestion. We conducted additional experiments by varying $\\lambda$ for the scRNA-seq dataset with values [0, 0.1, 0.3, 1, 3, 10, 30, 100]. On top of that, we implemented a preliminary pipeline for automatic $\\lambda$ selection that aggregate resulting subspaces using clustering and outputs aggregated results at different resolution, as implemented in Contrastive PCA. However, due to time constraints, we have yet to thoroughly test some design choices, such as the representation distance, which appears to be important in subspace clustering. To speed up training, we also introduced a mini-batch approach of Alg 1 (Appendix B.2) since the learnable projection matrix *U* is of size (n_feature, n_dim), independent of n_sample. We plan to implement additional features such as parallel training and better initialization search in the future.\n\nDespite these limitations, we confirmed the following results:\n\n#### **Impact on subspace similarity**\nSubspaces learned with similar $\\lambda$ values are generally more similar, as measured by Grassmann distance between representations, indicating gradual and continuous changes as a function of $\\lambda$. The infection subspaces for $\\lambda$ = [0.1, 0.3, 1] cluster together, as do the time subspaces for $\\lambda$ = [0, 0.1], while the other subspaces appear more detached.\n\n#### **Impact on gene selection and GO analysis** \nWe quantified overlaps of the top 100 genes list for the infection subspace (ranked by absolute PC1 weights) across different $\\lambda$ values. Results for $\\lambda$ = [0, 0.1, 0.3, 3, 10] are highly consistent, sharing more than 90 of the top 100 genes. The result also suggests that sPCA ($\\lambda$ = 0) is sufficient for the GO analysis in Figure 5 ($\\lambda$ =10). We examined sPCA representations and confirmed that the confounding effects from the temporal subspace (as visualized in Figure 9b and the new Fig S1c in **M2**) are mainly concentrated in PC2, which is not used here for feature selection. Notably, for reasons still under investigation, $\\lambda$ = 1 is an outlier, sharing only 68 of the top 100 genes with sPCA. A possible explanation is that Alg 1 (Appendix B.2) on large datasets may experience numerical issues with PyTorch’s SVD solver.\n\nWe will add visualizations of these results in the final version.\n\n### **W2. Connection with self-supervised learning**\nWe apologize for any ambiguity. Our reference of \"self-supervised learning\" is the connection to the auto-encoder (Appendix B.1), where target variables - features from the data - are used to disentangle the data itself. We chose the term mainly to reflect the difference with supervised learning. Here our goal isn't perfect target prediction, as labels are already known. Instead, we aim to better explain the data by reweighting the self-reconstruction loss with supervision guidance (as reflected in line 401). We recognize the potential confusion and will clarify the term in the final version. The connection to regularized linear regression is discussed in Appendix B.3.\n\n### **Minor points**\nThank you for pointing out the problems. \n- We will fix citation issues, and will increase the DPI of all figures to improve resolution.\n- We will provide more details on the preprocessing of the single-cell infection data. The data, preprocessed by the original authors, was subject to library size normalization followed by log1p transformation. A background correction step was also applied before normalization, where the mean expression in empty wells was subtracted. This makes the raw counts technically not integers, but the negative binomial likelihood in hsVAE-sc (**see M2**) seems to work just as well.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:
We thank the reviewers for their thoughtful comments. We first address the two main points raised by multiple reviewers, followed by responses to individual comments.

### **M1: Highlight the technical novelty (Reviewers PWKM and UeQu)**
We appreciate your acknowledgment of our paper’s technical contribution as novel and sensible (Reviewer debZ). To further clarify our design insights:
- **Linearity and interpretability by design.** PCA remains popular due to its interpretability. It is particularly valuable in data exploration, where the goal is to identify features that best represent the data in specific aspects. In this regard, linear models are indeed the state-of-the-art (See **M2.2**). We motivated our work as a multi-subspace extension of PCA to fully leverage its interpretability. It is the first linear supervised disentanglement model and is competitive with or even outperform non-linear counterparts (See **M2.3**). 
- **Theoretical advantages from linearity.** We show that sisPCA can be viewed both as linear auto-encoder (Remark 3.1 and Appendix B.1) and as regression for continuous targets (Remark 3.3 and Appendix B.3). Linearity also allows for more efficient and reliable optimization (Algorithm 1 and Conjecture 3.1).

### **M2: Comparison with non-linear SOTA methods (Reviewers UeQu and QSK9)**
We have now included additional baselines based on the HSIC-constrained VAE (HCV) ([1] of Reviewer QSK9).

### **M2.1 SOTA implementation**
We reimplemented the idea from [1] using the latest SCVI framework (scvi-tools v1.1.2) for variational inference training and designed three non-linear counterparts of PCA models:
1. **VAE**: Vanilla VAE with Gaussian likelihood.
2. **sVAE**: VAE with additional predictors for target variables.
3. **hsVAE**: sVAE with additional HSIC penalty.

Non-linear models generally have more hyperparameters. In our benchmark, we fixed the VAE architecture following the SCVI default (e.g. one layer of NN with 128 hidden units and batch normalization), or the scVIGenQCModel in [1] (e.g. equal weighting of prediction and reconstruction losses).

**Table S1: General model comparison**
| | Linear ||| Non-linear |||
|:--|:--:|:--:|:--:|:--:|:--:|:--:|
| | PCA | sPCA | sisPCA (this work) | VAE | sVAE | hsVAE |
| Supervision | - | HSIC | HSIC | - | NN prediction | NN prediction|
| Disentanglement | - | - | HSIC | - | - | HSIC |
| Interpretation | Linear projection *U* as feature importance ||| Blackbox ||
| Hyperparameters | 1) #dim | 1) Subspace #dim. 2) Target kernel choice | 1) and 2) of sPCA. 3) HSIC penalty | 1) #dim. 2) General NN design. | 1) and 2) of VAE. 3) Predictor design | 1), 2) and 3) of sVAE. 4) HSIC penalty|
| Optimization | Closed form | | Simple landscape | General limitations of NN and VI like robustness** ||

\** We were not able to run VAEs on the 6-dimensional simulated data in Figure 3 due to NaNs generated during variational training, which is not uncommon for SCVI models and is likely the result of exploding parameters.


### **M2.2 Interpretability**
We define and compare model interpretability based on the practical need to *learn how each feature contributes to specific aspects of the high-dimensional data*. 
Linear models, as demonstrated in Section 4.3.3, are inherently interpretable because the learned projection *U* directly shows how each subspace axis relates to original features. PCA-based models have the additional advantage of ordering subspace axes based on importance (variance explained). In Figure 5, we selected subspace-specific genes based on PC1 loading scores, which would not be possible without the mapping *U*.
In contrast, VAE models link features and subspace axes through a non-linear black box. While some interpretation approaches such as gradient-based saliency map are becoming standard practice, they are nowhere near as straightforward as the linear mapping. In addition, since gradient and loss can only be calculated per sample, aggregating them into a global feature importance score is non-trivial.

### **M2.3 Results on the single-cell infection data**
We have now included results from VAE-based models (**see attached Fig S1 for visualization**). We also provide below quantitative metrics in **Table S2**. While higher scores generally indicate better results, we caution that these metrics may not fully capture the biological relevance of the representations. As in Section 4.3.2, our analysis aims to isolate genes directly involved in parasite harboring from those associated with broader temporal changes post-infection. A better temporal subspace should tell us about subtle progression (pseudo-time) even of the same collection time beyond prediction accuracy (since we already knew the labels).

**Table S2: Quantitative evaluation of representation quality**

| Method | Linear ||| Non-linear ||||
|:--|:--:|:--:|:--:|:--:|:--:|:--:|:--:|
| | PCA | sPCA | sisPCA (this work) | VAE | sVAE | hsVAE | hsVAE-sc |
| Subspace separateness-Grassmann | 0 | 3.802 | 4.467 | 0 | 3.510 | 3.797 | 3.459 |
| Information density-Silhouette: (infection, time) | (0.068, 0.083) | (0.313, 0.279) | (0.294, 0.164) | (0.045, 0.137) | (0.058, 0.215) | (0.064, 0.197) | (0.235, 0.353) |
| Information density after UMAP-Silhouette: (infection, time) | (0.009, 0.238) | (0.153, 0.364) | (0.197, 0.228) | (-0.022, 0.335) | (0.008, 0.602) | (0.047, 0.586) | (0.258, 0.603) |

Overall, we've confirmed that:

- HSIC against target efficiently imposes supervision, achieving performance equivalent to neural networks and making sPCA indeed a strong baseline. 
- sisPCA outperforms its general-purpose VAE counterparts, particularly in the infection subspace. However, VAEs can gain further SOTA performance from domain-specific knowledge (such as count-based modeling).
- Disentanglement may come at the cost of weakening supervision, aligning with Section 4.2.3 and Table 2.

We will include the above results (**Fig S1, Table S1 and S2**) in the final version.

Meta-review: All four reviewers argued to accept the paper, albeit three of them only weakly. However, the paper has improved quite a bit relative to initial reviewer concerns during the discussion period, particularly noting the additional experiment which broadens the scope from the focus on single-cell RNA sequencing data; and the study into automatic settings of \lambda and the sensitivity to the parameter choice. Multiple reviewers raised their score after the author response.

Please be certain to incorporate everything discussed in the rebuttal phase when preparing the final version of the manuscript.
Decision: Accept (poster)