Factorized linear discriminant analysis for phenotype-guided representation learning of neuronal gene expression data
Authors: Mu Qiao
Markus Meister
Venue:  (2021)
Keywords: factorized linear discriminant analysis
phenotype
gene expression
representation learning

Abstract: A central goal in neurobiology is to relate the expression of genes to the structural and functional properties of neuronal types, collectively called their phenotypes. Single-cell RNA sequencing can measure the expression of thousands of genes in thousands of neurons. How to interpret the data in the context of neuronal phenotypes? We propose a supervised learning approach that factorizes the gene expression data into components corresponding to individual phenotypic characteristics and their interactions. This new method, which we call factorized linear discriminant analysis (FLDA), seeks a linear transformation of gene expressions that varies highly with only one phenotypic factor and minimally with the others. We further leverage our approach with a sparsity-based regularization algorithm, which selects a few genes important to a specific phenotypic feature or feature combination. We applied this approach to a single-cell RNA-Seq dataset of Drosophila T4/T5 neurons, focusing on their dendritic and axonal phenotypes. The analysis confirms results obtained by conventional methods but also points to new genes related to the phenotypes and an intriguing hierarchy in the genetic organization of these cells.

Under review as a conference paper at ICLR 2021
FACTORIZED LINEAR DISCRIMINANT ANALYSIS FOR PHENOTYPE-GUIDED REPRESENTATION LEARNING OF NEURONAL GENE EXPRESSION DATA
Anonymous authors Paper under double-blind review
ABSTRACT
A central goal in neurobiology is to relate the expression of genes to the structural and functional properties of neuronal types, collectively called their phenotypes. Single-cell RNA sequencing can measure the expression of thousands of genes in thousands of neurons. How to interpret the data in the context of neuronal phenotypes? We propose supervised learning approach that factorizes the gene expression data into components corresponding to individual phenotypic charac- teristics and their interactions. This new method, which we call factorized linear discriminant analysis (FLDA), seeks linear transformation of gene expressions that varies highly with only one phenotypic factor and minimally with the oth- ers. We further leverage our approach with sparsity-based regularization algo- rithm, which selects a few genes important to a specific phenotypic feature or feature combination. We applied this approach of Drosophila T4/T5 neurons, focusing on their dendritic and axonal phenotypes. The analysis confirms results obtained by conventional methods but also points to new genes related to the phenotypes and an intriguing hierarchy in the genetic organization of these cells.
single-cell RNA-Seq dataset
INTRODUCTION
The complexity of neural circuits is a result of many different types of neurons that specifically connect to each other. Each neuronal type has its own phenotypic traits, which together determine the role of the neuronal type in a neural circuit. Typical phenotypic descriptions of neurons include features such as dendritic and axonal laminations, clectrophysiologica properties, and connectiv- ity (Sanes & Masland, 2015 Zeng & Sanes 2017] Gouwens et al., 2019). However, the genetic programs behind these phenotypic characteristies still poorly understood.
Recent progress in characterizing neuronal cell types and investigating their gene expression, espe- cially with advances in high-throughput single-cell RNA-Seq (Zeng & Sanes 2017), provides an opportunity to address this challenge. With massive data generated from single-cel RNA-Seq. we now face computational problem: how to factorize the high-dimensional data into gene expression modules that are meaningful to neuronal phenotypes? Specifically, given phenotypic descriptions of neuronal types, such as their dendritic stratification and axonal termination, can one project the original data into a -iw-dimensional space corresponding to these phenotypic features and their interactions. and further extract genes critical to each of these components?
Here we propose a new analysis method named factorized linear discriminant analysis (FLDA). In- spired by multi-way analysis of variance (ANOVA) (Fisher 1918), this method factorizes data into components corresponding to phenotypic features and their interactions, and seeks a linear transfor- mation that varies highly with one specific factor but not with the others. The linear nature this approach makes it easy to interpret, the weight coefficients directly inform the relative importance of each gene to each factor. We further introduce a sparse variant of the method, which constrains the number of genes contributing each linear projection. We illustrate this approach by applying FLDA to a single-cell transcriptome dataset of T4/T5 neurons in Drosophila (Kurmangaliyev et al. 2019) focusing on two phenotypes: dendritic location and axonal lamination.

---

Under review as a conference paper at ICLR 2021
A
Complete table
B
Partial
representing the interaction of both feetures
Phenotypic feature
Phenotypic feature
j=C
j=b
j+2
[-3
J-0
j+1
j=2
j=3
1=0
Phenotopic festure
Phenotypic featuree
1=0
i=0
Rag
Tei
Pez
j=l
1-1
Ria
Pait
J=D
representing Fosture
y: representing
Figure I: Illustration of our approach. (A.B) the example, cell types are jointly represented by two phenotypic features, indexed with labels and respectively. If only some combinations of the two features are observed, one obtains a partial contingency table (B) instead of a complete one (A). (C) We seek linear projections of the data that separate the cell types in a factorized manner corresponding to the two features. Here u, v, and w are aligned with Feature 1, Feature 2, and the interaction of both features, with the projected coordinates y, 2, and respectively.
FACTORIZED LINEAR DISCRIMINANT ANALYSIS (FLDA)
Suppose that given gene expression data of single neurons which typically very high- dimensional. These cells are classified into cell types, as a result of clustering in the high- dimensional space and annotations based on prior knowledge or verification outcome (Macosko et al. 2015 Tasic et al. 2016 Shekhar et al. 2016 Tasic et al. 2018 Peng et al. 2019). We know the phenotypic traits each neuronal type, therefore each type can also jointed defined by the phenotypic features. We want to find interpretable low-dimensional embedding in which certain dimensions represent factors of phenotypic features or their interactions. This requires that variation along one of the axes in the embedding space causes the variation of only one factor. In reality, this is hard to satisfy due to noise the data, and we relax the constraint by letting data projected along axis vary largely with one factor while minimally with the others. In addition, we ask that cells classified as the same type are still close to each other in the embedding space, while cells of different types are far apart.
As a start, let us consider only two phenotypic features of neurons, dendritic stratification, and axonal termination, both of which be described with discrete categories, such different regions or layers in the brain (Oh et al., 2014) Buler et al. 2014) Sanes & Masland 2015 Kurmangaliyev al. 2019). Suppose that each cell type can be jointly represented by dendritic location indexed as and axonal lamination indexed as j, with the number of cells within each cell type naj. This representation can be described using a contingency table (Figure [TA.)). Note here that allow the table to be partially filled.
Let Ziyk 2, ...nn) represent the expression values of g genes in each cell (æxyk E R2)). How find linear projections Vijk = u" ziyk and 2iyk v" Zijk that are aligned with features and respectively (igure][)) We first asked whether we could factorize, for example, Bıyk. with respect to components depending on features i and j. Indeed, motivated by the linear factor models used in multi-way ANOVA and the idea of partitioning variance, we constructed an objective function as the following, and found 2² that maximizes the objective (see detailed analysis Appendix A);
W NAu
2u" arg max uER
u M.u
When we have a complete table, and there are a levels for the feature and levels for the feature j. we have
NA = MA A XyMAB
(2)
where MA MR. and MAB are the covariance matrices explained by the feature i, the feature j. and the interaction of them. X and X2 param-eas controlling the relative weights of M
2

---

Under review as a conference paper at ICLR 2021
and MAB with respect MA- Me is the residual covariance matrix representing noise in gene expressions. Formal definitions of these terms are the following:
MA
m, m mm.)m m..)
(3)
a
==
a
MB
m, mm)(m - m.)"
(4)
i=1 1=1
b
MAB
Mij m. m. + mm.)(iijj m, - -m. + m.
(5)
(a - () 1)
i=1 j=1
1
M
Bigh - m.:)( - -m;j)²
(6)
N ab
Thj k=1
1=1 /=1
where
1
m
mij
ab
j=1
m.
m
(8)
j=1
&
m.j
mij
(9)
a
in which
mij
Bijk
(10)
Tkjj
An analogous expression provides the linear projection o² for the feature j. and w" for the inter- action of both features and j. Similar arguments can be applied to the scenario of partial table find or ex the linear projection for the feature i or j (see Appendix B for mathematical details).
Note that N A is symmetric and M, positive definite. Therefore the optimization problem is generalized eigenvalue problem (Ghojogh et al. 2019). When M, is invertible, u² is the eigenvector associated with the largest cigenvalue M, "N.. In general, if we want to embed Zifk into a d- dimensional subspace aligned with the feature (d a). largest eigenvalues of M. N NA. which we call the top factorized linear discriminant components (FLDs). Since multi-way ANOVA handle contingency tables with more than two dimensions, our analysis can easily generalized to more than two features.
take the eigenvectors with the d
SPARSITY-BASED REGULARIZZTION OF FLDA
For this domain specific application in neurobiology. there is particular interest finding a small group of genes that best determine of the phenotypic features. This leads finding axes that have only a few non-zero components. To identify such sparse solution, we solved the following optimization problem:
3

---

Under review as a conference paper at ICLR 2021
u NAu u M,
u' arg max MER
subject to
||u|lo I
(11)
from which the number of non-zero elements of u² is less equal to I.
This is known as a sparse generalized eigenvalue problem, which has three unique challenges, as listed in Tan al. (2018): first, when the data are very disensiinniinnl,, Mc can be singular and non-invertible; second, because of the normalization term u² M.u, many solutions for sparse eigenvalue problems cannot be applied directly; finally, this problem involves maximizing convex objective over a nonconvex set, which is NP-hard.
To solve it, we used truncated Rayleigh flow (Rifle), method specifically developed to solve sparse generalized eigenvalue problems. The algorithm of Rifle is composed of two steps (Tan et al. 2018): first, to obtain an initial vector Ug that is close to u". We used the solution from the non-sparse FLDA as an initial estimate of un: second, iteratively, to perform a gradient ascent step on the objective function, and then execute truncation step that preserves the / entries of with the largest values and sets the remaining entries to O. Pseudo-code for this algorithm can be found in Appendix C.
As proved in Tan et al. 2018), if there is a unique sparse leading generalized eigenvector, Rifle will converge linearly to with the optimal statistical rate of convergence. The computational com- plexity of the second step Olgg + 9) for each iteration, therefore Rifle scales linearly with g. the dimensionality of the original data. Based on the theoretical proof, to guarantee convergence, the hy- perparameter 17 was selected to be sufficiently small such that 77maaz (M,) < 1, where Xmaz (Mc) is the largest eigenvalue of M.. In our case, the other hyperparameter 1, indicating how many genes to be preserved, was empirically selected based on the design of a follow-up experiment. As men- tioned later in Results, chose / be 20, reasonable number of candidate genes be tested in a biological study.
RELATED WORK
4.1 LINEAR DIMENSIONALITY REDUCTION
FLDA is a method for linear dimensionality reduction (Cunningham & Ghahramani 2015). For- mally, linear dimensionality reduction defined as the following: given n data points each of g dimensions, X = [x1 22 ] E R and choice of reduced dimensionality r < g. optimize an objective function f() to produce a linear projection P R and Y PX E Rxx the low-dimensional transformed data.
State-of-the-arr methods for linear dimensionality reduction include principal component analy- sis (PCA), factor analysis (FA), linear multidimensional scaling (MDS), linear discriminant anal- ysis (I.DA), canonical correlations analysis (CCA), maximum autocorrelation factors (MAF), slow feature analysis (SFA), sufficient dimensionality reduction (SDR), locality preserving projections (LPP), and independent component analysis (ICA) (Cunningham them are obviously unsuitable the problem want to solve, for example, MAF and SFA developed for data with temporal structures (Larsen 2002 Wiskott & Sejnowski, 2002), and LPP focuses the local structures of the data instead of their global organization (He & Niyogi). The remaining approaches can be roughly grouped for either unsupervised or supervised linear dimen- sionality reduction, and discuss them separately.
Ghahramani 2015). Some of
4.2 UNSUPERVISED METHODS FOR LINEAR DIMENSIONALITY REDUCTION
Unsupervised linear dimensionality reduction, including PCA (Joliife 2002), ICA (Hyvärinen) ct al. 2200], FA (Spearman) 1904) and more, project data into a low-dimensional space without using supervision labels. They are not suitable to solve our problem, due to the specific charac- teristics of the gene expression data. The dimensionality of gene expression data is usually very high, with tens of thousands of genes, and expressions of fair number of them can be very noisy. These noisy genes cause large variance among individual cells, but in an unstructured way. Without supervision signals from the phenotypic features, unsupervised methods tend to select these genes
4

---

Under review as a conference paper at ICLR 2021
to construct the low-dimensional space, which offers neither the desired alignment good sep- aration of cell type clusters. To illustrate this, we performed PCA on the gene expression data and compared it with FLDA. Briefly, we solved the following objective to find the linear projection:
u 6'XX"u 2
u arg max uER
(12)
u
The outcome of this comparison is shown in Results.
4.3 SUPERVISED METHODS FOR LINEAR DIMENSIONALITY REDUCTION
Supervised linear dimensionality reduction, represented by LDA Fiiherr 1936 Mc McLachlan 2004)
and CCA [Hotelling 1936] Wang et al.) (016). can overcome the above issue. By including super- vised signals of phenotypic features, we can devalue genes whose expressions are non-informative about the phenotypes.
4.3.1 LINEAR DISCRIMINANT ANALYSIS (LDA)
We name method FLDA because its objective function has similar format as that of LDA. LDA also models the difference among data organized in pre-determined classes. Formally, LDA solves the following optimization problem:
U Esu
u arg max
(13)
uER> U T E.Eu
where Eg and Eg are estimates of the between-class and within-class covariance matrices respec- tively.
Different from FLDA, the representation of these classes not explicitly formulated as contin- gency table composed of multiple features. The consequence is that, when applied to the example problem in which neuronal types are organized into a two-dimensional contingency table with phe- notypic features and j. in general, axes from LDA are not aligned with these two phenotypic features.
However, in the example above, the axes from each LDA to align with its specific feature. We call this approach "2LDAs". There two limitations of this approach: first, it discards information about the component depending on the interaction of the two features which cannot be explained by linear combination of them; second, it explicitly maximizes the segregation of cells with different feature levels which sometimes is not consistent with a good separation of cell type clusters. Detailed comparisons between LDA, "2LDAs" and FLDA can be found in Results.
perform two separate LDAs for the two features. This allows
4.3.2 CANONICAL CORRELATION ANALYSIS (CCA)
CCA projects two datasets X, E R and X, Rdxn the correlation between Y, and Y, is maximized. Formally, it tries to maximize this objective:
dxxn to Y. E R and Y, Rexn, such that
u u'(X,X.) 1 XX (X,X,)= tv
0²) arg max DERA,DER.
(14)
(ufu) - v)
To apply CCA to our problem, need to set X. to be the gene expression matrix, and X, to be the matrix of d phenotypic features (d = 2 for two features as illustrated later). In contrast with FLDA. CCA finds transformation of gene expressions aligned with linear combination of phenotypic features, instead of factorization of gene expressions corresponding to individual phenotypic features. This difference is quantified and shown in Results.
5

---

Under review as a conference paper at ICLR 2021
Feature
Genes
Ovent nettc
Overall mod mod.aattyy scom such
|+D
Typky Typta Typlos
2.0
685
Fiature
I+0
Typz
Type
158
123
i+1
Typz
Type
Calle
0.00
633
(2
RLA
Data Syntheus
635
5amm
Figure 2: Quantitative comparison between FLDA and other models. (A) Illustration of data syn- thesis. See Appendix D for implementation details. Color har indicates the expression values of the ten generated genes. (B) Normalized overall SNR metric of each analysis. The SNR values are normalized with respect to that of LDA. (C) Overall modularity score for each analysis. Error bars in (B,C) denote standard errors each calculated from 10 repeated simulations.
EXPERIMENTS
5.1 DATASETS
In order to quantitatively compare FLDA with other linear dimensionality reduction methods, such as PCA. CCA, LDA, and the "2LDAs" approach, we created synthetic datasets. Four types of cells, each containing 25 examples, were generated from Cartesian product of two features and j. organized in a 2x2 complete contingency table. Expressions of 10 genes were generated for these cells, in which the levels of Genes 1-8 were correlated with either the feature i, the feature j. the interactions them, and the levels the remaining genes were purely driven by noise (Figure A Details of generating the data be found in Appendix D.
To illustrate FLDA in analyzing single-cell RNA-Seq datasets for real problems of neurobiology, and demonstrate the merit of our approach in selecting few important genes for each phenotype, we used a dataset of Drosophila T4/T5 neurons (Kurmangaliyev et al. 2019). T4 and T5 neurons very similar in terms of general morphology and physiological properties, but they differ by the location of their dendrites in the medulla and lobula, two distinct brain regions. T4 and T5 neurons each contain four subtypes, with each pair of the four laminating their axons a specific layer in the lobula plate (Figure BA) Therefore, we can use two phenotypic features to describe these neurons: the feature indicates the dendritic location at the medulla or lobula; the feature ý describes the axonal lamination at one of the four layers (a/b/c/d) (Figure BB). In this experiment, we focused on the dataset containing expression data of 17492 genes from 3833 cells collected defined time during brain development.
5.2
DATA PREPROCESSING
The T4/T5 neuron datasct was preprocessed as previously reported (Shekhar et al.) 2016 Peng
al. 2019 Tran al. 2019). Briefly, transcript counts within each column of the count matrix igenesxcelly were normalized to sum to the median number transcripts per cell, resulting in the normalized counts Transcripts-per-media
TPM for Gene g in Cell We used the log- transformed expression data Ege = In PT + 1) for further analysis. We adopted a common approach in single-cell RNA-Seq studies that is based on fitting a relationship between mean and coefficient of variation (Chen et al. 2016:] Pandey et al.) 2018) to select highly variable genes, and performed FLDA on the expression data with only these genes. We preprocessed the data with PCA and kept principal components (PCs) explaining 999 of the total variance before running FL.DA but not the sparse version of the algorithm. the experiment below, we set the hyper-parameters λs in Equation to 1.
5.3 METRICS
We included the following metrics to evaluate our method: A signal-to- -o-is ratio (SNR) measures how well each discriminant axis separates cell types compared with noise estimated from the vari-

---

Under review as a conference paper at ICLR 2021
SuB matric sach
153
T4 merrons
T5 nncurons
Arcnol cemération
Devenne beatin
medala
Tile
1ib
Tdt:
Tid
HANE 200
ubula
15x
15b
15c
15d
7OU.3
FOR SOLESAARRRPPPER WOP.EEEE
150
192-
30.
F0.1
-15
-200
-150
-200
220
FLD
[LD
F1.
Figure 3: FLDA on the dataset of T4/T5 neurons. (A) T4/T5 neuronal types and their dendritic and axonal phenotypes. (B) T4/T5 neurons can be organized complete contingency table. Here i indicates the dendritic location and indicates the axonal termination. (C) SNR metric of each discriminant axis. (D) Projection of the data into the taree-dimensione space consisting of the discriminant axis for the feature (FLD,) and the first and second discriminant axes for the feature j (FLD;1 and FLD,2). (E-G) Projection of the data into the two-dimensional space made of FLD; and FLD511 (E), FL.D11 and FLDj2 (F), or FLD,2 and FLD (the third discriminant axis for the feature j) (G). Different cell types are indicated by different colors as in (A) and (D).
ance within cell type clusters. The explained variance (EV) for each discriminant axis measures how much variance the feature or j is explained among the total variance explained by that axis. The mutual information (MI) between each discriminant axis and each feature quantifies how "in- formative" an axis to specific feature. Built on the calculation of MI, we included the modularity score which measures whether each discriminant axis depends on at most one feature (Ridgeway & Mozer 2018). The implementation details of these metrics can be found in Appendix E.
RESULTS
To quantitatively compare the difference between FLDA and other alternative models including PCA, CCA, LDA, and "2LDAs"", we measured the proposed metrics from analyses of the synthe- sized datasets (Figure 2 Given that the synthesized data table, each LDA of the "2LDAs" approach could find only one dimension for the specific feature i or Therefore, as a fair comparison, we only included the corresponding dimensions in FLDA (FLD; and FLD,) and the top two components of PCA, CCA, and LDA. The overall SNR values normalized by that of LDA and the overall modularity scores were plotted for data generated with different noise levels (Figure BB. The performance of PCA is the worst among all these models because the unsupervised approach cannot prevent the noise from contaminating the signal. The supervised approaches in general have good SNRs, but LDA and CCA suffer from low modularity scores. This is expected because LDA maximizes the separation of cell type clusters but overlooks the alignment of the axes the feature or j. and CCA maximizes the correlation to linear combi- nation of phenotypic features instead of individual ones. By contrast, "2LDAs" achieves the highest modularity scores but has the worst SNR among the supervised approaches, because tries to max- imize the separation of cells with different feature levels, which is not necessarily consistent with maximizing the segregation of cell types. Both the SNR value and the modularity score of FLDD
organized in a 2x2 contingency
7

---

Under review as a conference paper at ICLR 2021
close to the optimal, as it not only considers the alignment of axes to different features but also constrains the variance within cell types. A representative plot of the EV and MI metrics of these models is shown in Figure []] reporting good alignment of axes to either the feature i or j in FLDA and '2LDAs", but not in the others.
A question of significance in neurobiology is whether the diverse phenotypes of neuronal cell types are generated by combinations of modular transcriptional programs, and if so, what is the gene sig- nature for each of the programs. To illustrate the ability of our approach in addressing this problem, we applied FLDA to the dataset of Drosophila T4/T5 neurons. The T4/T5 neurons could be orga- nized in a 2x4 contingency table, therefore, FLDA was able to project the expression data into a subspace of seven dimensions, with one FLD aligned with dendritic location i (FLD,). three FLDs aligned with axonal termination j (FI.D,1-3), and the remaining three representing the interaction of both phenotypes (FL.Di1-3). We ranked these axes based on their SNR metrics and found that FLD,1. FL.D,, and FLD)22 j2 have much higher SNRs than the rest (Figure (FC). Indeed, data repre- sentations in the subspace consisting of these three dimensions show a clear separation of the eight neuronal types (Figure BD). As expected, FLD; leases apart T4 from TS neurons, whose dendrites are located at different brain regions (FigureBE), Interestingly. FLD,1 separates T4/T5 neurons into two groups, a/b vs c/d, corresponding to the upper or lower lobula place, and FLDy2 divides them into another two, a/d b/c, indicating whether their axons laminate at the middle or lateral part of the lobula plate (Figure F... Unexpectedly. among these three dimensions, FL.D11 has a much higher SNR than FLD; and FLD whose SNR values are similar. This suggests a hicrarchical structure in the genetic organization of T4/T5 neurons: they are first separated into either a/b or c/d types, and subsequently divided into each of the eight subtypes. In fact, this exactly matches the sequence of their cell fate determination during development, as revealed in previous genetic study (Pinto- Teixeir et al.] 2018). Finally, the last discriminant axis of the axonal feature FLD separates the group a/c from b/d, suggesting its role in fine-tuning the axonal depth within the upper or lower lobula plate (Figure 3G).
To seek gene signatures for the discriminant components in FLDA, we applied the sparsity-based regularization to constrain the number of genes with non-zero weight coefficients. Here we set the number to 20, reasonable number of candidate genes that might be tested in follow-up biological study. We extracted a list of 20 genes each for the axis of FLD, or FLD,1. The relative importance of these genes to each axis is directly informed by their weight values (Figure C).)). Side-by bide-by-side, we plotted expression profiles of these genes in the eight neuronal types (Figure BDD)) For both axes, the genes critical in separating cells with different feature levels are differentially expressed in corresponding cell types. We compared our gene lists with those obtained using conventional methods which were reported in Kurmangaliyev et al. (219)) Consistent with the report, we found indicator genes for dendritic location such as TYAP-2, 2. 2pr2, CG34155, and CG12065, and those for axonal lamination such as klg, bi, pros. In addition, we found genes that were not reported in this previous study. For example, our results suggest that the genes Thor and pHC1-1 1 are impor- tant to the dendritic phenotype, and Lae and Mip are critical to the axonal phenotype. These are promising genetic targets to be tested in biological experiments. Lastly, FLDA allowed us to exam- ine the component that depends on the interaction of both features and identify its gene signature, which provides clues to transcriptional regulation of gene expressions in the T4/T5 neuronal types (Figures 6 and
As a supervised approach, FLDA depends correct phenotype labels to extract meaningful in- formation. But if the phenotypes are annotated incorrectly, can we use FLDA to raise a flag? We propose a perturbation analysis of FLDA to address this question built on the assumption that among possible phenotype annotations, the projection of gene expression data based on correct labels leads to better metric measurements than incorrect ones. As detailed in Appendix F, we generated three kinds of incorrect labels for the dataset of T4/T5 neurons, corresponding to three common scenarios of mislabeling: the phenotypes of cell type were mislabeled with those of another type; singular phenotypic level was incorrectly split into two; two phenotypic levels are incorrectly merged into one. FLDA was applied to gene expressions of T4/T5 neurons but with these perturbed annotations. Proposed metrics such as the SNR value and modularity score were plotted in Figure g Indeed, the projection of gene expressions with correct annotation leads to the best SNR value and modularity score compared with incorrect annotations. This implies that this type of perturbation analysis is useful practice in general: it raises the confidence that the original annotation is correct if FLDA on the perturbed annotations produces lower metric scores.
8

---

Under review as a conference paper at ICLR 2021
Her
Toor
Loe
Lac
Abb
0042505
CR4224
C0554
00964
C014664
0.8
w
C410001
0032065
P0p
200-
-02
0034355
0006715
0033522
Agupho
pannetins
POO
40
49
Selected for the dendritic phenotype
Selected genes for the axonal phenotype
Figure 4: Critical genes extracted from the sparse algorithm. (A) Weight vector of the 20 genes selected for the dendritic phenotype (FLD,). The weight value indicated in the color bar with color indicating direction (red: positive and green: negative) and saturation indicating magnitude. (B) Expression patterns of the 20 genes from (A) in eight types of T4/T5 neurons. Dot size indicates the percentage of cells in which the gene was expressed, and color represents average scaled expression. (C) Weight vector of the 20 genes selected for the axonal phenotype (FLD; Legend in (A). (D) Expression patterns of the 20 genes from (C) in eight types of T4/T5 neurons. Legend as in (B).
7 DISCUSSION
We developed FLDA, a novel supervised linear dimensionality reduction method for understanding the relationship between high-dimensional gene expression patterns and cellular phenotypes. We illustrate the power of FLDA by analyzing gene expression dataset from Drosophila T4/T5 neurons that are labeled by two phenotypic features, each with multiple levels. The approach allowed us to identify new genes for each of the phenotypic features that were not apparent under conventional found hicrarchical relationship in the genetic organization of these cells. These findings point the way for new biological experiments.
methods. Furthermore,
The approach is motivated by multi-way ANOVA, and thus it generalizes easily to more than two features. Future applications in neurobiology include the analysis of phenotypic characteristics such as electrophysiology and connectivity (Zeng & Sanes 2017 Gouwens al.) 2019 2020). More generally FLDA can be applied to any labeled data set for which the labels form a Cartesian product of multiple features. For example, this would include face images that can be jointly labeled by the age, gender, and other features of person (Moghaddam & Ming- -iung Yang 2002 Zhang et al., 2017).
FLDA factorizes gene expression data into features and their interactions, and finds a linear pro- jection of the data that varies with only one factor but not the others. This provides modular representation aligned with the factors (Bengio et al., 2012). Ridgeway & Mozer (2018) argued that modularity together with explicitness could define disentangled representations. Our approach is linear, which presents an explicit mapping between gene expressions and phenotypic features, therefore our approach can potentially serve as a supervised approach to disentanglement (Kingma
et al. 2014 Kulkarni et al. 2015 Karaletsos et al. 2016).
Compared with other non- linear embedding methods for cell types (Hinton & Salakhutdinov 2006 Feng et al. 2014; Gala et al. 2019), the linear nature of FLDA makes it extremely easy to interpret the low-dimensional representations, as the weight vector directly informs the relative importance of each gene. To allow the selection of small set of critical genes, we leveraged our approach with sparse regularization. This makes FLDA especially useful to experimentalists who can take the list of genes and test them in subsequent round of experiments.

---

Under review as a conference paper at ICLR 2021
8 APPENDIX
8.1 A. OBJECTIVE FUNCTIONS
Here we derive the objective functions used in our analysis. Again if æiys 1, ..,ig))) represents expression values of genes in each cell (Zijk E RR)), we seek to find a linear projection Viyk u² Tiyk that aligned with the feature
8.1.1 INSPIRATION FROM ANOVA
We asked what is the best way to factorize Yiyk- Inspired by multi-way ANOVA Fiisher 1918), identified three components: one depending on the feature i, another depending on the feature j. and the last one depending the interaction of both features. We therefore followed the procedures of ANOVA to partition sums of squares and factorize Vijk into these three components.
Let us first assume that all cell types defined by and contain the same number of cells. With cell types represented by a complete contingency table (FigureTA), Mijk can be linearly factorized using the model of two crossed factors. Formally, the linear factorization is the following:
Yijk a, + B3 (a:)) + Cigh
(15)
where Yigk represents the coordinate of the kth cell in the category defined by i and j; f is the average level of y: as is the component that depends on the feature i, and B, the component that depends on the feature j; (a:);; describes the component that depends on the interaction of both features i and ji Eijk ~ (0, o²) is the residual of this factorization.
Let say that the features and j fall into and b discrete categories respectively. Then without loss of generality, we can require:
0g = 0
(16)
Σ Bj =0
(17)
Σ (o iy
Σ (a:3)
(18)
Corresponding these, there are three null hypotheses:
Ho1 :a 0
(19)
Hoa Bg 0
(20)
Hos (aB)sz =
(21)
Here we want to reject Hoı while accepting Huz and Hos in order that Yifk is aligned with the feature i.
Next, we partition the total sum of squares. If the number cells within each cell type category is n, and the total number of cells is N. then we have
10

---

Under review as a conference paper at ICLR 2021
Wiik V...)² bn
(k.. y...) + (..j.
D...)
1=1
(22)
b
+ 71
Vig. Vi.. J. + F...)
Uijk
Dig. )²
km
where ÿ is the average of Vijk over the indices indicated by the dots. Equation
22 can be written as
SST SSA + SSp ASAB SS.
(23)
with each term having degrees of freedom - a-,, 6- (a-1)(b 1)(b-)), 1), Here SSA. SSB. SSAB, and SS, are partitioned sum of squares for the factors Pj. (a.5)ij. and the residual.
and - ab respectively.
ANOVA rejects or accepts null hypothesis by comparing its mean square (the partitioned sum of squares normalized by the degree of freedom) that of the residual. This is done by constructing F-statistics for each factor as shown below:
MSA
SS
FA
a-1
(24)
MS.
SS.
N-
8Sa
MSB
FB
b-1
(25)
MS.
58.
N ab
SSAB (a 14-1X6-11
MSAB
FAB
(26)
MS.
SS.
N-
Under the null hypotheses, the F-statistics follow the F-distribution. Therefore, null hypothesis rejected when we observe the value of F-statistic above a certain threshold calculated from the F-distribution. Here want FA to be large enough so that we can reject Hol. but FB and FAB to small enough for us to accept Hoz and Ho3. In other words, we want to maximize FA while minimizing F3 and FAB- Therefore, propose maximizing objective L:
= FA AB X2FAB
(27)
where A1 and A2 are with FA-
determining the relative weights of FB and FAB compared
8.1.2 OBJECTIVE FUNCTIONS UNDER COMPLETE CONTINGENCY TABLE
When the numbers of cells within categories defined by and j (nay) are not all the same, the total sum of squares cannot be partitioned as in Equation (22). However, we only care about distinctions between cell types instead of individual cells, we can use the mean value of each cell type cluster (Tig.) to estimate the overall average value (7...), and the average value each category i (Ti..) or
(. Therefore, Equation 22 can be modified the following:
b
6
(Bijk - J...)1-- > (Js. 72.)) + > (U.,. - )
1=1
j=1
a 6
Bag
Diy. Úr. 8.g. + )²
(yia Diy. Dig.) )21
7
(28)
11

---

Under review as a conference paper at ICLR 2021
where
Ek Vijk
Vij.
(29)
naj
Σg 2j=1 Dij.
Dr.
(30)
> -i Dij.
ÿ.j.
(31)
2 2im1 2j=1 Dij.
9..
(32)
ab
If we describe Equation 28 as:
58T = =S SSAB SS.
(33)
then following the same arguments, we want to maximize an objective function the following format:
$8
X1
S8n b-1
X2
SS (a-1)
e-1
(34)
58 ab
8.1.3 OBJECTIVE FUNCTIONS UNDER A PARTIAL CONTINGENCY TABLE
When we have a representation of a partial table, we can no longer separate out the component that depends on the interaction both features. Therefore, we use another model, linear model two nested factors, to factorize Vijk. which has the following format:
Diyk = 1 Ot + By(5) + Eiyk
(35)
Note that we now have By18) instead of B, + (or3) In this model, we identify a primary factor, for instance, the feature denoted by which falls into categories, and the other (indexed by j) becomes a secondary factor, the number of whose levels b, depends on the level of the primary factor. We merge the component depending on the interaction of both features into that of the secondary factor as Bz(t)-
Similarly, we have
D
myy
bg (Tt.. - 09..)9
(iikk - T...).
i=1 j=1
(36)
h,
6:
(Di.. 32.
Wjikk
Dzj. 11
j=1
which can be written as
SS, SSA + S8B + SS.
(37)
with degrees of freedom N 1, M a, and N M for each of the terms, where M is:
12

---

Under review as a conference paper at ICLR 2021
M
ba
(38)
Therefore, want maximize the following objective:
A
$8 M
L
a-1
(39)
SS
8.2 FLDA WITH A PARTIAL CONTINGENCY TABLE
Here we provide the mathematical details FLDA under the representation a partial table. When have partial table, if the feature i is the primary feature with levels, and the feature is the secondary feature with b, levels, then NA in Equation (T) defined as follows:
NA = MA AM
(40)
where
MA
m; m..)m, - -m..)
(41)
i=
by MB
mij m..)(mij m.
(42)
M
and M is defined in Equation
38
Correspondingly, M. in Equation is defined as:
1
M.
[iihh maj)( m.;j)²1
(43)
N
M
Thj)
and
b
m.
mij
(44)
M
(=1
1=1
b;
m.
m
(45)
b
The remaining mathematical arguments are the same those for the complete table. In this scenario, because we don't observe all possible combinations of features i and j. we cannot find the linear projection for the interaction of both features.
8.3 C. PSEUDO-CODE FOR RIFLE
We show pseudo-code for the Rifle algorithm as follows:
13

---

Under review as a conference paper at ICLR 2021
procedure RIFLE(NA,, A, Me, (.)
Y) is the step size
==1
indicates the iteration number Converge when u u
while not converge do
14 NA
Pt-1 u. M.
C I 1 Pt-1
(NA
Pc-1
Ou
Ut f 1[Cu--ill
Truncate ur by keeping the top entries of with the largest values and setting the rest
entries to 0
Ut 4 <- +-1+1 end while return ue end procedure
u1
8.4 D. IMPLEMENTATION DETAILS DATA SYNTHESIS
To quantitatively compare FLDA with alternative approaches, we synthesized data of four cell types, each of which contained 25 cells. The four cell types were generated from a Cartesian product of two features i and j. where i E {0. 1) and j E (0, 1]. Expressions of 10 genes were generated for each cell. The expression value of the hth gene in the kth cell of the cell type ij. g was defined the following:
Fisk = i Eijk
(46)
Ri Eijk
(47)
R3 inj + Eigh
(48)
9yA iVj. J + Eigh
(49)
Pijk 21 + Gifk
(50)
P3Kk = 2j + Eijk
(51)
0354 21 AJ J + Gigt
(52)
Lit 2i V j Eigk
(53)
R1 Eisk
(54)
Figh Cijk
(55)
where
[1, if 1 -
i iAj
(56)
0,
otherwise
14

---

Under review as a conference paper at ICLR 2021
and
iVJ = f0, if 0,j=0 = [1,
(57)
otherwise
were interactions of the two features. Here Etyk driven by Gaussian noise, namely,
Ejjk 0O,
(58)
We synthesized datasets of different a values (a E {0.2, 0.4, 0.6, 0.8, 1.0)). This was repeated 10 times and metrics for each a value were calculated as the average across the 10 repeats.
8.5
IMPLEMENTATION DETAILS OF THE METRICS USED IN THE STUDY
We measured the following metrics in experiments:
8.5.1
SIGNAL-TO-NOISE RATIO (SNR)
Because we care about the separation of cell types, we define the SNR metric as the ratio of the variance between cell types over the variance of the noise, which is estimated from within-cluster variance. For the entire embedding space, given cell types, if the coordinate of each cell is indicated by c, then we define the overall SNR metric as the following:
tr tr np(Cp. C.p. 2)
SN Roverall =
(59)
tr(Eg Eh (Cpk (Cpk )²)
where Cp. is the center of each cell type cluster, and é is the center of all data points. Let denote the embedded coordinate along a specific dimension. The SNR metric for that axis is therefore:
1 Mppp. - C2 5 " (Cpi Cp.)) p=
SNR
(60)
8.5.2 EXPLAINED VARIANCE (EV)
We want to know whether the variation of a specific dimension is strongly explained by that of specific feature. Therefore, measure, for each axis, how much of the total explained variance is explained by the variance of the feature or j. Formally, given the embedded coordinate yijks we calculate the EV as the following:
Ee Disl
Rij Vi. D...)
EV; =
(61)
E ===
(Vishh - 82..)
Ei= Σ 1 M4y (J.. -2..)
EV,
(62)
Σ Σ
(Bijk J...²
where y is the average of yifk over the indices indicated by the dots.
8.5.3 MUTUAL INFORMATION (MI)
The MI between a discriminant axis u and a feature quantifies how much information of the feature is obtained by observing data projected along that axis. It is calculated the MI between data representations along the axis y = u" and feature labels of the data f. where X the original gene expression matrix:
15

---

Under review as a conference paper at ICLR 2021
T(y. f) H(y) H(S) H(y. f)
=- Σ p(u) log2 p(u) - Σ p(f (og)((s)
ΣΣ BEY JEF
Ply, f) p(v.
(63)
yEY
JEF
Here H indicates entropy. To calculate (z) and H (y, f), we discretize y into 10 bins.
8.5.4 MODULARITY
Ridgeway Mozer (018) argued that in a modular representation, each axis should depend on at most single feature. Following the arguments in their paper, the modularity score is computed as follows: we first calculate the MI between each feature and each axis (M, denotes the MI between axis and one feature f)- If an axis perfectly modular, it will have high mutual information for only one feature and zeros for the others, we therefore compute a template Lyy as the following:
0i f arg max, Mig
4f = 0, otherwise
(64)
where 0, maxy mig. We then calculate the deviation from the template
[y - tzs)²
8:
(65)
O2(N
where N is the number of features. The modularity score for the axis is 1- -d. The mean of over is defined as the overall modularity score.
6; 1
8.6 F. IMPLEMENTATION DETAILS OF ANNOTATION PERTURBATION
To evaluate the effect mislabeling phenotypic levels, we made use of the dataset of T4/T5 neurons, and generated three kinds of perturbation the original labels:
First, we switched the phenotype labels of T4a neurons with one of the seven other types (T4b, T4c, T4d, T5a, T5b, T5c, T5d). In this scenario, phenotype labels of cell types were incorrect, but the number of cell type clusters was the same. We had two levels of the dendritic phenotypes (T4/T5), and four levels of the axonal phenotypes (a/b/c/d). Therefore we kept one dimension for the dendritic feature, and three dimensions for the axonal feature.
Second, we merged the axonal phenotypic level with another level (b/c/d), an incorrect new level (a+b/a+c/a+td) In this scenario, we had three axonal phenotypes, therefore we kept two dimensions for the axonal feature.
Third, we randomly split each of the four axonal lamination labels (a/b/c/d) into two levels. For instance, among neurons with the original axonal level a, some of them were labeled with level al, and the others were labeled with a level a2. In this scenario, we had eight axonal phenotypes and we kept seven dimensions for the axonal feature.
We performed FLDA the dataset of T4/T5 neurons but with these perturbed annotations. Metrics from each of the perturbed annotations were measured and compared with those from original annotation.
REFERENCES
Yoshua Bengio, Aaron Courville, and Pascal Vincent. Representation Learning: A Review and New Perspectives. arXiv:1206.55:: [cs], June 2012.
Hung-L. Harry Chen, Yufang Jin, Yufei Huang, and Yidong Chen. Detection of high variability in gene expression from single-cell RNA-seq profiling. BMC Genomics, 17 Suppl 7:508, August 2016. ISSN 1471-2164, doi: 10.1 1811186ss122 116//128644011-22977722
16

---

Under review as a conference paper at ICLR 2021
John P. Cunningham and Zoubin Ghahramani. Linear Dimensionality Reduction: Survey. Insights, and Generalizations. Journal of Machine Learning Research, (6(89):2859-2900 2015. ISSN 1533-7928.
Thomas Euler, Silke Haverkamp. Timm Schubert, and Tom Baden. Retinal bipolar cells: Elementary building blocks of vision. Nature Reviews Neuroscience, 15(8):507-519 August 2014. ISSN 1471-0048. doi: 10.10/8/nm3378.
Fangxiang Feng. Xiaojie Wang, and Ruifan Li. Cross-modal Retrieval with Correspondence Au- toencoder. In Proceedings of the 22nd ACM International Conference on Multimedia, MM '14, pp. 7-16, Orlando, Florida, USA, November 2014. Association for Computing Machinery. ISBN 978-1-4503-3363-3. doi: 10. 145/2147868..... 2654992.
R. A. Fisher. The Use of Multiple Measurements in Taxonomic Problems. Annals of Eugenics, 7 (2)1799-188, 1936. ISSN 2050-1430 doi: (0.1111/j.1469-180) 180..19999..96666
Ronald Aylmer Fisher. The Correlation between Relatives on the Supposition of Mendelian Inheri- tance. Royal Society of Edinburgh], 1918.
Rohan Gala, Nathan Gouwens, Zizhen Yao, Agata Budzillo, Osnat Penn, Bosiljka Tasic, Gabe Mur- phy, Hongkui Zeng, and Uygar Sümbül. A coupled autoencoder approach for multi-modal analy- sis of cell types. In H. Wallach, H. Larochelle, A. Beygelzimer, F. d' textquuoteeee Alché-Bus E. Fox, and R. Garnett (eds.), Advances in Neural Information Processing Systems 32, pp. 9267- 9276. Curran Associates, Inc., 2019.
Benyamin Ghojogh, Fakhri Karray, and Mark Crowley. Eigenvalue and Generalized Eigenvalue Problems: Tutorial. arXiv: 11240 [cs, stat), March 2019. Comment: pages, Tutorial paper.
Nathan W. Gouwens, Staci A. Sorensen, Jim Berg, Changkyu Lee, Tim Jarsky, Jonathan Ting, Su- san M. Sunkin, David Feng, Costas A. Anastassiou, Eliza Barkan, Kris Bickley, Nicole Ble- sie, Thomas Braun, Krissy Brouner, Agata Budzillo, Shiella Caldejon, Tamara Casper, Dan Castelli, Peter Chong, Kirsten Crichton, Christine Cuhaciyan, Tanya L. Daigle, Rachel Dalley, Nick Dee, Tsega Desta, Son Ding, Samuel Dingman, Alyse Doperalski, Nadezhda Dotson, Tom Egdorf, Michael Fisher, Rebecca A. de Frates. Emma Garren, Marissa Garwood, Amanda Gary, Nathalie Gaudreault, Keith Godfrey, Melissa Gorham, Hong Gu, Caroline Habel, Kristen Hadley, James Harrington, Julie A. Harris, Alex Henry, DiJon Hill, Sam Josephsen, Sara Kebede, Lisa Kim, Matthew Kroll, Brian Lee, Tracy Lemon, Katherine E. Link, Xiaoxiao Liu, Brian Long. Rusty Mann, Medea McGraw. Stefan Mihalas, Alice Mukora, Gabe J. Murphy. Lind- say Ng, Kiet Ngo, Thue Nghi Nguyen, Philip R. Nicovich, Aaron Oldre, Daniel Park, Sheana Parry, Jed Perkins, Lydia Potekhina, David Reid, Miranda Robertson, David Sandman, Martin Schroedter, Cliff Slaughterbeck, Gilberto Soler-Llavina, Josef Sulc, Aaron Szafer, Bosiljka Tasic, Naz Taskin, Corinne Teeter, Nivretta Thatra, Herman Tung, Wayne Wakeman, Grace Williams, Rob Young, Zhi Zhou, Colin Farrell, Hanchuan Peng. Michael J. Hawrylycz, Ed Lein, Lydia Ng, Anton Arkhipov, Amy Bernard, John W. Phillips, Hongkui Zeng, and Christof Koch. Classifi- cation of electrophysiologicall and morphological neuron types in the mouse visual cortex. Nat. Neurosci., 22(7): 182 1195, July 2019. ISSN 1546-1726, doi: 10.1038/s41593-0 019 009---99-111--000...
Nathan W. Gouwens, Staci A. Sorensen, Fahimeh Baftizadeh, Agata Budzillo, Brian R. Lee, Tim Jarsky, Lauren Alfiler, Katherine Baker, Eliza Barkan, Kyla Berry, Darren Bertagnolli, Kris Bick- ley, Jasmine Bomben, Thomas Braun, Krissy Brouner, Tamara Casper, Kirsten Crichton, Tanya L. Daigle. Rachel Dalley, Rebecca A. de Frates, Nick Dee, Tsega Desta, Samuel Dingman Lee, Nadezhda Dotson, Tom Egdorf, Lauren Ellingwood, Rachel Enstrom, Luke Esposito, Colin Far- rell, David Feng, Olivia Fong, Rohan Gala, Clare Gamlin, Amanda Gary, Alexandra Glandon, Jeff Goldy, Melissa Gorbam, Lucas Graybuck, Hong Gu, Kristen Hadley, Michael J. Hawrylyez, Alex M. Henry, DiJon Hill, Madie Hupp. Sara Kebede, Tae Kyung Kim, Lisa Kim, Matthew Kroll, Changkyu Lee, Katherine E. Link, Matthew Mallory, Rusty Mann, Michelle Maxwell, Medea McGraw, Delissa McMillen, Alice Mukora, Lindsay Ng, Lydia Ng, Kiet Ngo, Philip R. Nicovich, Aaron Oldre, Daniel Park, Hanchuan Peng, Osnat Penn, Thanh Pham, Alice Pom, Zo- ran Popović, Lydia Potekhina, Ramkumar Rajanbabu, Shea Ransford, David Reid, Christine Ri- morin, Miranda Robertson, Kara Ronellenfitch, Augustin Ruiz, David Sandman, Kimberly Smith,
17

---

Under review as a conference paper at ICLR 2021
Josef Sulc, Susan M. Sunkin, Aaron Szafer, Michael Tieu, Amy Torkelson, Jessica Trinh, Her- man Tung, Wayne Wakeman, Katelyn Ward, Grace Williams, Zhi Zhou, Jonathan T. Ting. Anton Arkhipov, Uygar Sümbüll Ed S. Lein, Christof Koch, Zizhen Yao, Bosiljka Tasie, Jim Berg, Gabe Murphy, and Hongkui Zeng. Integrated Morphoelectric and Transcriptomic Classifica- tion of Cortical GeBger Cells. Cell, 183(4):335-953.) November 2020. ISSN 0092-8674.
doi: 10.
Xiaofei He and Partha Niyogi. Locality Preserving Projections. pp. 8.
G. E. Hinton and R. R. Salakhutdinov. Reducing the Dimensionality of Data with Neural Networks. Science, 31345786):004-077, July 2006. ISSN 0036-8075, 1095-9203. doi: 10. 10.1126/science. 1127647.
Harold Hotelling. RELATIONS BETWEEN TWO SETS OF VARIATES. Biometrika, 28(3-4): 321-377, December 1936. ISSN 0006-3444, doi: 10.10933hiomee2883-43311.
Aapo Hyvärinen, Juha Karhunen, and Erkki Oja. Independent Component Analysis. Interscience, New York, 1st edition edition, May 2001. ISBN 978-0-471-40540-5.
Wiley-
I. Jolliffe. Principal Component Analysis. Springer Series in Statistics. Ppringer-Verlage New York, second edition, 2002. ISBN 978-0-387-95442-4.
Theofanis Karaletsos, Serge Belongie, and Gunnar Rätsch. Bayesian representation learning with oracle constraints.
1.0.1.5066111111 [cs, staf), March 2016. Comment: 16 pages, publishes in
ICLR 16.
Durk P Kingma, Shakir Mohamed, Danilo Jimenez Rezende, and Max Welling. Semi-supervised Learning with Deep Generative Models. In Z. Ghahramani, M. Welling, C. Cortes, N. D. Lawrence, and K. Q. Weinberger (eds.), Advances in Neural Information Processing Systems 27. pp. 3581-3589. Curran Associates, Inc., 2014.
Tejas D Kulkarni, William F. Whitney, Pushmeet Kohli, and Josh Tenenbaum. Deep Convolu- tional Inverse Graphics Network. In C. Cortes, N. D. Lawrence, D. D. Lee, M. Sugiyama, and R. Garnett (eds.), Advances in Neural Information Processing Systems 28, pp. 2539-25477. Curran Associates, Inc., 2015.
Yerbol Z Kurmangaliyev, Juyoun Yoo, Samuel A LoCascio, and S Lawrence Zipursky. Modu- lar transcriptional programs separately define axon and dendrite connectivity. elife, 8:e50822, November 2019. ISSN 2050-084X. doi: 10.7554/eLiff.50822.
Rasmus Larsen. Decomposition using maximum autocorrelation factors. Journal of Chemometrics, (18--10)427-4335 2002. ISSN 1099-128X. doi: 00.1002/ccm.744
Evan Z. Macosko, Anindita Basu, Rahul Satija, James Nemesh, Karthik Shekhar, Melissa Gold- man, Itay Tirosh, Allison R. Bialas, Nolan Kamitaki, Emily M. Martersteck. John J. Trombetta, David A. Weitz, Joshua R. Sanes, Alex K. Shalek, Aviv Regev, and Steven A. McCarroll, Highly parallel genome-wide expression profiling of individual cells using nanoliter droplets. Cell, 161 (5):1202-122 1214, May 2015. ISSN 0092-8674. doi:
Geoffrey McLachlan. Discriminani Analysis and Statistical Pattern Recognition. Interscience, Hoboken, N.J, August 2004. iSBN 978-0-471-691155M
Wiley-
B. Moghaddam and Ming-Hsuan Yang. Learning gender with support faces. IEEE Transactions on Pattern Analysis and Machine Intelligence, 24(5):707-711 May 2002. ISSN 1939-3539. doi: 109/344 1000244.
Seung Wook Oh, Julie A. Harris, Lydia Ng, Brent Winslow, Nicholas Cain, Stefan Mihalas, Quanxin Wang, Chris Lau, Leonard Kuan, Alex M. Henry. Marty T. Mortrud, Benjamin Ouel- lette, Thue Nghi Nguyen, Staci A. Sorensen, Clifford R. Slaughterbeck, Wayne Wakeman, Yang Li, David Feng, Anh Ho, Eric Nicholas, Karla E. Hirokawa, Phillip Bohn, Kevin M. Joines, Hanchuan Peng. Michael J. Hawrylyez, John W. Phillips, John G. Hohmann, Paul Wohnoutka, Charles R. Gerfen, Christof Koch, Amy Bernard, Chinh Dang, Allan R. Jones, and Hongkui Zeng. A mesoscale connectome of the mouse brain. Nature, 508(7495):207 214, April 2014. ISSN 1476-4687. doi: 10. (0.1038/nature) 13186.
18

---

Under review as a conference paper at ICLR 2021
Shristi Pandey, Karthik Shekhar, Aviv Regev, and Alexander F. Schier. Comprehensive Identification and Spatial Mapping of Habenular Neuronal Types Using Single-Cell RNA-Seq, Curr. Biol., 28 1011052-1065.ee April 2018. ISSN 1879-0445. doi: 10. 1t1.16/f.cub22018002.000
Yi-Rong Peng, Karthik Shekhar, Wenjun Yan. Dustin Herrmann, Anna Sappington, Gregory S. Bryman, Tavé van Zyl, Michael Tri. H. Do, Aviv Regev, and Joshua R. Sanes. Molecular Classi- fication and Comparative Taxonomics of Foveal and Peripheral Cells in Primate Retina. Cell. 176 1111222-1137..22 February 2019. ISSN 1092-8674, doi: 10.
Filipe Pinto- Teixeira, Clara Koo, Anthony Michael Rossi, Nathalie Neriec, Claire Bertet, Xin Li, Al- berto Del-Valle-Rodriguez, and Claude Desplan. Development of Concurrent Retinotopic Maps in the Fly Motion Detection Circuit. Cell, 173(2):485 e9,.e11,, April 2018. ISSN 0092-8674.
doi: 10. w..1066f..ell.208..2222
Karl Ridgeway and Michael C Mozer. Learning Deep Disentangled Embeddings With the F-Statistic Loss. In S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa-Bianchi, and R. Garnett (eds.), Advances in Neural Information Processing Systems 31, pp. 185-194. Curran Associates, Inc., 2018.
Joshua R. Sanes and Richard H. Masland. The types of retinal ganglion cells: Current status and implications for neuronal classification. Annu. Rev. Neurosci. 38:221-246, July 2015. ISSN 1545-4126. doi: 10.1 -orpammurey--euu--0777----------
Karthik Shekhar, Sylvain W. Lapan, Irene E. Whitney, Nicholas M. Tran, Evan Z. Macosko, Monika Kowalczyk, Xian Adiconis, Joshua Z. Levin, James Nemesh, Melissa Goldman, Steven A. Mc- Carroll, Constance I.. Cepko, Aviv Regev, and Joshua R. Sanes. COMPREHENSIVE CLASSI-- FICATION OF RETINAL BIPOLAR NEURONS BY SINGLE-CELL TRANSCRIPTOMICS. August 2016. ISSN 0092-8674. doi: 10. 10.1016f..c01.2011077..55
Cell, (3665)::308-12333
C. Spearman. "General Intelligence," " Objectively Determined and Measured. The American Journal of Psychology, 15(2):201-292, 1904. ISSN 0002-9556. doi: 10.2307/1412107.
Kean Ming Tan, Zhaoran Wang, Han Liu, and Tong Zhang. Sparse generalized eigenvalue problem: Optimal statistical rates via truncated Rayleigh flow. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 80(5):1057-1086, 2018. ISSN 1467-9868. doi: 12291.
Bosiljka Tasic, Vilas Menon, Thuc Nghi Nguyen, Tae Kyung Kim, Tim Jarsky, Zizhen Yao, Boaz Levi, Lucas T. Gray, Staci A. Sorensen, Tim Dolbeare, Darren Bertagnolli, Jeff Goldy, Nadiya Shapovalova, Sheana Parry, Changkyu Lee, Kimberly Smith, Amy Bernard, Linda Madisen, Su- san M. Sunkin, Michael Hawrylycz, Christof Koch, and Hongkui Zeng. Adult mouse cortical cell taxonomy revealed by single cell transcriptomics. Nat. Neurosci, 19())3335-346, February 2016. ISSN 1546- 1726. doi: 10.1038/nn.421
Bosiljka Tasic Zizhen Yao, Lucas T. Graybuck, Kimberly A. Smith, Thue Nghi Nguyen, Darren Bertagnolli, Jeff Goldy, Emma Garren, Michael N. Economo, Sarada Viswanathan, Osnat Penn, Trygve Bakken, Vilas Menon, Jeremy Miller, Olivia Fong, Karla E. Hirokawa, Kanan Lathia, Christine Rimorin, Michael Tieu, Rachael Larsen, Tamara Casper, Eliza Barkan, Matthew Kroll, Sheana Parry, Nadiya V. Shapovalova, Daniel Hirschstein, Julie Pendergraft, Heather A. Sullivan, Tae Kyung Kim, Aaron Szafer, Nick Dee, Peter Groblewski, Ian Wickersham, Ali Cetin, Julie A. Harris, Boaz P. Levi, Susan M. Sunkin, Linda Madisen, Tanya I.. Daigle, Loren Looger, Amy Bernard, John Phillips, Ed Lein, Michael Hawrylyez, Karel Svoboda, Allan R. Jones, Christof Koch, and Hongkui Zeng. Shared and distinct transcriptomic cell types across neocortical areas. Nature, 563(7729):72-78, November 2018. ISSN 1476-4687. doi: 10. 1038//41186-008-0555444
Nicholas M. Tran, Karthik Shekhar, Irene E. Whitney, Anne Jacobi, Inbal Benhar, Guosong Hong, Wenjun Yan, Xian Adiconis, McKinzie E. Amold, Jung Min Lee, Joshua Z. Levin, Dingchang Lin, Chen Wang, Charles M. Lieber, Aviv Regev, Zhigang He, and Joshua R. Sanes. Single-cell profiles of retinal neurons differing in resilience to injury reveal neuroprotective genes. bioRxiv, pp. 71 711762, July 2019. doi: 10. 10.1101/71176..
Weiran Wang, Raman Arora, Karen Livescu, and Jeff Bilmes. On Deep Multi- View Representation Learning: Objectives and Optimization. arXiv: lrrXivv/602.000 /602.01024 [cs], February 2016.
19

---

Under review as a conference paper at ICLR 2021
Laurenz Wiskott and Terrence J. Sejnowski. Slow Feature Analysis: Unsupervised Learning of Invariances. Neural Computation, 14(4):715-770, April 2002. ISSN 0899-7667. doi: 10.1162/
0899766023133189988
Hongkui Zeng and Joshua Sanes. Neuronal cell-type classification: Challenges, opportunities and the path forward. Nature Reviews Neuroxcience, 18(9):30-554, September 2017. ISSN 1471-0048, doi: 10.1038/mm.2017.555
Zhifei Zhang, Yang Song, and Hairong Qi. Age Progression/Regresionn by Conditional Adversarial Autoencoder. arXiv: 000.700244444.. 772.08423 [cs], March 2017. Comment: Accepted by The IEEE Confer- ence on Computer Vision and Pattern Recognition (CVPR 2017).
ADDITIONAL INFORMATION
0.83
0.056
0.84
0.022
0.30
0.0
0.057
0.83
0.022
0.85
0.5
0.35
0.49
0.36
0.48
0.33
0.2
FUD_
FLD.J
LD_1
LD_
H
MM
0.084
0.8
0.4
FLD_]
F.D_
LD_
LDJ 101 LD2
0CA1
CCA2
PC1
PC2
0.0
Figure 5: Representative plots (at a = 0.6) EV and MI metrics for FLDA and other models. (A,B) EV (A) and MI (B) metrics of FLDA. FLD; and FLDD, indicate the factorized linear discriminants for features and j. (C,D) EV (C) and MI (D) metrics of 2LDAs. LD, and LD; indicate the linear discriminant components for features and j. (E,F) EV (E) and MI (F) metrics of LDA. LD, and LD2 indicate the first two linear discriminant components. (G,H) EV (G) and MI (H) metrics of CCA. CCA, and CCA2 indicate the first two canonical correlation axes. (IJJ) EV (I) and MI (J) metrics of PCA. PC₁ and PC2 indicate the first two principal components. EV; and EV, are the explained variance of features and j along an axis, and MI, and MI, indicate the mutual inform between an axis and features and j respectively, Values of EV and MI metrics are also indicated by the color bars on the right side.
20

---

Under review as a conference paper at ICLR 2021
100
Tda
T4b
73
Tdc
30
50
T5a
T4d
25
T3b
T5r
T5d
-75
FUD_L)
AUO
Figure 6: Additional plots for FLDA the dataset of T4/T5 neurons. (A.B) Projection of the orig- inal gene expression data into the two-dimensional space made of the first and second (FLDij1 and FLD;j2) (A) the second and third (FLD,j2 and FLD (B) discriminant axes for the compo- nent that depends on the combination of both features i and j. Different cell types are indicated in different colors in (B).
CC34350
Di-eaa
Sipe
0.
CG42673
CORN79
0.0
gan
(9070)
CG2617
Peddis
Rkon 14
-2
G417209
CG37200
Pape2
PM2
notie
End
End
Acc
Figure 7: Additional plots for critical genes extracted from the sparse algorithm. (A) Weight vector of the 20 genes selected for the interaction of both dendritic and axonal features (FLD,y1). The weight value is indicated in the color bar with color indicating direction (red: positive and green: negative) and saturation indicating magnitude. (B) Expression patterns of the 20 genes from (A) in eight types of T4/T5 neurons. Dot size indicates the percentage of cells in which the gene was expressed, and color represents average scaled expression.
21

---

Under review as a conference paper at ICLR 2021
0.5
Nemateed DVR
Nomatted DVL
0.6
0.6
Neraatedd
6.2
6.2
0.2
gmal
orggnal
types
3.0
1.0
Modulaity STOP
3.6
0.6
Modulrtt
0.2
Téc
753 15c
750
750
onginal
3+0
34C
3-5
9.0
orggral
Spitt
Type 14u
Typee
Spri Types
Figure 8: Evaluation of the effect of incorrect phenotype annotation the dataset of T4/T5 rons. (A.B) Normalized overall SNR metric (A) and overall modularity score (B) of FLDA after switching labels of T4a type with another neuronal type. (C,D) Normalized overall SNR metric (C) and overall modularity score (D) of FLDA after merging the axonal phenotypic level with another phenotypic level (b/c/d). (E.F) Normalized overall SNR metric (E) and overall modularity score (F) of FLDA after splitting each axonal phenotypic level into two. Metrics under the original annotation are colored green, and the values are indicated by the dashed lines. Here the SNR values are normalized with respect to that of the original annotation.
22


Reviews:
[{"review_id": "cRXFjoW1WQF", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "*Summary*\nThis manuscript presents a novel dimensionality reduction method, called Factorized Linear Discriminant Analysis. The method starts from a real problem in neurobiology, and tries to link expression levels of neural genes to phenotypes. In particular, the main goal of the proposed technique is to find linear projections of the genes expression which vary maximally with one phenotypical aspect and minimally with the others. The approach is evaluated using a synthetic example and a real case study involving Drosophila T4/T5 cells.\n\n\n*Positive points:*\n-> The paper starts from a real and challenging bioinformatics problem, i.e. the analysis of single-cell RNA sequencing data for neurobiology; data analysis tools for these data are nowadays fundamental to unravel the high complexity of this bio-medical field.\n-> The idea is interesting, well motivated and well explained.\n-> The proposed approach is simple and understandable: interpretability of solutions and results is currently fundamental when dealing with bio-medical data.\n-> The method is tested using a real world application\n\n\n*Negative points / questions:*\n\nCOMMENT 1. \nThe main problem of this manuscript is that the proposed method is not well inserted into the state of the art. Actually the discussion of authors of the related works is very limited, with only one paragraph at the bottom of page 1 plus the description of Linear Discriminant Analysis in Section 4. Many different linear dimensionality reduction techniques have been proposed in the past, each one with different characteristics, goals and optimization techniques. Authors should discuss them, especially in relation with their approach. A good entry point is the survey of Cunningham and Ghahramani:\n\nJohn P. Cunningham, Zoubin Ghahramani: Linear Dimensionality Reduction: Survey, Insights, and Generalizations, Journal of Machine Learning Research 16 (2015) 2859-2900\n\nAlso for what strictly concerns Linear discriminant Analysis, I think that a deeper discussion is needed. Many different extensions of LDA have been proposed, some of them strictly related to the goals/methods of this paper, like 2D-LDA:\n\nMing Li, Baozong Yuan, 2D-LDA: A statistical linear discriminant analysis for image matrix,\nPattern Recognition Letters, Volume 26, Issue 5, 2005, Pages 527-532, \n\nThis contextualization with respect to the state of the art is fundamental: without this, it is very difficult to get the true contribution of the proposed approach. In the same spirit, I suggest the authors to include some more recent techniques in the experimental comparison. \n\n\n\nCOMMENT 2.\nThe sparisification approach has not been fully discussed and justified. Many different methods for sparsification have been proposed, why do authors choose this particular one? How does this method relate to alternatives? Some sparse algorithms have been introduced also for discriminant analysis (not strictly related to the medical field), such as:\n\nN. H. Ly, Q. Du and J. E. Fowler, \"Sparse Graph-Based Discriminant Analysis for Hyperspectral Imagery,\" in IEEE Transactions on Geoscience and Remote Sensing, vol. 52, no. 7, pp. 3872-3884, July 2014, \n\n\nCOMMENT 3 .\nIf I correctly understand, in the experiments authors apply a PCA to reduce the gene expressions before applying the proposed method (last two lines of page 5). Is this a reasonable choice? I know that this is commonly done in other scenarios (like in Fisher-faces for face recognition), but in these experiments genes are fundamental for the knowledge extraction. I guess this is not done for the sparsified version (otherwise genes would have not been extracted), can you comment on this?\n\n\n\n\nCOMMENT 4.\nSome authors argued that formulating the optimization of linear dimensionality reduction techniques as eigenvalue or generalized eigenvalue problems is not always an adequate choice:\n\nJohn P. Cunningham, Zoubin Ghahramani: Linear Dimensionality Reduction: Survey, Insights, and Generalizations, Journal of Machine Learning Research 16 (2015) 2859-2900\n\nIn such paper the authors also suggested an alternative. Can you provide a comment on this?\n\n\nCOMMENT  5.\nI found some difficulties in reading and understanding the presentation and the discussion of the results. One problem was definitely the fact that tables are put in the appendix, and that such tables contains many numbers. I suggest the authors to put a summarizing table inside the manuscript (if possible).\n\n \nCOMMENT 6.\nThe phenotypical features are assumed to be categorical. How strict is this assumption? How much information are we loosing in this context? And in other contexts? In other words, are there other scenarios/applications in which this method be applied?  Adding some other possible application scenarios would increase the value of the proposal.\n\n\nCOMMENT 7.\nThis comment is related to the previous one, but contains more an open suggestion rather than a comment. Did you consider to use non categorical features? I think that this would open the usage of approaches used in three-way data analysis. Even I’m not aware of dimensionality reduction techniques methods for three-way data, I think a relation exists; techniques to extract interesting relations between different directions of data have been proposed, especially in the context of Tri-clustering, see for example:\n\nHenriques, R., Madeira, S.C.: Triclustering algorithms for three-dimensional data analysis: a comprehensive survey. ACM Comput. Surv. 51(5), 95 (2019)\n\nI suggest the authors to take a look also at this field.\n", "extra_scores": {}, "extra_text": {}}, {"review_id": "DJhKsWre9LE", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This manuscript describes a generalization of ANOVA that is intended to be used in the interpretation of single-cell RNA-seq data. The method requires specification of an orthogonal, discrete categorization of cells, nominally by phenotype.  The method then linearly factorizes the observed gene expression values into features and their interactions, relative to the phenotypic categories.  The factors can then be used to help interpret the categories, especially in conjunction with a regularizer to reduce the number of genes involved in the factors.\n\nI found this paper frustrating to read.  The second paragraph does not make clear exactly what problem is being addressed.  It says that we are \"given phenotypic descriptions of neuronal types,\" but not where those descriptions come from.  So I skipped ahead to the methods section to try to figure it out.  But even after reading that section, I could not understand where the phenotype values come from.  It was only when I made it to Section 2 that I verified that, indeed, the phenotypes are not observed but only inferred.\n\nHaving understood the problem, it seems to me that the use case for this approach is quite specific.  We need to have RNA-seq data that can be clustered in such a way that we can assign clusters to pre-defined phenotypic categories.  The claim, in the discussion section, that \"our approach can be easily generalized to .. additional characteristics such as electrophysiology and connectivity\" was not clear to me, but I assume this still refers to phenotypes that are inferred from the scRNA-seq data.\n\nThe critique leveled against CCA is that \"this approach cannot factorize gene expressions according to individual features, making the result hard to interpret.\"  But this seems like it must also be true of the proposed method, since prior to any analysis the data is transformed via PCA (Section 5.2).  I am confused, therefore, about how the method can be used to select genes in Section 6.\n\nThe comparison to LDA is done using two metrics, based on signal-to-noise ratio and mutual information.  I would have liked to hear more about why these particular metrics are appropriate, and in particular how they relate to whatever use case the authors have in mind. Overall, I am still not convinced that this is a problem that needs to be solved.  \n", "extra_scores": {}, "extra_text": {}}, {"review_id": "B1Uvhq5nQZy", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "Pros:\nThe paper studies a very important problem in gene data analysis. The proposed method is technically sound. The method is intuitive in its idea and easy to implement. The results are interpretable. And according to the experimental evaluations, the proposed method is consistent to existing biological observations and could further identify unknown genetic targets. Therefore, it, potentially, has insightful scientific implications. \n\nCons: \n- Relevance and Generalizability are unclear: The paper is relevant to researchers in subareas only and it is best-suited to a bioinformatics or neurobiology community. The paper requires background in neurobiological data analysis to evaluate whether the proposed method brings meaningful scientific insights.  In terms of novelties in machine learning field, the improvement over existing works, algorithmically or computationally, seems relatively limited. It is also unclear whether the proposed method generalizes well outside the subdomain of gene expression. It seems FLDA could possibly be applied to any such tensor data but a discussion on its general applicability to other data will be nice. \n\n- Relation to prior work is insufficient: Related work has not been discussed adequately, thus making the significance of this paper unclear. The authors discusses CCA and autoencoders in introduction, and also mentions the differences between FLDA and LDA (2LDAs), which are compared with the proposed method in experiments. I find such discussions are inadequate for the readers to understand the baseline or state of the art in this field. Therefore, It is somewhat unclear how the work improves from existing works. \n\n- Experimental evaluation is incomprehensive: As a largely application work, the comprehensiveness or tricks in experiments should be explained more clearly. For example, one part of the experiment, the sparsity-based regularization of FLDA seems like the application of Rifle (Kean Ming Tan et al) on the data. Readers are unclear what are the unique challenges in current setting. Also, any computational details (convergence, scalability, etc.) or any guidance on the hyperparameter selection? And does this method also generate meaningful results on another gene data?\n", "extra_scores": {}, "extra_text": {}}, {"review_id": "7wdl2t-3gNu", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "## Summary\nThe paper provides and ANOVA-inspired method, called FLDA, for creating a low dimensional embedding of scRNA-seq data. Additionally they propose a sparsity based method to find gene signatures, which can be used for further biological validation. The authors extensively evaluated their method on a data set of real expression values in Drosophila neurons.\nThey compared FLDA to two simpler and similar approaches, namely linear discriminant analysis (LDA) and a more feature-aligned version 2LDA.\nOn all benchmarks and metrics the FLDA shows clear advantage over the other methods.\n​\n## Reasons for score\nOverall, I vote for accepting. I can imagine several use cases of this approach, which looks like it is computational very feasible and easy to apply.\nAdditionally it is theoretically well founded on existing statistical methods.\nHowever, my main concern is FLDAs dependency on correct feature annotation, which is not mentioned by the authors.\nI hope authors can address my concern and the other cons in the rebuttal period.\n​\n## Pros\n- The authors clearly justify their design principle of FLDA based on ANOVA.\n- Benchmarks using different scores show clear advantage of the FLDA method.\n- The gene signature found for the T4/T5 data set show known marker genes for neuronal development.\n- Extension to the >2 feature case and other possible extensions are clearly outlined.\n- Extensive appendix showing theoretical validation for most obstacles in scRNA-seq data.\n\n## Cons\n- No evaluation on how dependent the model is on correct feature annotation.\n- No comparison to other methods for estimating gene signatures.\n- Validation on just one biological/real data set.​\n\n## Questions during rebuttal period\nPlease address and clarify the cons above.\nMainly what happens if one level is incorrectly defined, e.g. the ground truth are two or more levels for this singular level. Or in other words how robust is the method to the annotation.", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: The paper introduces a linear projection method, inspired by ANOVA,
for finding a supervised low-dimensional embedding.

A positive aspect is that the method is straightforward, and it is
even slightly surprising that in the family of linear models, there
still was an uncovered "niche".

The paper was considered useful for the purpose studied in the
paper, single-cell RNA-seq data analysis. But to claim broader
usefulness, more evidence should be presented.

One particular detail which was brought up by all reviewers was the
PCA preprocessing. For ICA it is a sensible choice, as linear ICA is
essentially "just" a rotation of the PCA components. But the
justification is not as good for a supervised method. PCA may be
necessary in practice, but may lose important category-relevant
information.

The paper still needs a significant revision before publication.  Even
though the method is straightforward method, a lot of time and
discussion was required for expert reviewers to understand it.

Decision: Reject