Structured Low-Rank Tensors for Generalized Linear Models
Authors: Batoul Ahmad Taki
Anand Sarwate
Waheed U. Bajwa
Venue: Accepted by TMLR ()

Abstract: Recent works have shown that imposing tensor structures on the coefficient tensor in regression problems can lead to more reliable parameter estimation and lower sample complexity compared to vector-based methods. This work investigates a new low-rank tensor model, called Low Separation Rank (LSR), in Generalized Linear Model (GLM) problems. The LSR model – which generalizes the well-known Tucker and CANDECOMP/PARAFAC (CP) models, and is a special case of the Block Tensor Decomposition (BTD) model – is imposed onto the coefficient tensor in the GLM model. This work proposes a block coordinate descent algorithm for parameter estimation in LSR-structured tensor GLMs. Most importantly, it derives a minimax lower bound on the error threshold on estimating the coefficient tensor in LSR tensor GLM problems. The minimax bound is proportional to the intrinsic degrees of freedom in the LSR tensor GLM problem, suggesting that its sample complexity may be significantly lower than that of vectorized GLMs. This result can also be specialised to lower bound the estimation error in CP and Tucker-structured GLMs. The derived bounds are comparable to tight bounds in the literature for Tucker linear regression, and the tightness of the minimax lower bound is further assessed numerically. Finally, numerical experiments on synthetic datasets demonstrate the efficacy of the proposed LSR tensor model for three regression types (linear, logistic and Poisson). Experiments on a collection of medical imaging datasets demonstrate the usefulness of the LSR model over other tensor models (Tucker and CP) on real, imbalanced data with limited available samples.

Published in Transactions on Machine Learning Research (08/2023)
Structured Low-Rank Tensors for Generalized Linear Models
Batoul Taki
bulonl..aki@@ruugrrrs.
Department of Electrical and Computer Engineering Rutgers University- New Brunswick
Anand D. Sarwate
anand.saruate@rutgerr. edu
Rutgers University-New Brunswiic
Waheed U. Bajwa
waheed. bajun Bruttenrt edu
Rutgers University- New Brunswick
Reviewed on OpenReview: https: // openreview. net/ forum? id= qUaBs3Ln41
Abstract
Recent works have shown that imposing tensor structures on the coefficient tensor in regres- sion problems can lead to more reliable parameter estimation and lower sample complexity compared to vector- based methods. This work investigates a new low-rank tensor model, called Low Separation Rank (LSR), Generalized Linear Model (GLM) problems. The LSR model which generalizes the well-known Tucker and CANDECOMP /AMPPPPAAAAC (CP) models, and is special case of the Block Tensor Decomposition (BTD) model is imposed onto the coefficient tensor in the GLM model. This work proposes a block coordinate descent algorithm for parameter estimation in LSR-structured tensor GLMs. Most importantly, it derives a minimax lower bound on the error threshold on estimating the coefficient tensor in LSR tensor GLM problems. The minimax bound proportional to the intrinsic degrees of freedom the LSR tensor GLM problem, suggesting that its sample complexity may be significantly lower than that of vectorized GLMs. This result can also be specialised to lower bound the estimation error in CP and Tucker-structured GLMs. The derived bounds comparable to tight bounds in the literature for Tucker linear regression, and the tightness of the minimax lower bound is further assessed numerically. Finally, numerical experiments on synthetic datasets demonstrate the efficacy of the proposed LSR tensor model for three regression types (linear, logistic and Poisson). Experiments on a collection of medical imag- ing datasets demonstrate the usefulness of the LSR model over other lensor models (Tucker and CP) on real, imbalanced data with limited available samples.
1 Introduction
In machine learning, regression models are used to understand the relationship between a set of independent variables (also known as covariates) and a dependent outcome. More formally, given a vector of covariates, x, and outcome, y, jointly distributed according Pxys the goal to produce a function that will predict y when given x. Under the Mean Squared Error (MSE) criterion, the solution would be finding the conditional mean Eyyx [x])). which is obtained by modeling the conditional probability Pyjx and estimating the model class parameters. For example, in linear regression, the predictor is given by b where the model class parameters are denoted by vector b and estimated through set of training data samples fx,, Vi } where x; is the :th sample vector of covariates. Different regression models suit different types of prediction problems: some common examples are linear, logistic and Poisson regression, all which fall under a broader class of models callod Generalized Lincar Models (GLMs). GLMs were introduced encompass classes of models that cannot be appropriately modeled as a simple 'lincar-repposee model' (McCullagh & Nelder, 2019). In particular, GLMs refer to parametric statistical framework that models the conditional

---

Published in Transactions on Machine Learning Research (08/2023)
probability of scalar response variable ys that follows an exponential family distribution. Because the exponential family distribution encompasses a large range of widely-used distributions, GLMs allow one to study broader class of regression problems.
GLM regression models are used for wide array of datasets in a multitude of applications. However, modern- day technologies are creating data-intensive environments and collecting increasingly high-dimensional data. In particular, we often encounter data in the form of variegated and structured multi-dimensional arrays (tensors), where the mumber of available data samples is far smaller than the mumber of variables (the dimensionality of the data). Prominent examples of two-dimensional arrays include biological imaging data such as clectroencophalogaaphy (EEG) and fiber-bundle imaging (Dumas et al., 2019). Examples of three- dimensional arrays include functional Magnetic Resonance Images (IMRI) (Bellec et al., 2017) and Magnetic Resonance Angiography (MRA) (Yang et al., 2020). Though such data has been used in many instances throughout the literature (Li et al., 2018; Hung & Wang, 2013; Zhou et al., 2013; Dumas et al., 2019), classical parameter estimation methods assume voctor-structured covariates and estimate a corresponding coelficient vector. There are two major concerns associated with multidimensional (a.k.a. tensor) data and its vectorization. First, vectorizing data that was originally in tensor form destroys its underlying structure (which often contains rich information valuable for regression analysis). Secondly, the resulting vector model exhibits a very large mumber of parameters, in the sense that in the high-dimensional setting the GLM regression model becomes ill-posed. This is an instance of 'the curse of dimensionalityy (Hung & Wang, 2013; Zhang& & Jiang, 2018).
A common solution to the curse of dimensionality in tensor GLM problems is to impose structure on the model parameters: this is the focus of this work. If the covariates are tensor structured, we expect to estimate coefficient tensors whilst exploiting the multidinensional structure of the data and rich information lying in the correlation between tensor modes. Common structures include: the addition of a sparsity regularizer or a low-rank inducing regularizer to the regression problem (Abramovich & Grinshtein, 2018; Seber & Lee, 2003; Zhang & Jiang, 2018; An & Zhang, 2020; Raskutti et al., 2019); the imposition of some tensor factorisation on the model parameters (Ahmed et al., 2020; Li el al., 2018; Zhou el al., 2013; Zhang et al., 2020; Tan et al., 2013; Zhang & Jiang, 2016; Wu al., 2022; Taki et al., 2021); or both of the above. The imposition of such structures should ultimately lead to the estimation of fewer parameters, improving the computational complexity and performance of our parameter estimation problem.
Two commonly used tensor factorisations are the CANDECOMP CEEECOMPPPAAAAAEE /PARAFAC (CP) and Tucker decomposi- tions (Kolda & Bader, 2009). These decompositions impose a compact structure on the coefficient tensors, thereby restricting the class of possible solutions in the parameter estimation problem. Compared to sim- ple vector regression, these decompositions can decrease the number of training samples needed for reliable coefficient estimation (also referred to as 'sample complexity"). A smaller sample complexity can lower the variance of the model, yet the restrictive nature of these decompositions can also reduce the representation power of the coefficient tensors for many classes of tensors, causing non-favourable bias-variance trad-eoff. A less studied decomposition particularly in regression works - is the Block Tensor Decomposition (BTD) (De Lathauwer, 2008), which can be expressed as a Tucker decomposition with block diagonal core tensor. We further discuss BTD in relation to our work in Section 3.
In this work we impose a new docomposition on the coefficient tensors that we will promptly refer to as the Low Separation Rank (LSR) model. In fact we will show that the LSR decomposition is generalization of the Tucker decomposition, and a special case of the BTD model (De Lathauwer, 2008). The LSR model maintains lower sample complexity than vectorization-best based tensor GLM regression but is less restrictive than Tucker or CP models of same sized core tensor. Though estimating an LSR-structured tensor introduces a greater sample complexity than the aforementioned decompositions, this increase is compensated by a stronger representation power, leading to better parameter estimation performance. In other words, we show that the increase in variance of the LSR model is conquered by its stronger representation power for a larger class of tensors (ergo, its decrease in bias), leading to more favourable bias- variance trade-off (model compactness representation power).
2

---

Published in Transactions on Machine Learning Research (08/2023)
1.1 Contributions
In this paper we make the following contributions. We introduce the LSR-structured tensor problem under the GLM framework that we appropriately denote as LSR-TGLM. GLMs encompass various flavours of regression including linear, logistic and Poisson. We focus on the high-dimensional setting and discuss the various parameters of an LSR-structured lensor (such as 'rank' and 'separation rank', terms we will introduce in Section 2) that reduce the sample complexity of parameter estimation in GLMs. We also compare sample complexities between different tensor dccomposition models and the LSR model.
Additionally, we explore two problems at the core of this work. First, we propose a parameter estimation algorithm (which we name LSRTR) for the LSR-IGLM problem. The main idea is that parameter estimation in GLMs can be achieved through Maximum Likelihood Estimation (MLE) (McCullagh & Nelder, 2019). However, for structured lensor settings, such as estimating CP or pocker-structured tensors, the objective function of the MLE problem is highly non-convex (Li et al., 2018; Zhang & Jiang, 2016; Zhou et al., 2013; Tan et al., 2013). This is also true when the tensor is LSR-ttuuctued. To overcome this, we observe that the problem can be partitioned into several convex sub-problems that can then be solved alternately. On the basis thereof we propose block coordinate descent. algorithm to find the MLE of the LSR-structured cocfficient tensor. Secondly, and perhaps most importantly investigate the fundamental error threshold of the LSR-TGLM problem by deriving minimax lower bound on the estimation error. This minimax bound is useful for assessing the performance of the proposed algorithm and ascertaining the parameters that may affect the sample complexity of the parameter estimation problem. The bound is general and can be specialised to previously introduced regression types under the GLM framework, such as CP and Tucker tensor GLMss Obtaining the bound requires a special construction of a packing set of T.R-structured tensors. The methods we develop are systematic and can be appropriate in other works that consider similar topological properties of structured tensors (i.c., LSR-structured tensors). We also assess the tightness of our bounds in two ways: 1) Through numerical study where we show that the ratio of the empirical error through LSRTR and the minimax bound is approximately constant with increasing sample size, and 2) We specialise our minimax bound to the Tucker linear regression case and show that our bound matches the optimal error rates for Tucker linear regression found in recent works (Zhang et al., 2020).
Finally, we evaluate the performance of our algorithm through extensive numerical experiments on synthetic data. We also test the performance of imposing the LSR structure on several classification problems with multidimensional medical imaging datasets. We show that while the LSR model outperforms the vector model, its rich representation power also allows for enhanced performance over the Tucker (and CP) case.
1.2 Relation To Prior Work
Regression problems have been a major focus of high-dimensional statistics for many years (Giraud, 2021). Some works on linear and logistic regression impose sparsily on the model parameter in order to reduce the sample complexity of the vector-based regression problems (Abramovich & Grinshtein, 2018; Sun & Zhang, 2012). However, in very high-dimensional regimes such as when the data is tensor-structured, i.c., we have (X, J sparsity assumptions do not provide enough reduction in the sample complexity (Raskutti et al., 2019; Lee & Courtade, 2020). Several works overcome the limitations of sparse vector regression by extending regular regression to the high-dimensional and low-rank matrix settings. Low-rank assumptions on data are common throughout the literature and are used to reduce the sample complexity of estimation problems (Barnes & Övgür, 2019; Shi et al., 2014). Such works propose regularized matrix linear and logistic regression models to obtain low-rank and/or sparse estimates of the cocfficient matrix regression problems, such as those on inference on images or graph data (Hung & Wang, 2013; Zhang & Jiang, 2018; Shi et al., 2014; An & Zhang, 2020; Berthet & Baldin, 2020). Some works directly impose low-rank structures cocfficient matrices through the rank-r singular value decomposition (SVD) (Taki et al., 2021).
Additionally, though tensors and their decompositions have long since been introduced in the literature (Kolda& & Bader, 2009), their applications in regression analysis have recently become established. Analogous to the low-rank matrix regression works, a varicty of works have introduced low-rank structures on cocfficient tensors for tensor regression problens. For logistic regression, Tan et al. (2013) first introduced using low- (CP) decomposition. A more flexible generalization of this
rank and sparse CANDECOMP/PARAAN
3

---

Published in Transactions on Machine Learning Research (08/2023)
work is imposing the Tucker decomposition on the coefficient tensor in tensor logistic regression (Zhang & Jiang, 2016; Wu et al., 2022). The Tucker structure has also been introduced for tensor linear regression (Zhang et al., 2020; Ahmed et al., 2020; Wu et al., 2022). To the best of our knowledge, the BTD structure has yet to be introduced for regression and GLMs.
More works to our interest generalize tensor linear and logistic regression works by imposing the CP and Tucker decompositions in tensor GLMs (Li et al., 2018; Zhou et al., 2013). Both structures have been shown to significantly reduce the number of learnable parameters, leading to efficient estimation and prediction in variety of regression problems, particularly with medical imaging data. The aforementioned works develop efficient estimation algorithms and provide empirical results on their performance. The proposed approaches outperform vector-based methods (in terms of estimation and prediction accuracy) in the high-dimensional setting when the number of available samples is limited. However, these matrix and tensor structures are aimed at being compact (in the sense that they decrease the number of learnable parameters in given problem), and are therefore also quite restrictive in their representation power of the true coefficient tensor. A more general and flexible tensor model is required to achieve accurale and efficient. estimation while maintaining useful level of compactness.
In terms of theoretical guarantees, various regression works provide local identifiability guarantees of the proposed CP and Tucker tensor models for GLMs, and asymptotic consistency and normality results for the MLE estimator of the model parameter (Li et al., 2018; Zhang et al., 2020; Zhou et al., 2013). Some works on higl-dimensional regression also provide sample complexity bounds of the proposed model in the form of risk upper bounds (Zhang et al., 2020; Ahmed et al., 2020) or minimax lower bounds (Barnes & Özgür, Foster et al., 2018; Abramovich & Grinshtein, 2018; Lee & Courtade, 2020; Abramovich & Grinshtein, 2016; Raskutti et al., 2011); however, these works are specific to restor-based logistic regression or Tucker linear regression. Current works in tensor logistic regression or bensor-based GLMs do not provide any theoretical guarantees for sample complexity (upper or lower bounds).
2019; Zhang al., 2020; Raskutti et al., 2019;
In terms of the LSR model, its motivational roots are two fold. First, a special case of the LSR model was introduced by Tsilikkarids & Hero (2013) for covariance estimation problems. An extension of this model has only recently been used tensor data for dictionary learning (Ghassemi et al., 2020). Secondly, the LSR model can be rearranged into a specialised form of the BTD model, equipped with further constraints. In terms of the GLM framework, to the best of our knowledge, our work is the first to consider the LSR (or BTD) model in regression problems.
1.3 Organization
The organization of this paper is as follows. In Section 2 we establish a background on various tensor models, as well as the LSR tensor model. In Section 3 we formulate the LSR-TGLM model and introduce two objectives: parameter estimation and minimax risk. In Section 4 we discuss the estimation problem of LSR-structured coefficient tensors in GLMs and propose an efficient algorithm for parameter estimation. In Section 5 we provide a numerical study of the LSRTR algorithm with synthetic data and experiments on real data. In Section 6 we introduce a sample complexity bound in the form of a minimax lower bound on the estimation error of the low-rank LSR-GLM model and provide formal proof in Section 6.3. We conclude our work in Section 7. Proofs of leminas for the main theorem, and additional numerical results are provided in the appendix.
2 Preliminaries
This work is based on structured tensor decompositions. For more comprehensive tutorial on tensor decompositions, see the survey of Kolda & Bader (2009). We will now list some necessary preliminaries regarding lensors and tensor structures.
We use the following notation convention throughout the paper: x, x, X and X denote scalars, vectors, matrices and tensors, respectively. Specifically, X' is the tensor defined as the aggregation of n tensors: x X X,......... X.). Given fixed tensor X. x vec(X) is the column-wise vectorization of X. The
4

---

Published in Transactions on Machine Learning Research (08/2023)
tensor In is the ... m identity tensor, such that for any tensor S. the product I -
Given n K-mode tensors (X) ) of dimension x m2 x ... x m,; X is the combined (aggregated) tensor of n samples of dimension my x m2 x x.. X mK x n. For a positive integer K, the set [ 2..,...., .K) K} so that [EK)) is the ordered set (X1, X2........ Xx). The symbol -[K] for the reverse order, so that (Xx)Re-(K] = (Xx, x-......, X1). For matrix X, the vector xl6) is the jth column of X and the vector xTO) is its 5thh row. For a vector x, x(j) is the jth element of x. If z E then |ø] is the greatest integer less than equal to We the standard notation ||x for the p-norm (p 2 1) of a vector matrix X or tensor X, the Frobenius norms are |X and FF respectively. For two vectors X1 and X2, X1 X2 denotes their outer product. Similarly, for two matrices X1 and X2, X1 @ X2 denotes their Kronecker product. The inner product between two vectors, matrices or tensors is denoted as (,, For a set of K vectors {x: X1 o............................... oxk produces a k-dimensional rank-1 tensor. For matricos (X,) K X₁ - X @ X2 8 ....... & Xx produces 'K-order Kronecker-structured matrix'. We call a matrix a matrix it is product of K 2 2 matrices. The mode-k matricization of tensor X is (1) (Kolda & Bader, 2009); and given matrix B, X X B denotes the multiplication of X by B along mode k. Finally, for all & E [K] we have X [K] B X x1 B1 X2 X Bx.
For
K-order monnecker-structuree
We now formally define GLMs for vector-structured covariates, as discussed in the literature.
Definition 1 (Vector-structureed Generalized Linear Models). Consider on observation y, vector of co- variales x E IR"", and a bias and regression coefficient vector b, both to be estimated. Let be response variable generated from a distribution in the exponential family with probability mass/density function as follows:
Ply, n) = b(y) exp(nI(()) a(n)).
(1)
Here, nj is called the natural parameter, T④y) is the sufficient statistic and a(1) is the log-portition for cumulantt function. Consider a given regression problem of estimating y given x. This requires månimising the MSE as follows:
TMMSE = arg min Fx, - y)²(( - -)]
(2)
y
The Minimum MSE (MMSE) solution the expected value of y conditioned on x, or Ely|x] = u. Now, let H be the support of y and consider strictly increasing and invertible link function g()) H + R. The Generalized Linear Model is then defined as
s(4)= = 2 1/,bb,b (b.x) 2;
(3)
where the natural porameter v) is the linear predictor, = g (m), and is the bias.
The core idea behind GLMs is that though many regression problems are not linear, the distribution of the observation y is only affected by the linear combination (b,x) z. If the distribution of falls under the exponential family (which it often does), then y is related to the linear predictor via g(-). In the case of linear regression where E IR is assumed to be a contimuuus response variable with Gaussian distribution, we have μ = (b, x) + z and g(-) is just the identity function. In logistic regression is Bernoulli distributed, taking values y E {0, 1}, and we have / = (xxxp(((bx)+2) and g(v) log 1 In Poisson regression y E N is a Poisson distributed random variable that expresses the count of event occurring in fixed time interval, while μ exp( (b, x) + z)and g(14 = log (1)).
1-
Definition 2 (CP Decomposition (Kolda & Bader, 2009; Zhou et al., 2013)). Consider a K K-mode tensor B E XX-XXXX The rank-r CANDECOMP/PARATT (CP) decomposition decomposes B inlo sum of rank-1 tensors follows:
B Σ <E(r]
b1 o... bg.
(4)
5

---

Published in Transactions on Machine Learning Research (08/2023)
where Da I e [K], i E [r] is column vector. Equivalently, (4) can be expressed in vector form
follows:
vec ()) a b
Σ b bb, iE(r)
@ b
(5)
Definition 3 (Tucker decomposition (Kolda & Bader, 2009; Li et al., 2018)). Consider K K-mode tensor B e R Xmx The
Tucker decomposition decomposes B as follows:
B G Bi X BR:
(6)
where G E xrxxx.xrx denotes the core tensor and (B, R x ]EE[K] denote the factor matrices. Equiva- lently, by defining g a veo(G), (6) can be expressed in vector form as follows:
vec B b [Bx @ @ BB.g.
(7)
Ba.
B,,
Hiu
ma x Ps
B
LSR Decompostion
Decompositio
mag
mig x
B;
Mag
(a)
(b)
Figure 1: (a): A third-order tensor under the Tucker model. Tensor B is decomposed into core tensor G and factor matrices B₁ multiplied along the 4th mode of G for k [3]. The CP model appears as a special case of the Tucker model where G is a diagonal tensor equal dimension along each mode. (b): A third-order tensor under the LSR docomposition. Tensor B is comprised of a sum Tucker-structured lensors, with core tensor G fixed across sunmands.
A visual depiction of the Tucker model is Figure 1a. The Tucker model decomposes a tensor into a core tensor G of dimension r1 X ... x Fи which is then multiplied by a factor matrix BN along cach mode E [K]. Assuming G is small (i.e., г & my Vk E [K], the factor matrices B& Vk E [K] are tall, rank Ts matrices, Additionally, the CP model appears as a specialised case of the Tucker model as follows: Fix the number of basis vectors along all modes to some R (i.e., fix the rank of all factor matrices to r) and
impose G E G X as diagonal tensor. The CP model lends a desirable compactness as the number of learnable parameters in a CP-structured tensor can be far fewer than that of an unstructured tensor. Despite this, however, the restrictive nature of the CP model (specifically its rank restriction, where each factor matrix must have equal rank) renders unfavourable against the more 'rank-flexible" Tucker model.
Our second observation refers to (7). The vectorization of fncker-structured tensor shows the Kronecker- structured matrix composed of the K factor matrices. However this 'Kronecker' structure which implies that the cocfficient vector is composed of scparable sub-matrices weighted by some vector = vee(G) is also quite restrictive, when considering different. tensor structures having the same rank. To overcome this, the Block Tensor Decomposition (BTD) was introduced and studied in recent works (De Lathauwer, 2008; Rontogiannis et al., 2021; Fu et al., 2020). alternative way of viewing the BTD structure is as a summation of S Tucker-structured tensors. The LSR decomposition we define next is special case of a BTD that uses the concept of matrix separation rank.

---

Published in Transactions on Machine Learning Research (08/2023)
Definition 4 (Matrix Separation Rank (Tsiligkaridis & Hero, 2013)). Fir m (m1, ma. ........ mk) N and r r= (r.., Fkk) E NK, and set mn = [[ME[K] my and == = Tk- Then for a matric E mmxrr its separation mank 6 m.) is the minimum number S of K-order Kronecker-struchure matrices such that
B
Σ sE(S]
((.,)) 8 ... 6 B(K,))
where B(k,s) E X
The LSR model in (8) generalizes (7) by replacing the single boonecker-structured matrix in (7) with a sum of the form in (8). It poses that a matrix can be expressed as a sum of Kronecker-structured matrices, with S being a spocialised case and what we observe in the Tucker model when the tensor vectorizod, shown in (7)- Having introduced the concept of separation rank, we now ready to define Low Separation Rank (LSR) tensors.
Definition 5 (Low Separation Rank (LSR) Tensor Decomposition). Consider K-mode tensor B E
Rm Xmx, The rank-(r1,.... .Kk) LSR decomposition with separation rank S decomposes B as follows:
B = Σ G B(,)) X2 ... X B(Kss),
(9)
SEISI
where G E RTIX denotes the core tensor and B (K,ss) E muxk k E [K], s [S] denote the Kronecker- structured factor matrices. Equivalently, by defining g A vec(G), (9) can be expressed in vector form follows:
(B(r,s) ... (...))(I.)))) g. vec B b = Σ SE(S]
(10)
A visual depiction of the LSR model Figure 1b. For the purposes of this work, we pose some constraints on the LSR decomposition. First, we assume that the KS factor matrices B(k,s) are 'tall' and rank rh (i.c., Fk & my vk [K]) and have orthonormal columns. Secondly, the LSR model corresponds to tensor B having a separation rank that is relatively small so that 1 < S < min (Ine[K)) ma: [IKE(K] rk). Thus, defining mmr as the m x Sticfel manifold, and for a fixed tensor rank (f1, rK) and separation rank S, the LSR structured tensor B belongs to the following parameter space P P(rx).S
P()))
B'
Σ G' x B((,)) x XK B(x,,) R G' E R X."XXx sE|S]
ran (r (ri...., rk), B(sss) E OmaXn, Ke [K] SE [S]
(11)
The orthonormality assumption on the columns of B ((k,e) is common assumption made in the GLM litera- ture (Zhang & Jiang, 2016). Similar assumptions such as unit-norm columns and columns with fixed entries are also common (Zhou et al., 2013; Li et al., 2018).
We have described how the LSR structure (9) is a generalization of the Tucker (and therefore of CP) decomposition. The form in (9) also special case of BTD, or a summation of S Tucker-structured tensors, equipped with orthogonality constraints and common core tensor among all summands. The BTD can also be rearranged into specialised Ttukerrstructuree tensor with block-diagonal core tensor of dimensions Kr; x ... x Krk and factor matrices of size my x Srt. In our work, however, we compare the number of learnable parameters between tensor decompositions (CP, Tucker and LSR) of same-sized core tensor (of same rank). When comparing tensor decompositions for a fixed rank (ri,..., rK), LSR generalizes the Tucker dccomposition. We also remark that while the differences between LSR and a general BTD might appear nuanced, they are significant, particularly the critical requirement for common core tensor in the LSR structure, which allows for summing the Kronecker-structure factor matrices. This distinction

---

Published in Transactions on Machine Learning Research (08/2023)
is fundamental for retaining the LSR matrix structure, thereby offering a number of benefits, which we will elaborate in Sections 3, 4, and 6.
Moving forward, Table reviews the number of learnable parameters for the three models of CP, Tucker and LSR, for a fixed rank TK). The CP model contains the least number of parameters, especially if small. The LSR model is the most complex of the three models in that it has more parameters than the Tucker model. The working hypothesis throughout this work is that the price we are likely pay for increase in sample complexity (by using the LSR model for tensor GLM problems) is worth the gain we are likely to achieve in representation power and estimation accuracy.
CP
Tucker
LSR
Parameters
Σ k= (mar) +r
k (m=(() + (mar))
(E==(mmrr) + I P'k
Table 1: Number of learnable paramcters in the three tensor models (CP, Tucker and LSR).
3
Problem Statement
We are now ready to propose the Low Separation Rank Tensor Generalized Linear Model (LSR-TGLM). Consider a response variable with probability distribution belonging to the exponential family with pa- rameter TJ, as in (1). Consider also tensor-structured covariates X E ImmIx-Xxxx and a coefficient tensor B E R X that assumes a low-rank LSR structure as shown in (9). Given a link function g()), the
LSR-TGLM model assumes vj is given by
((u) =
Σ G X2... XK B ((K,)), X 1
x/+ +
(12)
It is important to note here that we assume that the tensor rank (r1, , Tx) and LSR rank are known for reasons we will discuss in Section 3.1. Additionally for algebraic simplicity, from this point forward
will consider the standard case where 2 0 in (3) and without loss of generality will express the LSR-TGLM model as
s(4)=
G B(,e) ' XK B BK,,)), X
(13)
[E(S)
The LSR-TGIM model in (13) can be written standard GLM by vectorizing the parameters as follows
g(n) (vec B vec(X))
Σ (B(K,e) .... B( g.x SE[S]
(14)
The sum of Kronecker-structuree matrices in (14) is due to the LSR-TGLM model in (13). We can now formally define the two goals that are the core of this work: 1) Parameter estimation for -SRM and 2) Minimax lower bound on the estimation error. One of the benefits of the requirement of common core tensor is an intuitive reasoning pertaining to the role of the core tensor in the LSR-structured tensor. If one was to view G as a 'weight lensor', the same G implies that the S groups of factor matrices must be given the same weight. We also emphasise here that the conceptual novelty our work resides not introducing the concept of L.S-structured tensors, but in our unique application of the LSR dccomposition the GLM model and our comprehensive analysis of the resulting GLM, which is non-trivial.
3.1 Parameter Estimation for LSR-
The objective in regression theory to predict an outcome y based on X, which is achieved by first estimating the model parameter. For T.S.-GGMMs, we wish to find estimate of B that best fits the model in (13).
8

---

Published in Transactions on Machine Learning Research (08/2023)
The underlying (true) B is low-rank LSR-structured tensor belonging to constraint set C:
C _ e T
XMRK B E PI S Rg. S min II m II Th , [K] R
(15)
This the set of all LSR-structured tensors with constrained tensor rank tuple (ri, 22........ "2.,,,, TK) and con- strained separation rank S. Now consider the space of tensor-structured covariates X x and space of scalar observations Y C I.. There exists a probability measure, denoted as Pxy, that allows a learning procedure to randomly draw points {X E X, E Y} from the product space X x Y. We do not know Px,,, but we have access to independently sampled observations (X, Mi H= Therefore, we find an estimate of through Maximum Likelihood Estimation (MLE), making the the objective function the negative log-likel:ooo::
a
C. B = Σ
log(b (bg()(()))
Σ (B)
(16)
To find the MLE, minimise (16) over the constraint set C:
arg min log(b(u))) + Σ ((B (( Bec
(17)
The vectorized LSR-structured tensor uses sum of S Kronecker-structured matrices. Lemma 1 from Ghassemi et al. (2020) shows that there is a 1-to-1 mapping between Kronecker-structung matrices with separation rank S and set of rank S lensors. Finding the rank of a lensor is NP-hard, (Hāstad, 1990), and thus is finding the separation rank of a Krooecker-structuree matrix. Therefore, in the context of our work on LSR-TGLM, finding the LSR rank S is NP-hard, and the problem in (17) is therefore intractable. To mitigate this issue first assume that the tensor rank (r1, "2,..,,,, rK) and separation rank S are known and we solve the following factorised problem for parameter estimation for LSR-TGIMM
8 B (
)) (18)
arg min
(bg()()))) Σ
Σ G
B X X
T(y)) -
(B
subject to ((k,s) E 0maE
The expression in (18) provides a tractable relaxation for (17), where the coefficient tensor is explicitly written in terms of the core tensor G and factor matrices ((k,.) of the low-rank LSR structure. In this work, we will study the problem in (18).
3.2 Minimax Lower Bound for LSR-TGLM
Our second goal is to derive lower bound on the minimax risk of estimating LSR-structured cocfficient tensors for the LSR-TGLM problem in (13). Minimax bounds can be useful tools in developing insight into the parameters on which an achievable error of a given problem might depend and shed light on the benefits of imposing tensor structures in regression problems. They also provide a mcans of quantifying the performance existing algorithms. Previous studies of tensor-structured GLMs (Zhou et al., 2013; Li et al., 2018) fall short of providing a sample complexity analysis. The analysis this work thus instrumental in revealing the potential benefits of imposing structure on the coefficient tensor within the GLM context. We will show that the LSR model exhibits lower sample complexity compared to the vector case, and provides an expressive representation of lensor data.
We adopt a local analysis and assume that the LSR-TGLM's underlying (true) B resides within a neighbour- hood with known radius around a fixed point. However, for a sufficiently large neighborhood, the minimax lower bounds derived in this work effectively become independent of the radius. Therefore we assume that for fixed tensor rank (r1, .2.....,K)) and scparation rank S, the underlying belongs to the set
(0
B' E Pin), (BB d),
(19)

---

Published in Transactions on Machine Learning Research (08/2023)
the ball of radius d with distance metric p 11-11p; which resides in the parameter space P(rx). defined in 2 (11). Thus Oa( C Prr... and B has energy bounded by |B F < d?. Note that we fix the reference point as the tensor of all zero-elements without loss of generality. Indeed, any neighbourhood Ba(A) around point A E )r..).s is just a translation from (o( with known distance ||A|] Additionally, we note that the point 0 also belongs to the parameter space (rr)).s- The minimax risk is defined the minimum worst-case behaviour for any estimator. Mathematically, it is expressed as follows.
e" = inf sup Ey.x
0 BBB
(20)
B (D)
Here, B denotes an estimator of B, y [V U............................... .Un]: == :,) and & is function with (()) If define o = 11-11> _ R. + R, with 6(0) 0 then the minimax risk simply the worst-case Mean Squared Error (MSE) for the best estimator, ie.,
-. = inf sup Ey, B ()
| B-B
B B
(21)
Proving a lower bound e* > Eo the minimax risk shows that any estimator must have risk lower bounded by Eo- Existing minimax bounds the parameter cstimation problem in GLMs regression models provided in the literalure (Abramovich & Grinshtein, 2016; Lee & Courtade, 2020; Raskutti et al., 2011) cannot be applied here for two reasons. Primarily, these bounds do not account for the impact of the structural assumptions we make on our model. In fact, we require bounds that accurately reflect the sample complexity of estimation algorithms for LSR-structured coefficient tensors. Secondly, the link function in GLMs also involves the analysis of the space of T.SR-structured coefficients, making this part of our analysis non-trivial and fundamentally different such works. We elaborate this point further in Section 6.3. The derived minimax bound in this section conveniently generalizes the CP and Tucker tensor structures and can be specialised to existing bounds in the literature such as that for Tucker-structured linear regression (Zhang et al., 2020).
We now introduce some standard lemmas and assumptions used in this work. Assumption 1 (Covariate Distribution). For X fpm1--xma define vec A and m II mk- kE(K]
Then, x VO,,E,) and thus E(x] = 0 and B[xx"]= Ez.
Lemma 1 ((McCullagh & Nelder, 2019)). Any observation generated according to a distribution from the exponential family has mean a' (n),, i.e., the first derivalive of a(n), and variance a" dn)), i.e., the second derivative of a(n).
Assumption 2. The first derivative of the cumalant function, a(r), with respect to n is bounded unáformly by a constant M 2 0: d'(1)) S M.
Lemma 2. Consider the slandard GLM problem in Dessnition with negalive log-likelihood function
log(b():)) + Σ ((b, x;) TT(::) - a((b)x:)))
(22)
The gradient of (22) with respect to is Ei= (T(n:) g (((b, x,,
The proof of Lemma 2 follows the same steps explained in McCullagh & Nelder (2019). Assumption 1 states that vec(X) 3 x is a zero-mean Gaussian random variable with covariance matrix Σ.. We do not place any further assumptions on Σ,. Assumption 2 implies that the mean of any GLM observation is bounded. This is a common assumption made in the literature (Lee & Courtade, 2020). For binary logistic regression, e.g., this assumption is satisfied with M 1. For linear or Poisson regression, this assumption implics that the energy of y is boundd.. Though the range of y is and N for linear and Poisson regression, respectively, any outcome can be bounded by positive constant M with high probability in most instances.
10

---

Published in Transactions on Machine Learning Research (08/2023)
4 Estimation Problem and Algorithm
To solve (18) we propose approach similar to those found in prior works (Li et al., 2018; Zhou et al., 2013; Zhang & Jiang, 2016; Tan et al., 2013). The main idea behind this approach is to recognise that although the objective function in (18) is non-convex in all elements G and (B(r.s) (ee.s))/[(KK]eee [e([ jointly, it is convex with respect to each element separately. In this work, we propose a Block Coordinate Descent (BCD) algorithm that estimates each element of the LSR structured lensor B alone while holding all other elements constant. We note for the reader that BCD is simply Alternating Minimisation (AM) when there are only two factors to be estimated. Additionally, one may element; however, choose to solve it via gradient steps in this work. First, for every k' E [K] and s' E [S], we estimate the factor matrix B while keeping all other factor matrices (B(r,s) (B( (e( and core tensor G fixed. Secondly, we estimate the core tensor G while keeping all factor matrices (B, [ fixed. We also point out that though our BCD approach is similar to some existing works, in this work we also show how each sub-problem reduces to a smaller scale GLM problem with lewer learnable paramelers. We also evaluate the computational complexity of each sub-problem. We derive the optimisation problem for each step of the algorithm below.
any solver to solve each convex sub-problem to estimate each
4.1 Estimating Factor Matrices
To estimate the factor matrices we perform gradient descent over our objective function and project each iterate onto the Stiefel manifold, with the objective function corresponding to the block B(wyy)) given by:
Cn (Bi)))))
(23)
Bir Xias
B(k)dd)
Gimy)
+
Σ <E(S))s
GG [[K] B(kss): X.)
T(y:)
-a (BBi.o)) Xiee
Bik
GT
GG [ B(k.,), Bix Xi)
KKTUK'
In order to derive the expression in (23) for every Buyrarr for k [K] and s' [S], notice that:
B.X G 11 (1,ss) '.... x B B(K.)). X
Σ
G B (1.)) X2 ... X B(K,)), X)
#E[S]
(Basss) (BiR,)) @ B
B - 0 8... B (l.,.)) G(m))
(24)
Σ (Gx X XK B (K.)),
Now define the following notations in order to keep the expressions in (24) concise. From the first summand in (24) we define
(w) vec (X(N) (B(r,et) ... B(++,, @ ((v- .... (....BB))))
Ba,,))) GiTT)):
(25)
and from the second summand in (24) we define
(K())) = Σ
(Gx,Brl. X2 ... X B(K,)).)
(26)
We can then express (24) as
B. X (vec (B (())))))))))))))) +00 ([vec(B (
(27)
11

---

Published in Transactions on Machine Learning Research (08/2023)
By defining x a ( (ar,)) we rewrite (23) as
arg min
Σ [[ 1), x:)) T(v:) ([
(28)
B aty 60m0 Xee
The most important insight here is that the resulting problem can be viewed as a parameter estimation prob- lem for GLMs with B( as the 'paramcter' and x; as the 'predictor' or structured covariates. Estimating s( alone in particular results a low-dimensional problem with myTK parameters. We solve (23) via projected gradient descent. Using the above defined notations, the gradient of Lu (())))))) with respect ((,,)) is
8cn Oveo (3)
(T(w)) /- (((ve(B(x.s)))1)).)) &,,
(29)
and the projection operator H. IR XTW of orthogonal matrices is defined as
Xmw that projects the obtained iterate onto the manifold
HB arg min Beamaexyy
B B-BBrra
(us .z) |
(30)
The solution to the problem in (30) is simply obtained using the QR decomposition of B(ass))
4.2 Estimating the Core Tensor
For core tensor G, we are only required to perform gradient descent. The optimisation problem is
arg min La (G). where
(31)
G
)
Ln G
g,
Σ
B
(k,,)
xi)
T(w:)
g. sE
Σ
B
(k.o)
x
i=1
<E(s]
[
[A]
In order to see the derivative of the expression in (31), notice that,
T
B. X
(Vec(B), vee (X)
Σ
B
((saa)
g.x
g: Σ
(
(32)
+E[s]
[K]
eE(S]
[K]
With a slight overload of notation we also define the following:
x
Σ
(B(K,s) 8 ... @ Baa,)) x,
(33)
SE[S]
and further express (32) and (31)
[B (g.
(34)
and
arg min Σ (g. a x,)),
(35)
respectively. Once again, the resulting problem can be viewed as parameter estimation problem for GLMs with G as the 'parameter' and &, as the 'predictor' or structured covariates. Estimating G alone results in a low-dimensional problem with II Fa parameters. We solve (31) via gradient descent, where the gradient k
of (G) with respoct G
OCn Övec(G)
= ΣΣ (T(s)) --- ((g, x,
(36)
1=1
12

---

Published in Transactions on Machine Learning Research (08/2023)
4.3 Final Algorithm: LSRTR
We summarise the procedure discussed above in Algorithm 1 and name algorithı Low Separation Rank Tensor Regression (LSRTR). We also show the prediction procedure performed using the estimated cocfficient tensor in Algorithm 2, where we calculate the mean E|y[X] of the postcrior distribution basod on the estimated coefficient tensor B from Algorithm 1. The posterior mean allows us to make predictions for an observation y and report confidence probabilities. Note that the convergence of BCD on non-convex problems is a function of the initialisation. In this work, we initialise LSRTR randomly on the constraint set. That is: an LSR-structured tensor with (KS) orthogonal factor matrices, randomly generated on the Stiefel manifold, and core tensor with random Gaussian entries. We also note that though BCD algorithms such as LSRTR are popular amongst tensor-structured regression works and have been shown to be effective in practice (Zhou et al., 2013; Tan et al., 2013; Li et al., 2018; Zhang & Jiang, 2016), BCD algorithms in GLM works do not explicitly exploit the coefficient tensor's Kronecker structure that appears upon vectorization. Keeping common core tensor and maintaining the LSR matrix structure (14) allows to explore other parameter estimation algorithms that exploit the Kronecker matrix structure, similar to those in existing dictionary learning works (Ghassemi al., 2020). However, leave this for future work.
prior
We next provide brief discussion on the per-iteration computational complexity of each sub-problem of LSRTR. For factor matrix B(k,)), each gradient step has complexity of order ((mrrk), and the QR projection step has complexity of order ((mmrk))) Thus the per-iteration computational complexity of estimating a factor matrix is ((m)r))))) Since there is no projection when estimating the core tensor G, its per-itcration sample complexity is just O rk).
Algorithm LSRTR: A block coordinate descent algorithm for LSR-TGLMs
1: Input: training samples (X,, 5: Mi h step size a, separation rank S, tensor rank (ri, (2k..,,K))
2: Initialise: Factor matrices B A l E [K], s E [S], core tensor Go and t t- 0.
(k.s)
3: repeat:
4: for s' E [S]] do
5: for k' E [K] do
(1)
( (kt,a))
Σn (T(y:) - g- (((ee)BB().)))) ,x,)), x;
(K) vec
7:
B ((+1)
(K's'))
H
(t))
(K'
end for
end for
10:
G
[1+1)
t (G)) Ei=l (T(y:) -g -(g
x:)))
11:
t +-t+1
12:
until convergence
return B Evelst (())
X1 (() (1,8) X2 X B (K,
13:
Algorithm 2 Posterior prediction for LSR-TGLMs
Input Estimate B E mi xxxxxx and Nhe test data points (X, 1 Output Expectation ğa = E[yn 1,), Elya X Elyn,. X m
1: Define: X' X.,.,, X
Compute fü for input as: ça = g ( return ī
In this work, we do not provide convergence guarantees for LSRTR; this is non-trivial task that we defer to future work. The main challenge in proving convergence is the constraint set in LSRTR, which is a prod- nct space of Stiefel manifolds and is neither closed nor convex. The non-convexity of the constraint space makes a convergence analysis for LSRIR difficult since we cannot apply general results for the convergence of BCD (Bertsekas, 2016). Therefore, one at best expect L.RTR to converge to local minimum. The
13

---

Published in Transactions on Machine Learning Research (08/2023)
non-convexity also presents challenges for analysing projections, in particular when using the QR decompo- sition. More specifically, projection operators in projected gradient methods must be non-expansive in order to prevent potential amplification of the estimation error. Since the general method of proving the non- expansiveness projection also assumes a closed and convex constraint set - which the Stiefel manifold is not (Bertsekas, 2009) demonstrating the non-expansiveness of QR projection the general sense also non-trivial. In many works, however, the QR decompasition a common tool for projecting onto and/or optimising over the Stiefel manifold (Absil et al., 2008), yet proving stronger guarantees on optimisation over Stiefel manifold is larger and non-trivial problem. Therefore, we leave the theoretical analysis of Algorithm 1 as significant direction for future works. However, we make some promising remarks on its practical performance, convergence, and projection behaviour. Specifically in Section we assess the perfor- mance of our proposed algorithm on various GLM problems for 2-dimensional and 3-dimensional synthetic data. Moreover, in Section 5.2, we use the LSR-TGLM model and the proposed LSRTR algorithm to solve classification problems on several real-world medical imaging datasots.
5 Numerical Study and Experiments on Medical Imaging Data
In this section we provide a comprehensive numerical study in order to assess the efficacy of the L.RR- TGLM regression model and the corresponding proposed LSRTR algorithm. The objectives of this study are three-fold. First, we investigate the performance of our algorithm, particularly: 1) The performance of our proposed approach (parameter estimation and prediction) against substantial growth in number of parameters, and 2) The performance of our algorithm, particularly when faced with a large increase in sample size. Secondly, we compare the estimation and prediction accuracy gains of the LSR model to current regression models and algorithms in the literature on synthetic data (Zhou et al., 2013; Li et al., 2018; Zhang et al., 2020; McCullagh Nelder, 2019; Sober & Lee, 2003). Thirdly, we assess the performance of the LSR-TGLM model on classification problems for medical imaging datasets.
5.1 Experiments on Synthetic Data
We begin with experiments on synthetic data for three GLM problems: linear regression, logistic regression and Poisson regression. The purpose of the experiments is to answer the following questions:
1. When the underlying coefficient tensor has LSR structure, does the proposed LSRTR algorithm lead to reliable parameter estimation? Specifically, does our algorithm achieve a reduction in sample complexity and estimation error compared to pector-based methods (for a fixed sample size)?
Does our proposed algorithm: 1) provide reliable parameter estimation and prediction against a substantial growth in number of parameters and 2) have per-iteration computation time comparable to other tensor-based methods, when faced with a large increase in sample size?
Empiriaal convergence behaviour: does the proposed algorithm converge to a stationary point?
In our experiments we compute and report the following:
1. A coefficient tensor B estimated through the following learning methods: (i) The LSRTR method (Algorithm 1), (iz) Low-rank Tucker model methods, specifically, bllokk relaxation' algorithm for parameter estimation of Tucker-structured GLMs proposed by Li et al. (2018) (we will call this procedure TTR), and similar iterative procedure for logistic Tucker regression (LTuR) with Frobenius norm regularization proposed by Zhang & Jiang (2016), (iii) unstructured (vector-based) methods for regression. In the case of (iiz) we use Least Squares (LS) for linear regression, a first- order method we name LR (Seber Lee, 2003) for logistic regression, and GLM fitting algorithm we name PR, with log link function from McCullagh & Nelder (2019) for Poisson regression. The
performance of each method is evaluated via the normalised estimation error defined as
B
14

---

Published in Transactions on Machine Learning Research (08/2023)
2. A predicted response vector y. The prediction accuracy for linear, logistic and Poisson regression is evaluated via the normalised squared error defined as - the Mean Absolute Error (MAE) layll
defined as اا>او| for The testing samples, and the normalised squared logarithmic error defined as
respectively.
[[og
3. The magnitude of the computed gradient of the loss function, denoted as [VL(())l at the final iteration of every sub-problem in LSRTR.
Since we are reporting normalised estimat ion error, we only consider errors below 1. Indeed, the normalised error compares the performance of an estimator relative to the 'trivial estimator' that always outputs 0. A normalised error of 1 occurs when the estimate is O, yet an algorithm can potentially perform worse than the trivial estimator, allowing an error greater than 1. Such an error is insubstantial and we disregard when analysing the performance of models and algorithms.
5.1.1 Experimental Setup
LSR-structured coelficient lensor B E R X- with separation rank S in (9) as follows. The entries of the core tensor G are sampled from a Gaussian(O, pizg..rx) distribution. Each factor matrix B(s..) is constructed by first constructing matrix A(k,s) E mmxx with independent standard normal entries. The QR. decomposition is then applied to A(kss) in order to obtain the orthonormal matrix A(k,s) E 0 The first тІ columns are then ex- tracted from ((tas) in order to obtain B ((k,s)- With the above construction, we have BI G (we derive more precise bound on [BI in Section 6.3). We also generate n data samples {X, y} used for training and estimation, and Nite samples used for testing and prediction. We generate independent tensor-st ructurre covariates (X,)=1 and (X, with zero-mean Gaussian entries and covariance Σ,, and obsorvations {y and {yi (yr) which are randomly generated according to the probabilistic model in (13) and with appropriate link function g(u). Additionally, is assumed to follow standard Normal distribution. Thus, for linear regression, the observation yi conditioned on data X, follows Gaussian distribution, i.e., yi ~ W (B..,)) ,), 1). In logistic regression, the observation Yi conditioned on data X, follows a Bernoulli
In our experiments we generate an (r,, " Fxr.a.
distribution, i.e., y/ ~ Bernoulli [+exp(-(B.x.)) (B.,,, Finally, in Poisson regression, the observation con-
ditioned on data X; follows a Poisson distribution, i.e., Yi ~ Poisson (exp( (B, X)))). The experiment for each GLM problem is performed for various model sizes (m=)<e(k)) tensor ranks (r*)Re(K) and scparation rank S. The coefficient lensor B and data samples 1/ and 1x,, are generated as described above and each experiment is repeated for incrcasing value of n. We then compute the estimate B and Finally, for each GLM problem, unless otherwise stated, we conduct 50 repetitions of each experiment (for fixed coefficient tensor B) and average the error over the repetitions.
prediction ÿ [01. 22................. ... In... 1i, ... Pn... where j. - ((B.)((
5.1.2 2D Synthetic Data
We begin with the two-dimensional (matrix) case with K = as many medical imaging data fall under this case (x-rays, EEG, fiber-bundle images.). In this set of experiments the underlying cocfficient matrix B has dimensions m m and rank r. The dimensions and rank of B are thus represented as the tuple rm, For all GLM problems we fix the separation rank $ = thus the number of learnable parameters is S(mr + mr) + p2 4mr + r2. We now proceed with the specific experimental setups of the three GLM problems under study.
For linear regression we consider various model sizes and ranks, i.e., m E (64, 128, 256), r E {4,8). We range sample size n from (roughly) the degrees of freedom (number of Icarnable parameters) of the smallest model size to no smaller than the degrees of freedom of the largest model size. The smallest model size is for the tuple (m- 64, 4) and the largest model size is for the tuple (m 256, 8), which, under the LSR
15

---

Published in Transactions on Machine Learning Research (08/2023)
structure, have 4(64 4) 4 = 1040 and 4(256 x 8) + 8 = 8256 learnable parameters, respectively. We compare the performance of LSRTR to TTR and T.S. Figure 2 reports the mean normalised estimation and prediction accuracy, respectively, across 50 repetitions. The shaded regions are one standard deviation of mcan cstimation and prediction accuracy, based on 50 replications.
For logistic and Poisson regression we consider m E {32, 64, 128) and r E (4,8). Since the models are non- linear, we follow the heuristic rule of ranging sample size n from (roughly) 5 times the degrees of freedom of the smallest model size to smaller than that of the largest model size (Li et al., 2018). The smallest model size is for the tuple (m = 32, = 4) and the largest model size is for the tuple (m = 128, r= 8), which, under the LSR structure, have 4(32 x 4) + = 528 and 4(128 x 8) + 8 = 4160 Iearnable parameters, respectively. We compare the performance of the LSRTR method to i) LTuR and LR and i) TTR and PR for logistic and Poisson regression, respectively. Figure and Figure 10 in the appendix report the estimation and prediction accuracy for logistic and Poisson regression, respectively, with shaded regions depicting one standard deviation of mcan cstimation and prediction accuracy, based on 50 replications.
We pause to make a few remarks. Firstly, in the case of synthetic data, the rank r and separation rank S are assumed to be known and other algorithmic parameters (the step size a) are set using separate validation experiments. Secondly, though the generated coefficient B is matrix, the numerical study by Li et al. (2018), shows that the low-rank Tucker model effectively estimates two-dimensional parameters and the TTR algorithm performs well in the matrix setting even with huge increases in model and sample sizes. The study by Li et al. (2018) also shows how TTR generalizes CP tensor regression algorithms and improves upon the performance of several other methods (such as PCA and Bayesian regression methods) in terms of estimation and prediction accuracy in experiments with synthetic data and for several regression types (linear and logistic) (Li et al., 2018). For these reasons we limit our scope of comparative algorithms for matrix-structured regression to TTR.
We also note that since we report normalised errors, we do not investigate any errors greater than 1 in Figures 2, 3 and 10. The plots in Figures 2, 3 and 10 show that, for all algorithms in all GLM problems, the estimation and prediction accuracy improves and the shaded regions (which characterise the slandard deviation of the error points over 50 repetitions) decay with an increase in observations. In particular, on the whole, the tensor-basod methods LSRTR, TTR and LTuR outperform the vestor-based methods LR, LS and PR, particularly with larger model size m. This is because the number of observations provided in the experiments is typically much lower than the sample complexity of the vector methods. For example, consider linear regression with m = 256. LS requires at least 256 x 256 = 65,536 observations for reliable parameter estimation, while if we set r = and S = 2, TTR and LSRTR require only 2064 and 4112 observations for reliable parameter estimation, respectively. Additionally, for TTR, LTuR and LSRTR, as the core tensor dimensions grow (i.e., as we increase the Tucker rank of our tensor models), we require a relatively larger sample size to achieve better accuracy. This is not surprising as a larger core tensor increases the number of parameters to be learned. Figures 2, 3 and 10 also show that for all GLM problems under study (namely lincar, logistic and Poisson regression), the LSRTR algorithm shows 'consistent' performance for estimation and prediction. What we mean by this is that, even with a substantial growth in model size, the estimation error induced by LSRTR decreases with increasing sample size, and the LSRTR algorithm can still achieve acceptable estimation and prediction errors as long as the number of observations provided is large enough. Additionally, we can see that LSRTR algorithm can avoid over-fitting in the prediction step even with very large number of observations.
We also gain some additional insights from the results. First, if the underlying parameter model B is LSR-structured, then LSRTR can recover B with a relatively low number of observations, in the sense that the number of observations required is proportional to the intrinsic degrees of freedom in the LSR model (rather than the extrinsic dimensionality of B). This implies that though we do not provide any theoretical guarantees for Algorithm 1, it can be employed successfully in practice. Second and more interestingly, if the underlying parameter model B is sttuttuu then imposing the Tucker model on B in the parameter estimation procedure- as in TTR and LTuR algorithms results in less accurate estimation compared to LSRTR. Particularly, we see that for relatively large model sizes, such as in Figures 2e and 2f, or for non- linear models, such as in Figures 3e, and 10c, TTR and LTuR algorithms' performance plateaus, even when the number of observations has far exceeded the degrees of freedom of a structred coefficient matrix
16

---

Published in Transactions on Machine Learning Research (08/2023)
LS: m 64
LS: m 64
TTR: (m,
TTR: (m, (64,4)
Nomalised Estimalion Error
TTR: (m, = (64,8)
(m.r) (64.4)
Normalised Squared Emor
2.5
LSRTR: (m.r) (64.4)
TTR: (mn,7 = (64.8)
LSRTR: (m..r) (64.8)
LSRTR: (m.rr) (64.8)
1.5
1.5
0.5%
1000
Number of observations
Number observations (n)
(a)
(b)
10
10
TTR: (m,r) (128,4)
TTR: (m,r) (128,
TTR: (m.r) (128.8)
TTR: (m. (128,
10%
LSRTA: (m..r (128.
Normalised Squared Error
10
LSRTR: (m.r) (128, 4)
Nomallsed Estimation
LSRTR: (m..r) (128, 8)
LSITT:: (m.rr (128, 8)
10
10
n0-2
10
1000
5000
8000
1000
3000
8000
Number of observations (n)
Number of observations (n)
(c)
(d)
10
10
LS:
L.S 256
TTR: rm, (256,
TTR: (m, (256. 4)
TTR: (m.r) (256,
(m.r) (256.
10%
LSRTR: (m.rr) (256, 4)
10%
LSRTR: (m,r) (256, 4)
Normalised Estimation
LSRTR: (m..r) (256, 8)
LSRTR: (m.r) (256.
Normalised Squared
10
1012
1000
2000
3000
4000
5000
6000
7000
8000
8000
1000
2000
3000
4000
5000
6000
7000
9000
Number observations (n)
Number observations (n)
(e)
Figure 2: Comparison of LSRTR with L.S. and TTR for two-ddmennsinall synthetic data when (64, 128, 256), {4,8) and 5 2. Normalised estimation error for 64, 128, and 256 is shown in (a), (c), and (e), respectively. Normalised prediction error for 64, 128, and 256 is shown in (b), (d), and (f), respectively. Each marker represents the mean normalised estimation iiimaation//pedictioo preediction errors, over 50 repetitions. The shaded regions correspond to one standard deviation the mean normalised errors.
17

---

Published in Transactions on Machine Learning Research (08/2023)
LR:
LR: ma
LTuR: (m.r) (32.
0.9
LTait: (m, (32,
Nomalised Estimation Error
LSRTR: (m.rr (32,4)
LSRTR: (m.r) (32.4)
0.7
Mean Absolute Error
0.7
0.6
0.8
0.5
0.5
0.4
0.4
0.1
0.1
Number of observations
Number observations (n)
(a)
(b)
1.2
LTuR: (m, r) = (64.4)
0.9
LTui:: (m, = (64, 4)
Nomallsad Eslimation Errox
LSRTR: (m.r) (64.4)
LSRTTR: (m.r) = (64.4)
Meen Absolute Error
0.8
0.8
0.6
0.5
0.4
0.4
0.3
0.2
0.2
0.1
7500
15000
17500
Number of observations (n)
Number of observations (n)
(c)
(d)
1.8
LR: 128
LR: 128
LTaR: (m, = (128, 4)
0.9
LTuR: (m, = (128, 4)
Nomallised Estimaaiion Error
1.4
LSRTR: (m..r) (128, 4)
0.8
(m.r) (128,
Meen Absolute Ergr
0.7
1.2
0.5
0.8
0.4
0.3
0.2
0.2
2500
5000
7500
10000
12500
15000
17500
2500
5000
7500
10060
12500
15000
17500
Number observations (n)
Number observations (n)
(e)
(f)
Figure 3: Comparison of LSRTR with LR and LTuR for two-dimennional synthetic data when E {32, 64, 128), 4 and 5 = 2. Normalised estimation error for 32, 64, and 128 is shown in (a), (c), and (e), respectively. Normalised prediction error for 32, 64, and 128 is shown in (b), (d), and (f), respectively. Each marker represents the mean normalised estimation repetitions. The shaded regions correspond to one standard deviation the mean normalised errors.
crimaationn/redicttoonn ppredictionn errors, over 50
18

---

Published in Transactions on Machine Learning Research (08/2023)
B. In fact, even LS and LR outperform TTR and LTuR when the model size is small enough and enough observations have been provided. Though we expect LSRTR to perform better than TTR or L.Tu when B is originally LSR-structured, the results suggest that if regression problems with real data have model parameters that are approximately-LSR structured, then TTR or LTuR may not be suitable algorithms and we can motivate the use of algorithms such as LSRTR that allow for the consideration of richer class of tensors that may lead to more reliable estimation.
Notice also that in Figures 2e and 2 (as well as in Figures 3e and 3f), when = 8, the normalised errors are greater than 1. Here, LSRTR is operating in the under-sampling regime, where the model size is very large and the number of available samples (9000) is just about the number of Icarnable paramcters of the model. In such a setting we expect the estimation error to be large and to decrease with increasing samples. We also make specific remarks on Figures 3c and 3e, In the case of m = 64, we see that LSRTR shows lower estimation error at n = 17500 (0.11 vs 0.29), yet has a larger standard deviation (shaded region) over 50 experiment replications (0.04 vs 0.006). This is an acceptable standard deviation, especially considering the reduction in estimation error compared to the Tucker-based LTuR algorithm. In the case of m 128, LTuR shows lower estimation error than LSRTR until about be duc to the fact that LSRTR is operating in the undersampling regime for n = 7500 (as LSR has more learnable parameters than Tucker). Additionally, we see in Figure 3e that the estimation error for LTuR increases after n = 7500, and in Figure 3f that LSRTR outperforms LTuR, which suggests that LTuR may be overfitting the GLM model.
7500 (errors of 0.6 vs 0.4). This could
We also make a specific remark on Figure 10a in the appendix, where we witness a slight increase in estimation error for rm, = (32,8) 8) between 5000 and n 7500 (0.025 vs 0.08). However, we see that the estimation error continues to decrease after n 7500, and we do not witness this phenomena in any other experiment in Figure 10.
5.1.3 3D Synthetic Data
We have evaluated the performance of the LSRTR algorithm against an increase in model and sample size for K 2. 2. Now, we wish to further explore the performance of our method compared to the other state-of-the-art methods with 3-dimensional data (K 3). In these experiments, the underlying coefficient tensor B has dimensions mi X ma my and Tucker rank fr. r). The dimensions and rank of B are thus represented as the tuple (m,r) where m [m1, m2, ma]- We study the linear regression problem and we fix the separation rank S 1, i.c., B is Tucker structured. Thus B is 3-dimensional Tucker tensor of dimension mi m2 my and rank r along each mode, and the number of learnable parameters is S(mir + m2r + m3r) +r3 = r(m: + m2 + m3 + r²). We consider the model size m = [16, 32, 64], r = 4. The model size under the Tucker structure is 4(16 + 32 + 64 + 42) = 512. We range sample size n from (roughly) two times the degrees of freedom to roughly eight times the degrees of freedom. We compare the performance of L.RTR to TTR and L.S. Figure 4 reports the estimation accuracy and prediction accuracy, respectively, with shaded regions depicting one standard deviation of mean estimation accuracy, based on 50 replications. For LSRTR the rank r is assumed to be known and other parameters (the separation rank S and the step size o:) are set using separate validation experiments. The separation rank was thus chosen as 2.
We observe from Figure 4 the same trends as in the previous set of experiments. We also conclude that even with the underlying coefficient tensor being Tucker-structured, our proposed LSRTR method in this experi- ment cxhibits better statistical performance with fewer observations than the 'Tucker specific' algorithm TTR as well as LS. This suggests that the LSR model induces a richer class of solutions with greater representation power that is able to capture meaningful information in the data that the more compact state-of-the- tensor models do not - even with increased sample size. The results also suggest that the projection step in LSRTR (imposing orthonormal factor matrices) may contribute to the enhanced performance of LSRTR in terms of estimation opredicctooo error (as TTR does not impose such a constraint).
Finally, we numerically investigate the per-iteration computational complexity of LSRTR vs TTR for three-dimensional synthetic data, for varying sample size n. We repeat the previous experiment, i.c., [16,2, 32, 64), r 4, S 2, over varying sample size (n 100, 300, n 500), and report the
19

---

Published in Transactions on Machine Learning Research (08/2023)
10
10
IS (116, 32, 64))
LS: (116,32,64) 64
TTR: (m, r) = (116, 32, 64)
TTR: (m. r) = ([16, 32, 64,44
Nonmalised Estimation Errr
10%
LSRTR: (m.r) (11,,22,44.
Normalised Squared Error
10
LSRTR: (m.r) - (116, 64.
1
10x1
10
1000
2000
3000
1000
2000
3000
4000
number of observations (n)
Number of observations (n)
(a)
(b)
Figure 4: Comparison of LSRTR with LS and TTR for three-dimensional synthetic data when m [16, 32, 64], 4, S 1. For LSRTR, S 2 was chosen. Normalised estimation error is shown in (a). Normalised prediction error is shown in (b). Each marker represents the mean normalised estima- tion/ /on//eediction errors, over 50 repetitions. The shaded regions correspond to one standard deviation of the mean normalised errors.
n = 100
n 300
n = 500
LSRTR
0.1 (7e 5)
0.37 (2e 4)
0.65 (5e
TTR 0.048 (3e 5)
0.18 (1e 4)
0.31 (2e- 4)
Table 2: Per-iteration computation time (in seconds) of LSRTR and TTR over varying sample size (n = 100, n = 300, n n= 500). We report the mean (and variance) over 50 repetitions.
mean per-iteration computation time (in seconds) over 50 repetitions of the experiment, for LSRTR and TTR. Table shows the results. We can see that LSRTR and TTR have comparable per-iteration compu- fation times. In fact LSRTR has computation time of roughly = times that of TTR.
5.1.4 Convergence Behaviour of LSRTR
We have previously discussed the challenges of providing theoretical guarantees for the proposed algorithm, specifically an analysis of its convergence and the non-expansiveness of the QR projection. Nonetheless, we are still interested in developing an understanding of the practical behaviour of LSRTR. For two-dimensional synthetic data, for fixed sample size n, we repeat the experiment from Section 5.1.2, i..., for linear regression of model sizes (64, 4), (128, 4) and (256, 4), with n 7000, and logistic regression of model sizes (32, 4), (64, 4) and (128, 4), with n 7500, and with S error for LSRTR and the magnitude of the computed gradient at the final iteration of every sub-problem in Algorithm 1. Figure 5 depicts the convergence behaviour of LSRTR. We see that the per-iteration normalised estimation error decreases with increasing number of iterations. In regards to the non-expansiveness of the QR projection, these results suggest that the error does not amplify, and the QR step does not prohibit LSRTR from finding a 'good estimate" of the underlying coefficient tensor. The decrease in the estimation error also suggests that LSRTR converges. Morcover, Table depicts the magnitude of the computed gradient of the loss function the final iteration of every sub-problem in LSRTR. We can see in all cases that the magnitudes of the gradients are very close to 0, suggesting the convergence of LSRTR to stationary point.
We report the per-iteration normalised estimation
20

---

Published in Transactions on Machine Learning Research (08/2023)
Normalised Estimation Error (Linear Regression)
2
Normalsed Estimation Error (Logistic Regression)
LSKIR: (m,r) (64,4)
LSKUR: (m,r)
LSRTH: (m,r) (128.4)
1.8
LSKTR: (m, (61.4)
(m,r) (256.4)
LSRTR: (m.r)
1.8
1.6
1.4
1.2
1.2
0.8
0.8
0.6
0.6
0.2
0.2
Number of iterations
Number of iterations
(a)
(b)
Figure 5: Convergence behaviour of LSRTR for two-dimensional synthetic data for linear regression when m E {64, 128, 256], and r 4, for n = 7000, and logistic regression when m E {32, 64, 128}, and r 4, for 7500. For LSRTR, is shown in (a) and (b), respectively.
was chosen. Normalised estimation error for linear and logistic regression
Linear Regression (m,r)
Logistic Regression (m,
gradient norm
(64, 4)
(128, 4)
(256, 4)
(32, 4)
(64, 4)
(128, 4)
(v))(l.)))ll
3.1 x 10-3
5.2 x 10-3
4.3 x 10-3
0.24
0.38
0.21
(VL(B(..))
7.8 x 10-3
1.6 x 10-3
9.3 x 10-3
0.2
0.32
0.3
VL
1.7 x 10-2
4.7 x 10-3
8.73 x 10-4
0.4
0.1
0.18
VL
3.5 x 10-2
2.1 x 10
4 x 10-3
0.29
0.15
0.23
[V l2
7.42 x 10-4 3.8 x 10-3
1.2 x 10-3
0.37
0.28
0.19
Table 3: Magnitude of the final iteration gradient, per sub-problem, lincar and logistic regression.
5.2 Experiments on Medical Imaging Data
We move on to investigating our approach medical imaging data. Regression analysis of medical imaging data can be a useful tool in medical decision making. We chose this type of data for several reasons. Medical images are usually multi-dimensional,, and since data acquisition in medical sciences is expensive, the regression problems under study with such data is very high dimensional. Another major reason is that medical images model complex biological structures (such as brain maps). One expects the model parameter the corresponding regression model model this complexity. Lastly, the of observations in medical imaging data is usually severely imbalanced. In other words, positive medical diagnoses are relatively rare, and thus given set of observations will contain far fewer positive diagnoses than negative diagnoses. The high-dimensionallty complexity and imbalanced obscrvation set of medical imaging data serve as an excellent assessment of the elliciency of our proposed model and algorithm. That is, we contend that due to the LSR model's favourable bias-variance trade-off and augmented representational power, LSR-TGLM may be well-suited handling such data and may perform more robustly than some other compact tensor- structured GLMs. Here we study three datasets; ABIDE Autism (Craddock al., 2013; Lodhi & Bajwa, 2020), ADHD200 (Bellec et al., 2017) and Vessel MNIST 3D (Yang et al., 2021b).
21

---

Published in Transactions on Machine Learning Research (08/2023)
5.2.1 ABIDE Autism
The Autism Brain Imaging Data Exchange (ABIDE) dataset contains the resting state fMRI data collected from 98 subjects. The data has already been preprocessed for motion realignment and correction, slice timing correction and image normalisation. Each data sample corresponds to 111 cortical and sub-cortical brain regions scanned over 116 time periods (Craddock al., 2013). Therefore each sample is 111 x 116 matrix of fMRI data from one subject. Each observation is binary response variable y E {0, 1} depicting subject's 's diagnosis of either having autism (y 1) not (u = 0). The goal is to classify the subjects as either autistic or not. We perform a single train-test split procedure of 80 and 14 training and testing samples respectively. As for the autistic-contro) split, we choose 40 autistic and 40 control samples uniformly at random for training, and 14 test subjects the same way. The case-control ratio is this dataset is thus 1 11. ABIDE Autism is the only balanced datasot this set of experiments.
5.2.2 ADHD200
ADHD200 is repository of resting state fMRI images of subjects from 8 research sites: Peking University, Brown University, Kennedy Krieger Institute, Donders Institute, NYU Child Study Center, Oregon Health and Science Universily, University of Pittsburgh and Washington University in St. Louis. The data includes brain maps of fractional Amplitude of Low-Frequency Fluctuations (fALFF) of 959 child subjects. fALFF brain maps can be useful in predicting Attention Deficit Hyperactivity Disorder (ADHD) in children. The data has been preprocessed, (ADHD-200 200 Consortium, 2012; Bellec al., 2017), and we also perform standard preprocessing of the data such as removing missing observations images with poor quality (a list. of images with poor quality has been made public, (ADHD-200 Consortium, 2012)). Each data sample corresponds to 3-dimensional patient IMRI brain scan (T1-weighted image) of size 121 x 145 x 121. Figure 6 shows the axial, sagittal and coronal views (views across each tensor mode) of a fMRI data sample. Each observation is a binary response variable y E (0. } depicting a child subject's diagnosis of either having ADHD (y = 1) or being typically developing (y = 0). The objective here is classification of children subjects with ADHD. The dataset has already undergone train-test split. procedure. The training dataset consists of 762 samples, 280 of which are labelled as ADHD crrperactive/mmpllsii labelled as typical (control). The testing dataset consists of 197 samples, where 76 and 93 are labelled ADHD and control, respectively. The data is slightly imbalanced where the case-control ratio in the training and testing sets is 3 : 5 and 4: 5, respectively.
inattentive and combined) and 482 of which are
(a)
(b)
(c)
Figure 6: Axial (a), coronal (b), and sagittal (c) views of an (MRT data sample.
5.2.3 Vessel MNIST 3D
McdMNIST (Yang et al., 2021b) is a collection of medical imaging datasets suitable for various tasks such as classification and regression. One of the datasets the Vessel 3D dataset (Yang et al., 2020) containing 1909 reconstructions of vessel segments from brain Magnetic Resonance Angiography (MRA) data of different
22

---

Published in Transactions on Machine Learning Research (08/2023)
Dataset
Covariate Dimensions
Train Size
Test Size
ABIDE Autism
(111, 116)
80
14
ADHD200
(121, 145, 121)
762
169
Vessel MNIST 3D
(28, 28, 28)
1335
382
Table 4: Description of three medical imaging datascts (ABIDE ADHD, ADHD 200, and Vessel MNIST 3D). Described are the model dimensions and train and lest sizes for all three datasets.
subjects. Each data point corresponds to 3-dimensional vessel segment of dimension 28 x 28 x 28. Each vessel segment may or may not exhibit an intracranial ancurysm and the datasct is diverse in the sense that it includes a variely of shapes and scales of intracranial aneurysms. The data has also been preprocessed to restore incomplete scans made by neurosurgenss and to remove duplicated data. Each observation is a binary response variable y E {0, 1} depicting a vessel segment diagnosis of cither having an ancurysm (y = 1) or being healthy (y = 0). The goal is to therefore classify each vessel segment as either exhibiting aneurysm or being healthy. The dataset has already undergone train-test split procedure. The training dataset consists of 1335 samples, 150 of which are labelled as 'aneurysm' and 1185 of which are labelled as 'healthy'. The testing dataset consists of 382 samples, where 43 and 339 are labelled as 'ancurysm' and healthy', respectively. The data is severely imbalanced where the case-control ratio in both the training and testing sets is 3 25.
5.2.4 Results of Experiments on Medical Imaging Data
All three medical imaging datasets support logistic regression problem. We estimate the unknown coefficient tensor B through our LSRTR algorithm and use our estimate to predict the responses of the test data through a logistic regression model, i.c., we calculate the posterior probability using the sigmoid function = (1/(1 + exp( (B. X,) + +))) for every test subject i and use a cut-off threshold of 0.5 to classify each response as or O, where the positive class is chosen if s > 0.5. Parameters including the Tucker rank (чa)иe(K). the scparation rank S and the step size a are set for each dataset using separate validation experiments. On all three datasets we compare our method with four other methods, two of which have been previously introduced, namely, LR and LTuR. The third method is low-rank CP model method (which we name LCPR), specifically a block coordinate descent approach for parameter estimation of CP structured logistic regression models proposed by Tan et al. (2013). The fourth method is Support Vector Machine (SVM) with Gaussian Kernel (Hearst et al., 1998). In the case of the Vessel MNIST dataset, we also compare our mothod with an established bascline in the literature for Medical MNIST datasets: ResNet 50, which is residual neural network that is 50 layers deep. This neural network has also been augmented with a 3D convolutional layer designed to handle 3-dimensional data.
To evaluate the performance of each method we report the following scores: Specificity defined as the true negative rate, sensitivity defined as the true positive rate, the F1 score, the average accuracy score defined as the MAE, and the Area Under the Curve (AUC). We note that we penalise false positives and false negatives with the same severity and compute the MAE to be consistent with the current literature (Li et al., 2018; Zhou et al., 2013; Zhang & Jiang, 2016). We also note that it is important to use the F1 score and AUC, as they are more appropriate metries in the case of class imbalance, compared to average accuracy. Table summarises the description of each dataset and Tables 5a, 5b, and 5c summarise the results. The chosen rank (Fx) (eK) for ABIDE Autism, ADHD200 and Vessel MNIST 3D are (6, 6), (6, 6, 6) and (3, 3,3), respectively. The chosen LSR rank for LSRTR for ABIDE Autism, ADHD200 and Vessel MNIST 3D are S 2, S = and S = 2, respectively. The results show that the LSRTR method performs well in terms of scnsitivity and average accuracy for all datasets and the vector based-method LR has poorer performance with all datasets, where for some datasets all observations are predicted as positive. Particularly with the ADHD200 dataset, we are faced with the challenge of efficient estimation as the sample size is very small with respect to the model size (762 vs. 2,097,152). The relatively good performance of LSRTR suggests that LSRTR performs well in the high-dimensional setting (i.e., when the sample size is significantly smaller than the number of
23

---

Published in Transactions on Machine Learning Research (08/2023)
covariates). We notice that LR performs very poorly with ADHD200, there over 2 million parameters
to estimate. With the chosen rank, the LSR model decreases the dimensionality to 2322 x S + 216 which enables significantly more efficient classification. LTuR and LTR decrease the dimensionality to 2322 + 216 and 2322, respectively, but cxhibit poorer performance than LSRTR (and even SVM), suggesting that CP and Tucker models cause a loss of representation power that is essential to parameter estimation for some datasets. With the Vessel MNIST dataset we have more training samples than ADHD200 and lower- dimensional problem (only 21,952). For this reason all methods perform relatively well, with LSRTR having one of the best performances. When compared to ResNet 50, method approaches its performance in terms of F1 score and beats its performance in terms of average accuracy.
We make a remark on Table 5a. The sensitivity for LSRTR is 1, which is likely due to the fact that we use only 14 test samples from the ABIDE Autism dataset. In order to provide better interpretability, Figure 7 shows a histogram of the posterior probabilities for the 14 test data samples. We see that all of the samples labelled as y = 1 are classified as such since their posterior probabilities are all over 0.5. Most samples of Class 1 are classified with high confidence (in the sense that posterior probabilities are well above the 0.5 threshold). Only one sample of Class 0 was mis-classified.
To provide a more comprchensive comparison between LSRTR and other state-of-the-art neural networks, we also compare the performance of LSRTR on vessel MNIST dataset to a benehmark analysis (Yang et al., 2021b). We compare three ResNet 18 methods, with various additional convolution layers (namely 2.5D, 3D and ACS, or, "axiallcoronnal-ssiittallll" two additional ResNet 50 methods, with additional convolution layers (namely 2.5D and ACS) (Yang et al., 2021a)), and AutoKeras. The results are shown in Table 6. The results show that, although LSRTR may be outperformed by some complex methods, it exhibits comparable and occasionally superior performance to several others. Additionally, we note that while neural networks often outperform regression models in certain fields like computer vision due to their increased flexibility, their potential overfitting in small sample regime (limited sample sizes are often associated with medical datasots from clinical studics) remains a concorn. Morcover, regression models are frequently termed 'white- box models' owing to the transparency and interpretability of their parameters. This interpretability an indispensable attribute in numerous applications allows for a detailed analysis of covariates' statistical significance and the generation of confidence intervals (not just accuracy scores) as posterior probabilities qualities that are challenging to obtain with neural networks (Dreiseitl & Ohno-Machado, 2002; Issitt et al., 2022). Such confidence probabilities are important information for medical professionals to consider in medical studies, as they provide valuable insights into the classification results. Additionally, unlike the neural network methods, LSRTR is casily extendable beyond the 3D case, and is trivially applied to other regression types (such as linear and Poisson), and does not require fine tuning for each case.
0.8
0.7
Posterior Probablii
0.6
0.5
0.2
Class 0
Class
Figure 7: Histogram of posterior probabilities for test data in ABIDE Autism Dataset for the case of LSRTR.
24

---

Published in Transactions on Machine Learning Research (08/2023)
SVM
LR
LCPR
LTuR
LSRTR
Sensitivity
0.71
0.71
0.71
0.71
Specificity
0.14
0.71
0.85
0.85
0.85
Fl score
0.55
0.71
0.77
0.77
0.93
AUC
0.42
0.51
0.84
0.84
0.9
Average Accuracy
0.43
0.71
0.78
0.78
0.92
(a) ABIDE Autism
SVM
LR
LCPR
LTuR
LSRTR
Sensitivity
0.41
0.62
0.51
0.83
Specificity
0.78
0
0.34
0.52
0.4
Fl score
0.5
0.62
0.35
0.48
0.65
AUC
0.62
0.5
0.56
0.62
0.67
Average Accuracy
0.62
0.45
0.47
0.51
0.59
(b)
ADHD200
SVM
LR
LCPR
LTuR
LSRTR
ResNet 50 + 3D
Sensitivity
0.39
0.53
0.26
0.32
0.47
0.85
Specificity
0.95
0.55
0.946
0.94
0.96
0.86
score
0.44
0.21
0.3
0.37
0.55
0.57
AUC
0.84
0.52
0.6
0.66
0.81
0.9
Average Accuracy
0.89
0.55
0.869
0.87
0.91
0.85
(c) Vessel MNIST 3D
Table 5: Comparison of LSRTR with LR, LCPR and LTuR for diagnosis of test subjects in datasets:
(a)ABIDE Autism (b) ADHD200 and (c) Vessel MNIST 3D. In the case of Vessel MNIST 3D, LSRTR also compared to ResNet. 50.
RN18+2.5D RN18+3D RNIS+ACS
RN50+2.5D
RN50+ACS
Autokeras
AUC
0.74
0.87
0.93
0.75
0.9
0.773
Table 6: Neural network benchmark analysis for classification diagnosis of test subjects in dataset Vessel
MNIST 3D. RN stands for ResNet.
6 Minimax Lower Bounds for Tensor-Structured GLMs
In this section we derive lower bounds on the minimax risk of the LSR-TGLM problem'. We numerically assess the tightness of the derived bound in Section 6.1. Tn the results we will derive, the minimax risk depends explicitly on the parameters of the LSR-struetured tensor model, namely {πa]ke(k), {mg the distribution of X and the number of samples Specifically, we derive lower bounds ex using an argument for estimation problems based Fano's inequality described by Yu (1997). This approach relates the minimax risk of coefficient tensor B to multiple hypothesis testing problem. It states that if there exists an estimator (let us call it B X)), and which can be any parameter estimation algorithm),
Preliminary results appeared conference paper taki2021 minimmxx...
25

---

Published in Transactions on Machine Learning Research (08/2023)
with error matching the minimax risk ex (i.c., B X) B E*), then this estimator 2() can be used to solve multiple hypothesis testing problem (MHTP) (Khas'minskii, 1979). In this MHTP, the hypothesis set consists of collection, B2, of L LSR-structured coefficient tensors, where the goal of the estimator B(y, X),, is to detect the correct 'generating' L.Rt-structured tensor. Fano's inequality provides a fundamental limit/bound on the error probability for the multiple hypothesis testing problem. This limit in turn provides a lower bound on the minimax risk e*. As mentioned in Section 3, our analysis is local and our approach to deriving a lower bound on the minimax risk is information- and thus involves the analysis of the mutual information (defined as I)) between observations y = ... yn) and hypothesis E [L]:
F
I(y:1) & Ey, log f(y,1) f(y)f(0]
(37)
where f(:.)) and f()) are joint and marginal probability distribution functions respectively. We are now ready to state our main result.
Theorem 6. Consider the rank-(r1.... rkk) and separation rank S LSR- TGLM problem in (13) with n ii.d observations, (K, where vec(X,) ~ (0, Ex), and the true coefficient tensor 2 < d² The minimar risk gx is then lower bounded by
6² (1t) ) SKE(K mm 1)rs + IIEE (Tx - -1--1 -1
(38)
128M 12
where M is defined in Assumplion 2.
Theorem can be specialised to the Tucker and CP decomposition regression problems in the existing literature. That is, when 5 1 we have the Tucker model, and when rı 12 тк and G is constructed a diagonal tensor, we have the CP model. With these specialisations we obtain the following corollaries:
Corollary 1. Consider the rank-(ri..... (Tk Tucker-TGLM problem (Li et al., 2018; Zhang et al., 2020; Zhang& 8 Jiang 2016) with n i.i.d observations, (X, Us J where vec(X,) ~ NoO, Σx) and the true coeffi-
cient tensor B < d. The minimax risk for this problem is lower bounded by
2 (10 ExEjK((m: - 1)ra + IIEIKK (rk - 1) -1
(39)
128M Ezllgn
Corollary 2. Consider the rank-r CP-TGLM problem (Zhou et al., 2013; Ton et al., 2013) with n i.i.d observations, (X,,3i"" where vec(X,) (0, Ex) and the true coefficient tensor B < d. The minimax risk for this problem is lower bounded by
(10 E KEEK (m+ 1)r + (r 11
ex N
(40)
128M |E
Table 7 provides a summary of the minimax lower bounds from this work and the relevant current literature. Specifically, it shows the order-wise lower bounds on the minimax risk for logistic regression, linear regression and GLMs for several tensor structures on the model parameter. We make three remarks prior our discussion of Table 7. First, define the terms m [[Iε[K] my and F [kę[[]] T&. Additionally, oß is the variance of observation y, D is fixed constant and an upper bound on the second derivative of the cumulant function (a" (n)), and we place a dash (-) in the cells where there are no specific results in the existing literature. However, note that our derived minimax bounds for tensor-structured GLMs cover these gaps in the literature (namely the case CP-sttr structuud linear and logistic regression, and Tucker-structured logistic regression). Secondly, the order-wisc lower bounds reported for the unstructured case are in fact bounds on the minimax risk for prediction error, rather than estimation error. Nonetheless we include them in our comparison the sample complexities of estimation and prediction in given problem tend to be
26

---

Published in Transactions on Machine Learning Research (08/2023)
Structure of B
Regression
Unstructured
CP
Tucker
LSR
Linear
y Y
oy
E ke[k]
MATK
+
n
(Raskutti et al., 2011)
(Zhang et al., 2020)
Logistic
m (Abramovich & Griushtein, 2016)
L Dn
mar +r KE
Σ
MATK +F kE[K]
Σ
8 Σ MiTK + kE[K]
GLM
E2
M Ex 12
M || |2
(Lee & Courtade, 2020)
Corollary
Corollary
Thcorem 6
Table 7: Summary of order-wise lower bounds on the minimax risk for linear regression, logistic regression and GLM settings and for various tensor structures (unstructured, CP, Tucker and LSR). For unstructured problems, in [IEE[K] my is the total number of elements of the tensor. For Tucker and LSR models, IIE Tk-
proportional. Thirdly, here we have reported lower bounds on the minimax risk for non-sparse regression, as we are not studying sparsity in the scope of this work.
In the first two rows of the table, we report the order-wise minimax lower bounds for linear and logistic regression, In all cells of the first two rows, we see that the minimax risk decreases proportionally with increasing sample size. In the unstructured case, the minimax risk is proportional to the number of learn- able parameters in the model m.). For linear regression, however, this only applies if the number of parameters is less than n. This condition puts this result at disadvantage for the high-dimensiona setting under study. We can see that the minimax risk of linear regression in the Tucker-structured tensor case is proportional to the number of parameters in the Tucker model (Ex m +II rk). The numerators (Il, my) and (Ex mark + II rk) give insights into the parameters on which an achievable minimax risk might. depend, and thus give insights into the sample complexity of the parameter estimation problem. We can see that imposing low-rank Tucker structure may significantly decrcase the number of learnable paramcters.
The third row is the row of interest as it summarises results for GLMss Similar to the previous rows, the minimax risk decreases with an increase in sample size. We also have dependence on of and # Intuitively speaking, this means that the minimax risk increases with the variance of y. In the unstructured case the sample complexity is (II. my); and does not account for any consor-structured GLM settings. Theorem 6 and Corollaries and 2 make up the last three columns and address this gap in the literature. Therefore, our proposed minimax lower bound is unified bound, as it encompasses minimax lower bounds for the Tucker and CP-structured GLM models that were introduced in prior works (Zhou et al,, 2013; Li al., 2018), feature that has yet to be examined in existing literature. The minimax risk is proportional to the number of parameters for the CP, Tucker and LSR case. This is intuitively plcasing as it suggests that by imposing either the CP, Tucker or LSR structure on B, one can significantly reduce the sample complexity from the unstructured case to Σg m&r + r. Et MATK + II Tk, and S 'E mArk + II, Tk: respectively, and that we can develop algorithms (such as LSRTR Algorithm 1) that meet these bounds. As previously mentioned, the CP model induces the lowest number of parameters and the LSR model has (S - 1) Ex mark more parameters than the Tucker model. Additionally, our results also imply an inverse relationship between the
27

---

Published in Transactions on Machine Learning Research (08/2023)
minimax risk and |E and M, which is an intuitively pleasing outcome of our analysis. The variance of the covariates x (and thus [15xl1)) symbolises the signal power held by the GLM,, in the sense that the signal- to-noise ratio increases with increasing variance of x. To illustrate perspicuously, we can take the specific case of binary logistic regression: In this sctting, an increased variance of x induces a more varied natural parameter, VJ, which in turn causes easier distinction bet ween the two response classes. In a similar manner, smaller variance of induces a more difficult classification problem (the classes are harder to distinguish). In the extreme case, a variance of 0 causes all observations classes indistinguishabl. Therefore with increased variance of x (and thus ||Exll)). the estimation error should decrease. This argument is also consistent with various other minimax lower bounds in the literature (Barnes & Özgür, 2019; Shakeri et al., 2016).
collapse onto single point, making the
Lastly, though we do not provide any theoretical guarantees for the tightness of our bounds in Theorem 6 and Corollaries 1 and 2, we can see that our minimax bound for the Tucker-GLM meets that of the Tucker lincar regression works by Zhang et al. (2020), (we do not consider the r term in (Zhang et al., 2020) as that only accounts for the non-singular transformation indeterminacy, which do not discuss in the scope of this work). Zhang et al. (2020) have shown that their bounds are indeed optimal, which suggests that the minimax lower bounds derived in Theorem 6 and Corollaries and 2 are also tight. Morcover, we now provide a numerical analysis to assess the tightness of our bounds.
6.1 Tightness of Theorem
We utilise the results of our experiments synthetic data in Section 5 to investigate the tightness of the minimax lower bound in Theorem 6. We have shown that our result matches an existing bound in the literature for the specific case of Tucker linear regression. Now, we perform an empirical assessment of our result. Figure shows the ratio of the mean empirical error (over 50 repetitions) of LSRTR from our experiments in Section 5.1.2 to our obtained lower bound in Theorem 6. We do this for lincar and logistic regression and we can see the ratio is approximately constant as function of sample size for both regression types. Such a result suggests that bound may be achievable in the sense that we can develop algorithms that meet the minimax lower bound and take advanlage of the LSR tensor structure to lower the sample complexity of GLM problems.
1.25
Emipirical error to minimax bound ratl
LSRTR: (m,r) = (64.4)
Emipirical error to minimax bound ratlo
LSRTR: (m.r) = (64,4)
LSRTR: (m,r) (64,8)
1.2
LSRTR: (m.r) (128. 4)
LSRTR (m,r) = (128.4)
LSRTR: (m.r) (128.8)
1.15
1.1
1.05
0.95
0.85
0.8
7000
7500
10000
12800
15000
17500
Number of observations
Number of observations (n)
(a)
(b)
Figure 8: Comparison of the empirical error in LSRTR to the minimax lower bound for linear and logistic regression. (a) shows the ratio of the empirical error in linear regression the minimax lower bound for m 64, 128 and 4,8, (b) shows the ratio of the empirical error in logistic regression the minimax bound for m 64, 128 and 4.
28

---

Published in Transactions on Machine Learning Research (08/2023)
6.2 Roadmap for Proof of Theorem 6
We begin by introducing lew important concepts. First, in order to set up the multiple hypothesis testing problem consider the constructed packing set
B [B, E P(rntss: E [4] C (a(O), LEN,
(41)
from which the true index corresponding to the truc coefficient tensor is generated uniformly at random.
In order to ensure a tight bound on our minimax risk, B2 must possess the following properties: (i) The minimum packing distance between any two distinct hypotheses By E Br. 1, E [T], 1 + Specifically, for a strictly positive parameter 8 we require construction such that | B, F v88. (ii) The hypothesis testing problem is hard in the sense that it is difficult to detect the true index I (and thus the true coelficient tensor) based on observation y. (iii) The construction of B3 must induce lower bounds on the minimax risk that reflects the reduction in the sample complexity of the LSR tensor structure, i.e., we require a bound that exhibits a relation between the intrinsic degrees of freedom of the LSR tensor {mm*kee( (ri]ke k): and the number of samples (n). Secondly, since the objective of the multiple hypothesis testing problem is detect the index 1, it is solved through the minimum distance decoder
is large.
i(y)
arg min
B B F
442)
rE[L]
where B is estimated through any learning algorithm, such as Algorithm proposed in this work. In the LSR-TGLM problem we have the following hypothesis detection criteria:
If ||| B,
S V28: P((()) + O.
If
B B v 25: A detection error might occur, P((()) +1) 2 0.
The minimum distance decoder detects the true hypothesis if the estimate B lies within the open ball B(B,, V28), and detection error can only occur if bounded by
B B > V25. Thus the probability of error is
P((y)) +1) TP
B B F 2 V26
(43)
In this approach we interpret the coefficient estimation problem as a communication problem. From the set in (41), the source selects the truc hypothesis B, (the true coefficient tensor) uniformly at random. The hypothesis B, then generates the channel output according to 'channel model' as in (13) with a chosen link function g()). Specifically, the channel outputs y: such that BEyy #4 = و ((B X,)). The minimum distance decoder then successfully recovers the true hypothesis if the estimator B is within ball of radius 288 around B. Thirdly, we recognise that the problem under study is supervised and the response variables Vi Vi E [n] are conditioned on input data X, vi E [m]. On the basis thereof we opt to evaluate the conditional mutual information 1(y;4/8 2 ()y:
We use the aforementioned concepts to achieve our result in (38) through the accomplishment of the following tasks: 1) We must construct the set B2 as an exponentially large (with respect to the dimensions of an L.SR- structured tensor) family of L distinct LSR structured tensors, and satisfying the properties described above. This involves constructing KS + individual sets (for the KS factor matrices (B(R,.) kε(k), *E[s] and core tensor G) and deriving conditions under which all sets can exist simultaneously (with high probability). 2) We must find tight bounds on the conditional mutual information I(y; As explained in Section 3, 'Off-the- shelf packing sets and bounds on the minimax risk for the unstructured (vector-based), low-rank matrix, or low-rank Tucker-struetured parameter estimation problems are not useful in our application as these results do not capture the LSR tensor structure and would yield suboptimal lower bounds. Thus designing the right geometric insight into LSR tensors by constructing a packing that respects the factorisation structure and constraints (e-g., orthogonality) on the factor matrices is novel, our goal is to understand
29

---

Published in Transactions on Machine Learning Research (08/2023)
B3
28
d
B
B
B
V 86
B
BL
Ba(B,)
u.a.r
B
B
14 9 (4=gB((BB
(B,
min B(y, [E(C]
B,ll i
Source
Channel
Decoder
Figure 9: Information theoretic approach for deriving bounds on the minimax risk. Consider a packing set By. = (B; E [4]} C Ba(D) with minimum distance between any two elements as v8ō. The true coefficient tensor selected uniformly at random (u.a.r) and the LSR-IGLM model generates the channel output. We detect the true coefficient tensor using the minimum distance decoder. Here, the minimum distance decoder will detect B, as the true hypothesis.
the benefits of imposing the LSR structure on the problem's sample complexity. The novelty our work is that we explicitly leverage this structure leading a hypothesis B, with a structure that cannot be achieved through current methods. Another challenge in using this technique is actually bounding the relevant information-theoreet quantities. To that end, we derive divergence in tensor-structured GLMs. Our bound involves explicitly accounting for the link function with an LSR-structured coefficient tensor. This is a non-trivial derivation, as analysing the link function entails the careful construction of packing set and the derivation tight bounds on the packing distance as in Lemma The results Lemma 5 and the bounds on the KT--divereence can be used in other problems involving LSR tensors and distinguish our work from previous studies. Lastly, though our analysis is local and depends on a fixed tensor in a neighbourhood with known radius, our minimax risk becomes independent of this radius if it becomes sufficiently large.
tight bound the Kullback-Leibler
6.3 Proof of Theorem,
Four lemmas lay the foundation to the formal proof of Theorem 6. Through Lemmas 4, 5, and 7 we achieve the following: We construct the hypothesis set B,, containing distinet tensors. Since each B, E B, is LSR-structured, as defined (9), we construct B2 by constructing KS + individual sets that correspond [S]) and a set for G. Our construction of B, results in tight upper and lower bounds on the distance [B, B for any two distinct 1,1' E [T.], sampled uniformly at random (u.a.r), which will be used to derive bounds on the mutual information II(y; () and minimax risk e*. Complete proofs of these lemmas are in the appendix. Lemma (and Corollary 4) introduce the sets of 'generating" vectors (and matrices) from which we will later construct G (and ((...) for all [k], E [S]), that are in turn used to construct each B, Specifically, for some E N, aim is to construct set of (a)-dimensional
to KS sets for B(kks) (k E [K],
vectors with entries sampled uniformly at random from { [ 11
Ja with some designated minimum packing distance the Hamming metric between any two vectors. Topologically speaking, this can also be viewed as the construction of subset of an (a)-dimensional hypercube with required minimum Hamming distance between any two distinct elements in the set. Lemma and its subsequent corollary utilise the Gilbert- Varshamov bound constructing binary codes with minimum distance in order to construet the 'generating'
30

---

Published in Transactions on Machine Learning Research (08/2023)
sets we have just described. The following can be found in the book by Tsybakov (Isybakov, 2009, Lemma 2.9), where we restate the bound here in the interest of keeping this work self-containd..
Lemma 3 (Gilbert- Varshamoa Bound, Lemma 2.9 (Tsybakov, 2009)). Let d 2 8. Then there exists a subset C C {0, 1] of size
CC| 2 24/8
(44)
such that for any pair v, v E C, we have |V- 1-v11 Corollary 3. For any d 2 8, there enists a set of binary vectors C S a- of size |C| 2 24/8 such that for any pair v, v' EC, we have ||v v'lla 2 d g.
8.
Lemma 4 and Corollary 4 introduce a set of 'generating' binary vectors and matrices, respectively, with minimum distance.
Lemma 4. Let Tk > 0 Ak E [K] and = = IIEEK]] Tkk For any integers (Ti,T2, -. ,T)) such that 72 8, there erists collection of F 11/8 vectors {sy E --11-111--1 f [F]} such that for all f. E [F]
have ||s/ - sy 1l 2 (F - 1)/8.
Corollary 4. Let my > 0 and n > 0 4. E [K]. For any inlegers (mz,, such that (mi - 1)rs 2 8 there earists a collection of P = 2 11)r//8 matrices {Sp E [(my - 1)rk
[-1,1] (m 1)ra;
p [P]} such that for all p.p' E [P] we have ||Sp Sp l1o 2 (m² 1)ra/8.
The sets generated in Lemma 4 and Corollary 4 can be used to construct the set BL of L tensors in (41). We now introduce Lemma 5, which derives conditions on L such that the sets in Lemma 4 and Corollary 4 exist simultaneously. In Lemma 5 we construct B, with a certain set of properties. Note that every element in B, must have an LSR structure. That is, every clement is comprised of KS low-rank and orthogonal factor matrices (of dimensions my X Thy E [K]), and a core tensor (ol dimensions rι X .... TK), that we carefully construct from the set of "generating" matrices and vectors defined Corollary 4 and Lemma 4, respectively, according to (9). Additionally, any two distinet elements B, B, E Bt must be a suitable distance apart. Lemma 5 derives upper and lower bound this distance; the lower bound determines the minimum packing distance, while the upper bound is used to derive upper bound on the conditional mutual information in Lemma 6. The LSR tensor structure induced by having a common core tensor G is critical to such an analysis, specifically in deriving the lower bound in (46), which is a result of the algebraic steps shown between (74) and (75), in the proof of Lemma 5 in the appendix. Not only would these steps be considerably more challenging without the convenience of the LSR matrix structure of (14), but the effect of not achieving tight bounds on the packing distance may ultimately lead to minimax bounds that do not reffect the intrinsic degrees of freedom of the chosen GLM model.
Lemma 5. Define F Ππ rk- Let S > O, my > 0, rk > 0, and (mx - 1)rk 28 vk [K]. There exists collection of L tensors B2 e [B, :1 [L]] c Ba(0) for some d of cardinality
I. > 2²
(45)
such that for any
1 32 - 1)
d F 1- 3 7
<85
(46)
5
we have
F
F 7
82
2 <
B, - BB < 452
2
(47)
7
The bounds in (47) match (up to a constant) and are tight, which help ensure that the upper bound on the conditional mutual information is also tight. Until now, we have completed the tasks of constructing
31

---

Published in Transactions on Machine Learning Research (08/2023)
the set B,, with the packing distance 86 52
T 7- and with packing number (cardinality) parameterised
by the parameters frx (ra)Re(K), f ) S. We transition to the next task of deriving upper bounds the conditional mutual information, I(y: Upper bounds on L(y; ()2) involve the evaluation of the Kullback- Leibler (KL) divergence. For this we employ well-cstablished results in the literature (Cover & Thomas, 2012; Jung et al., 2016). Specifically, define fi(y|X') and Jr (y( as the conditional probability distribution of given X with any two distinct coefficient tensors B and By respectively. Denote DKL as the relative entropy, then
1 Σ T 2,1"
I(y;1 < D
(48)
We note here that we require the hypothesis test to be hard, which corresponds to relatively small KL divergence. Therefore, our aim is to derive tight upper bound on the KI. divergence in (48), by relying on the tight bounds in (47). Lemma thus derives an upper bound on [ly;; 4xx).
Lemma 6. Consider the LSR-TGLM problem given by the model in (13) such that B E Ba(0) for some d and consider the set By constructed in Lemma 5. Consider a i.id observations yi following probability distribution belonging to the exponential family when conditioned on vec(X,) -N(,,E)) Then we have:
F 7-
I(y; Mn Eull2
(49)
The steps between (77) and (80) in the proof of Lemma in the appendix elucidate how the analysis of the KT-ddivegencc involves both the link function and the bounds appearing in (47). Specifically, the KT- divergence between two distributions fi, fr in the exponential family is a function of the link functions m = (B,, X) and T (Br,X). X). Additionally, a tight upper bound on the KL-divergence requires a tight upper bound on 7t- n or alternatively, [B, Brll Since any B, is LSR-structured, we require the result in (47). The final step is to lower bound T(y; using Fano's inequality (Yu, 1997), stated below:
I(y: I(y: 2 (1 1 -((y +1)) logz - 1.
(50)
Evaluating (50) is simple. For this we refer the reader back to (43) in Section 6.2, which bounds the error probability P((()) A 1) for the recovery of hypothesis B, using the minimum distance decoder. Lemma 7 bounds this probability, with proof that follows exactly that for Lemma 8 by Jung et al. (2016), and thus is omitted in this work.
Lemma 7. Consider the minimumm distance decoder in (42). Consider also the LSR-TGLM regression problem in (13) with minimaa risk e*. Suppose e* 5 for some non-negalive scalor 8. Tf the L lensors distributed in Lemma 5 satisfy |B Brll 2 86 then the detection error of the minómum distance decoder is upper bounded by P(((y) # S
Proof of Theorem 6. Fix (f.. rk). From Lemma 5 for any > 0 satisfying (46) there exists packing
set Br. C (0) with cardinality T. and packing distance |B, Bp satisfying (45) and (47), respectively.
Also, from Lemma 6, with the set B2, I[y: satisfics (49). Now suppose there exists estimator B which guarantees a risk e
S2, - z If we set 86 = 4 p-1
s2 g2 then from Lemma for the set B₁ we have
(P((y)) f1)5 S J. Combining Lemmas 6 and have
1 log, L -
< 128MM |EZl,E",,
(51)
2
which be rearranged to achieve the result in (6).
7
32

---

Published in Transactions on Machine Learning Research (08/2023)
7 Conclusion
In this work, we investigate the LSR model tensor-structured GLM problems. Specifically, imposed low-rank LSR structure on the coefficient tensor in GLMs. The parameter estimation problem is highly non-convex, and we propose a block coordinate descent algorithm, called LSRTR, for this purpose. Fach convex sub-problem in LSRTR estimates a separate element/ /componen of the LSR-structured coefficient tensor. In our thcoretical analysis, we provide a minimax lower bound on the estimation error of the parameter estimation problem of LSR-structured tensors and specialise these bounds for Tucker-structured and CP-structured tensors. These bounds show that the LSR, structure reduces the sample complexity compared to that of the vector case. We evaluate the tightness of these bounds numerically and through the comparison of the specific case of Tucker regression to tight bounds in the literature. The methods we use may be of interest to readers as they can also be utilised for deriving minimax bounds on other LSR-structured estimation problems. Furthermore, we evaluate the LSRTR algorithm and the LSR-model on several synthetic and real datasets. These experiments demonstrate that the LSR model is less restrictive than other tensor models (such as Tucker or CP), and can be effectively employed for analysis on balanced and imbalanced medical imaging data. Some possible future work include thcoretical analysis of the tightnoss of the minimax bounds, and theoretical guarantees of the LSRTR algorithm.
Acknowledgments
This work was supported by the US National Science Foundation under awards CCF-1910110, CCF-1907658 and OAC-1940074, the US National Institutes of Health under award 2R01DA040487, and the Army Research Office under award W911NF-21-1-0801.
References
Felix Abramovich and Vadim Grinshtein. Model selection and minimax estimation generalized linear mod- cls. IEEE Transactions on Information Theory, 3721-337,,,, 6 2016. doi:10.1109/ TIT 2016.2555812. /ITI
Felix Abramovich and Vadim Grinshtein. High-dimensional classification by sparse logistic regression. IEEE Transactions on Information Theory, 3365):3068-300 3079, 5 2018. Nol:10.119/TTTT:200 201...201.220 2884499..
P.-. Absil, R. Mahony, and Rodolphe Sepulchre. Optimization algorithms on matriz manifolds. Princeton University Press, 2008.
ADHD-200 Consortium.
The ADHD-200 Consortium: A model to advance the translational po- tential of neuroimaging in clinical neuroscience. Frontiers in Systems Newroscience, 6:62, 2012. doi::0.3389/fnns ht0.3389/hhnss.2220000222. 2012.00062.
Talal Ahmed, Haroon Raja, and Wahcod U. Bajwa. Tensor regression using low-rank and sparse Tucker decompositions. SIAM Journal on Malhemalics of Dala Science, 2(4):944 966, 2020. Noi:10.1137119112990
Baiguo An and Beibei Zhang. Logistic regression with image covariates via the combination & and Sobolev regularizations. PLOS One, 156)):e0234975, 2020. doi: 1o:10..77//journal/ pon. 0234975.
Leighton Pate Barnes and Ayfer Özgür. Minimax bounds for distributed logistic regression. Technical Report arXiv:1910.01625 [cs.IT], ArXiV, 10 2019. URL pettps://rrxiv..om www.xi..or//bb//1/// 01625.
Pierre Bellec, Carlton Chu, Francois Chouinard-Decorts Yassine Benhajali, Daniel S. Margulies, and R. Cameron Craddock. The Neuro Bureau ADHD-200 preprocessed repository. Neuroimage, 144(Part B):275 286, 2017. 10t.1101....ioiiii. w0w.:oe...........0.
Quentin Berthet and Nicolai Baldin. Statistical and computational rates in graph logistic regression. In Silvia Chiappa and Roberto Calandra (ods.), Proceedings of the Twenty Third International Conference on Artificial Intelligence and Statistics, volume 108 of Proceedings of Machine Learning Research, pp. 2719 2730. PMLR, 8 2020. URL https://procediigg..... mlr grees//v10//erreeeeeeeeeeeeeeeee html. http:://ggrcccccnnnnoeeeeeee//er
33

---

Published in Transactions on Machine Learning Research (08/2023)
Dimitri P. Bertsckas. Conver Optimization Theory, volume 1. Athena Scientific, 2009.
Dimitri P. Bertsekas. Nonlinear Programming: 3rd Edition. Athena Scientific, Nashua NH, USA, 2016. URL
http: http://wwwww www.waaheeascc.ommm com/nonlinbook html.
Thomas M. Cover and Joy A. Thomas. Elements of Information Theory. John Wiley & Sons, 2012.
Cameron Craddock, Yassine Bemhajali, Carlton Chu, Francois Chouinard, Alan Evans, Andrés Jakab, Bud- hachandra S. Khundrakpam, John D. Lewis, Qingyang Li, Yan Milham, Michael, Chaogan, and Pierre Bellec. The Neuro Bureau preprocessing initiative: Open sharing of preprocessed neuroimaging data and derivatives. In Frontiers in Neuroinformatics: Conference Abstract: Neuroinformaties 2013, volume 7, pp. 27, Stockholm, Sweden, 8 2013. doi:10.338)) (o0.0///oo...n../2000010 2.1..09000....
Lieven De Lathauwer. Decompositions of a higher-order tensor in block terms Part 1: Lemmas for partitioned matrices. STAM Journal on Matrix Analysis and Applications, 30(3):1022-1032, 2008.
(bb:100.133700606666555
Stephan Dreiscitl and Lucila Ohno-Machado. Logistic regression and artificial neural network classi- fication models: A methodology review. Journal of Biomedical Informalics, 35(5 6):352 359, 2002. 11110000166/11522446608333333
John P. Dumas, Muhammad A. Lodhi, Batoul A. Taki, Wahced U. Bajwa, and Mark C. Pierce. Computa- tional endoscopy framework for improving spatial resolution in fiber bundle imaging. Optics Letters, 44(16)33968 3971, 2019. 1ed:10.1334/000 1...666.......00
Dylan J. Foster, Satyen Kale, Haipeng Luo, Mehryar Mohri, and Karthik Sridharan. Logistic regression: The importance of being improper. In Sébastien Bubeck, Vianney Perchet, and Philippe Rigollet (eds.), Proceedings of the 31st Conference On Learning Theory, volume 75 of Proceedings of Machine Learning Research, pp. 167-208. PMLR, 2018. URL actps://proceedingg .mr.............mm.............. mlr. wtt.:/ddwwn......peeeeee////eeee /ttp:///o//ofsorre.aaaa.. html.
Xiao Fu, Nico Vervliet, Lieven De Lathauwer, Kejun Huang, and Nicolas Gillis. Computing large-scale matrix and tensor decomposition with structured factors: A unified nonconvex optimization perspective. IEEE Signal Processing Magazine, 37(5):78-9) 94, 9 2020. doi:10.1109/MM93 2tt...222..033....
Mohsen Ghassemi, Zahra Shakeri, Anand D. Sarwate, and Waheed U. Bajwa. Learning mixtures of separable dictionaries for tensor data: Analysis and algorithms. IEEE Transactions Signal Processing, 68:33-48, 2020. doi:10.1109/ 1w1.10..19/TTSP... 2..9..99
Christophe Giraud. Introduction to High- figh-dimensional Statistics. Number 139 in Monographs on Statistics and Applied Probability. Chapman and Hall/CRC. New York, NY, USA, 2021.
Johan Hástad. Tensor rank is NP-complete. Journal of Algorithms, 11(4):644 654, 1990. doi: 101:10.001/0196 6774(90)90014-6.
Marti A. Hearst, Susan T. Dumais, Edgar Osuna, John Platt, and Bernhard Scholkopf. Support vector machines. IEEE Intelligent Systems and their Applications, 13(4):18 28, 1998. 4b::10.1109//5554.008222.
Hung Hung and Chen-Chier Wang. Matrix variate logistic regression model with application to EEG data. Bioslatistics, 14(1):189-209) 2013. doi:10.1093/ 1..10.1193//iossaassiiic. biostatistiics/kxs022
Richard W. Issitt, Mario Cortina- Cortina-Borja, William Bryant, Stuart Bowyer, Andrew M. Taylor, Neil Sebire, and Stuart A. Bowyer. Classification performance of neural networks versus logistic regression models: Evidence from healthcare practice. Cureus, 14(2):e22443, 2022. foi:10.7759/cureus.222 22443.
Alexander Jung, Yonina C. Eldar, and Norbert Görtz. On the minimax risk of dictionary learning. IEEE Transactions on Information Theory, 62(3):1501 1111110011111,,,-11 3 2016. doi:10.1109/THI /T00099/TTT.00 TTt.//TT.21161160 206.55177066
Rafail Z Khas' minkki.. lower bound on the risks of nor-maraettric estimates of densities in the uniform
metric. Theory of Probability & Its Applications, 23(4)7794 -948 1979. doi:10.1137/1 131111137/11 1123095.
34

---

Published in Transactions on Machine Learning Research (08/2023)
Tamara G. Kolda and Brett W. Bader. Tensor decompositions and applications. SIAM Review, 51(3):
155-500, 2009. doi:10. 1311111177070001111
Kuan- Yun Lee and Thomas Courtade. Minimax bounds for generalized linear models. In H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin (eds.), Aduances in Neural Information Processing Systems, volume 33, pp. 9372-9382. Curran Associates, Inc., 2020.
Xiaoshan Li, Da Xu, Hua Zhou, and Lexin Li. Tucker tensor regression and neuroimaging analysis. Statistics in Biosciences, 10(3):520-544 545, 3 2018. 661:10.0007//2255--00---111333
Muhammad Asad Lodhi and Wahood U. Bajwa. Learning product graphs underlying smooth graph sig- nals. Technical Report arXiv:2002.12777 [eess.SP), ArXiV, 10 2020. URL attps://ant //t::i///xivvvv/ia org/ /bs//////////
2002.11277,
Peter McCullagh and John A. Nelder. Generalized Linear Models. Routledge, 2019.
Frank Nielsen. Statistical divergences between densities of truncated exponential families with nested sup- ports: Duo Bregman and duo Jensen divergences. Entropy, 24(3):421, 2022. (ot:10.3390/e2403333
K. B. Petersen and M. S. Pedersen. The Matrix Cookbook, 11 2012. URL http://www2. compute dtu. dk/ (tbdb//pub//374-fullll full. html. Version 20121115.
Garvesh Raskutti, Martin J. Wainwright, and Bin Yu. Minimax rates of estimation for high-dimensionall linear regression over ba-balls. IEEE Transactions on Information Theory, 5710))6666--6994 10 2011.
11111.1199/TII..200 2011.22111165799
Garvesh Raskutti, Ming Yuan, and Han Chen. Convex regularization for high-dimenssonal multiresponse tensor regression. The Annals of Statstics, 47(3):1554-1584, 2019. 15:10.214/18-A0S1722
Athanasios A. Rontogiannis, Eleftherios Kofidis, and Paris V. Giampouras. Block-term tensor decomposition: Model selection and computation. JEFF Journal of Selected Topics in Signal Processing, 15(3):464-475, 4 2021. doi:
George A. F. Seber and Alan J. Lec. Linear Regression Analysis, volume 330. John Wiley & Sons, 2003.
Zahra Shakeri, Waheed U. Bajwa, and Anand D) Sarwate. Minimax lower bounds for Kronecker-structured dictionary learning. In Proceedings of the 2016 IEEE Int. Symposium Information Theory (ISIT), pp. 1148-1152, Barcelona, Spain, 7 2016. IEEE. doi:10.1109/ wtw.1109/IIIT.000 2222001161117777
Jianing V Shi, Yangyang Xu, and Richard G. Baraniuk. Sparse bilincar logistic regression. Tochnical Report
arXiv:l MrXiv:1404.4100 404444104 (math.OC), ArXiV, 2014. URL https htttps///rxxv...e wargg/b////aaaeeeeeeeeee 0bbs/1404.44.
Tingni Sun and Cun-Hui Zhang. Scaled sparse linear regression. Biometrika, 99(4):879-88) 898, 12 2012.
109:10.1093/ biomet/ biomet/ass043.. ass5043.
Batoul Taki, Mohsen Ghassemi, Anand D. Sarwate, and Waheed U. Bajwa. A minimax lower bound for low-rank matrix-variate logistic regression. In Proceedings of the 2021 55th Asilomar Confer- ence on Signals, Systems, and Computers, pp. 477-484, Pacific Grove, CA, USA, 10 2021. IEFE.
1oi::0.1109/1EEBBE 1LLL1///EEEEOEE3335502222222 2021193934452211972339999
Xu Tan, Yin Zhang, Siliang Tang, Jian Shao, Fei Wu, and Yueting Zhuang. Logistic tensor regression for classification. In Jian Yang, Fang Fang, and Changyin Sun (eds.), Intelligent Science and Intelligent Data Engineering, pp. 573-581, Berlin, Heidelberg, 2013. Springer Berlin Heidelberg. ISBN 978-3-642-36669-77 010-100077777---4-----6
Theodoros Tsiligkaridis and Alfred O. Hero. Covariance estimation in high dimensions via Kro- necker product expansions. IEEE Transactions on Signal Processing, 61(21):5347-5360, 11 2013.
11010.1100//STP22132227883
Alexandre B. Tsybakov. Introduction to Nonparametric Estimation. Springer New York, NY, 2009.
35

---

Published in Transactions on Machine Learning Research (08/2023)
Martin J. Wainwright. Information- limits sparsity recovery in the high-dimensional and noisy setting. TEFE Transactions on Information Theory, 55(12):5728-51 5741, 12 2009. 1tt.00.1109/IIA T0//TIIT220000 20..
Ying Wu, Dan Chen, Chaoqian Li, and Niansheng Tang. Bayesian tensor logistic regression with applications to neuroinaging data analysis of Alzheimer's disease. Statistical Methods in Medical Research, 31 (12):
2368-2382, 9 2022. 011100117700922202222111
Jiancheng Yang, Xiaoyang Huang, Yi Hc, Jingwei Xu, Canqian Yang, Guozheng Xu, and Bingbing Ni. Reinventing 2D convolutions for 3D images. IEEE Journal of Biomedical and Health Informalics, 25(8):
8009-30183 8 2021a. 1oi:10.1119/.JJ JBHII 2021.344424
Jiancheng Yang, Rui Shi, and Bingbing Ni. MedMNIST classífication decathlon: A light AutoML benchmark for medical image analysis. In Proceedings of the 2021 IEEE 18th Inter- national Symposium Biomedical Imaging (ISBI), pp. 191-195, Nice, France, 4 2021b. IEEE. doi:10.1109/ /100////SS12221122118822222
Xi Yang, Ding Xia, Taichi Kin, and Takeo Igarashi. IntrA: 3D intracranial aneurysm dataset for deep learning. In Proceedings of the 2020 IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), pp. 2656-2666, Seattle, WA, USA, 6 2020. 1ww.11/0///VR22200002202
Bin Yu. Assouad, Fano, and Le Cam. In Festschrift for Lucien Le Cam, pp. 423 435. Springer, 1997.
Anru R Zhang, Yuetian Luo, Garvesh Raskutti, and Ming Yuan. ISLET: Fast and optimal low-rank tensor regression via importance sketching. SIAM Journal on Mathematics of Data Science, 2(2):444 479, 2020.
Jianguang Zhang and Jianmin Jiang. Decomposition- tensor learning regression for improved classi- fication of multimedia. Journal of Visual Commanication and Image Represenlation, 41:260-277 2016. doi:10. 1oi:10.1016/j 2016..001...
Jianguang Zhang and Jianmin Jiang- Rank-optimized logistic matrix regression toward improved matrix data classification. Neural Computation, 30(2):505 -5555550--525 2018. 11111..11//encc a 01038.
Hua Zhou, Lexin Li, and Hongtu Zhu. Tensor regression with applications in neuroimag- ing data analysis. Journal of the American Statistical Association, 108(502))540-552 2013. dbi::0011800001211555 45.11...1....1..77969999
A Supporting Results
A.1 Proofs
Proof of Corollary 3. By Lemma 3 there exists packing with minimum distance & in the & norm containing 24/8 binary vectors. Mapping 0 + -a and 1 for all vectors in this packing shows that there exists a packing of at least 24/8 vectors in (-a, a3d with minimum distance da in the l norm. For any pair v, v in that packing we have ||v - v 2 $ since every entry of v v'v' is {-2a, 0, 2a).
Proof of Lemma 4. This a direct consequence of Corollary
The proof of Lemina contains the derivation a light upper bound on [BB By The following lemna is a component needed to achieve such an upper bound; specifically, we will see that it proves useful in dcriving an upper bound on || for any I E [Z].
Lemma 8. Let x E RM be any m-dimensional real vector and O mmxm be any orthogonal basis of T Define u, cos 0zl, where O, is the ongle between any possible and basis veclor os The function
f= f= (14 + 1)2 is minimized when u; = vm or when x is equiangular to all basis vectors o,, Vi E [m]]
-1
36

---

Published in Transactions on Machine Learning Research (08/2023)
m Σ (u 1)2. Additionally, for basis O Ommx and
Proof of Lemma 8. Consider the function f
E Rm 2cos(cos 1, thus have the equality constraint g Σ u = 0. Denote À the Lagrange multiplier, and thus the Lagrange function is defined as
Lg Σ (4 + 1)² X Σ u 1).
(52)
The partial derivative of bg with respect to ui, Vi [m]] is
2(u; + 2\us, vi E [m].
(53)
au;
The partial derivative of bg with respect to
aTg ig 02
m
u 1.
(54)
1=1
By setting (53) and (54) to zero we get
-1
"; =
Ei [m],
(55)
1 A
and
m
Σ
4 + 1 = 0,
(56)
respectively. What we have in (55) and (56) is a system of m + equations and m + unknowns. We solve the system to find the critical points of f()), which are u, = and u, 7m Since uy 2 0, the solution is u,
Vi [m] minimizes the function f(-). Vm
Proof of Lemma 5. Fix the following arbitrary real orthonormal bases: Q of RR,, and K sets of ræ bases, {Ux.) TT | of TM vk [K].
Next, consider the following hypercubes subsets thereof: 1) The set of F vectors {sy} from Lemma 4:
Sy
(57)
where f E [F], and 2) KS sets of P(c,s) matrices Ak E [K], Vs E [S], from Lemma 4:
Spona E { (ime - 1)rk (m= 1)r
-1
+1
(m² (mx-1)xrk
Wk E [K], S [S].
(58)
where P(k.s) E (P(k.))]-
We proceed with the following steps order to construct the final set B, of coefficient tensors from the sets in (57) and (58). Since By C (2))), we know that the energy of any B, is upper bounded by d? We will construct Gy; and matrices with hothonlllal columns, namely B Bp... Vk [K], E [S), all of which will be used to construct every B, e Br. Specifically, due to our LSR model, any tensor B will have rank (r1, P2,,2,..., LSR structure.
We use the notation (f,i) to denote the 4th step in constructing the fth element of G,- Hence in the first step, we construct vectors S(.1)) E R² for f E [F], using Q and sy, as follows:
B(..)) Q
vf [F].
(59)
Sf
37

---

Published in Transactions on Machine Learning Research (08/2023)
From (59), since ||yll| we have:
7
B(.
Q
r-1
F 1
5s
Now, we define each vector gy as
81 = Vi BU.))
(60)
for some positive number
Similarly, the notation (P(r,)):) i) to denote the 8th step in constructing the P(k.s) element of B Bpa.. Hence, we construct matrices B ()) E Immxx for Pik.s) E (P(x,)),, Vk E [K] 8 E [S]. Define B (Po..) as the jth column of B 1)), constructed as follows:
1), vk E [K], [S], and S the jth column of Spc.." Let the columns be
B (P(c.)):)) UKi
(Din.a) (P(..a)), Ke [K], E [S].
(61)
From (61) we have:
B ())0))))))0
== Th Px+1
wk E [K]-
rk
We now construct matrices B a). E mR for P(k.s) E [P(k,)) k E [K] E [S]. The construction of each B Poo.: follows the same procedure for all E [K] and s E [S]. Define B(J E Rm as the jth colurnn of Pik..)
B Bpuny for j E [rx]. We set
B (Pa..):))
(
(62)
P(kk))
B
(1)
and define
3
a (-11)
BB(++)
(++)
(P
Σ B(++1)
(P))o)))) P(L.,) Ple.)
B(5))
B B
(63)
j'ml
and
a++1
B (++1) P(k.a)
(64)
||
The steps in (62), (63) and (64) constitute the well-known Gram-Schmids process. Thus, the set of vectors
|| = and B(i) (P)).):1)) for j E [rw]s P(k.s) E [P(,.)), E [K] and s E [S] are orthonormal, ie., P((k.) 22
B(d) 1 LB')) P for any two distinct j.j E [rx]- Consoquently, (Bra...) (B (Bpt.... = I Now by defining
P(k.a
the set
C -
se(s) (F).)( E (P E [K], E [5]]
(65)
as the set of all tuples, (f.P(1.))... P(11.)..... .(K.,)) ... ,P(K.))) we have
L [C] (0)
2 (1/8)) ( (>-1) (m=
(66)
38

---

Published in Transactions on Machine Learning Research (08/2023)
where (a) follows from Lemma 4 and Corollary 4. We define the set of coefficient tensors, By as,
B1 B Σ G X1 B po..)
B Bpp...)
X B Bp(c,) :lE [4], P(k,), (P(,,)) e [K] [s]]
(67)
and we restrict such that
32 1) F < < ØV
F-1
(68)
S
F
We make the final note that, due to the Kronecker product, we can express vec(B,) as:
vec B,
Buxx. & Bpox-1.. 8.... B Bp )8f
(69)
We have the following remaining tasks at hand: 1) We must show that the energy of any B is less than between any two distinct d. 2) We must derive an upper and lower bound on the distance [ -B.-B.l13 11)))
tensors B, Bp, BL. We begin by showing 1 d::
B p
Σ G f x1 B Prr.)) X2 B Bpca. X3 ... X B p(x..)
B PIx..) & B Poc-110 ..... B BP)a)) gi
B B) 8 ... & B po.e)
|grr112
(70)
2
B Bpok.)
X B
Ba. 1F
||g/112
(71)
S2rE2
(c) d?
S2
B pn 1r |Br122
(72)
where (b) follows from the fact that for any matrix and any vector a, ||Aa|2 S ||A ||al|2 and the fact that || ||A|| (Petersen & Pedersen, 2012). Additionally, (c) follows the triangle inequality and (d) from the fact that the matrix norm of the Kronecker product is the product of the matrix norms. Additionally, (e) holds due to (68).
We proceed with deriving lower and upper bounds on B - Brl for any two distinct B, By Br. We first denote the square matrix Bpo.. E Rm= as the completed orthonormal matrix of each low-rank matrix B Bp mRr Also, G, E Rmix.xmm has entries G)(+) defined follows:
G (1: ri..., rg) G1: (1l::r,-.... (., l:rk) G (r(r + 1: mı..... .r. + 1: mx) = G G, FI. :tl.,,..................... rx).
(73)
1:
Also define B, S G X B Pa... x k B p...) for any 1 E [T... With these definitions, we have the
A=1
equality B, B = B B
Defining x
B P(n.a)
$ B Bpin. 8 ... x B Bpa..) for any e [L] and we
K
39

---

Published in Transactions on Machine Learning Research (08/2023)
have the following:
E²
8
k=1
B
- B
Σ
B k)
8(1.1) -
Σ
B Bli
((
r
km K
==1
A K
2
B pik.a
B(f.1)
+
B Br/nas
B(F
#=1
A= K
k=K
BUr..)) ))
2
B Bptk.))
B(..),
Σ
K
s=1
k=K
Define T2 = 8
k=1
B
p(r., V. =
8
B
for any E [S], and To,, v! as the jhh column of T, and V.:
k=K
k=
respectively, then we have
B, B 1/
8²
Σ
T
Tgg(111)
+
Σ
S
V ((((11)
T g (111)s Σ V V.B(r (.BBA.A)
F
2
8
Σ
T B((.1)
+...
T
B(S.1)
+
VTW B(''))
+ +...+
V (m)) (()
2
Σ
T
(())
B(r.1)
vi
V. B()r1)) + +...
Σ
T B(11))
((11))
VT'(m) B((),1)
We group every
(IS= TT (TBBBr)) + (Ix- v. v.T) 2
T
TTU (...))
(Is- vi B(r.))
for E [m]. We get
2
8 T
B(r .1) )r
B, B
(1)
m 8(r1)
Σ
V (1)
(74)
The expression (74) contains inner products. Specifically,
ΣS T (()
B(..))
= (ES T B
Denote A = cos <
T T())
B(..)) , then we have
===1
2
S
s Σ
2
B
B
2
A
(3)
1B
(0)
| [B1.
i=1
5=1
2
2
(g)
TE2 2
Σ
A
8 Σ 1=1
T
(2)
T())
2
7(F
1=1
a=1
(h)
762
m
52(), + +1)
7(7
(1) 7
S2
Σ S2
2 >
(75)
P(F
mm
where (f) is due to applying Cauchy-Schwartz Schwartz inequality to ES
B (g) due to (60), (h) is due TV(( are lower and upper bounded by S and S, respectively, Σg
to the fact that |>> T T())
and
40

---

Published in Transactions on Machine Learning Research (08/2023)
and (i) is from the result in Lemma 8. Finally, for finding upper bounds
B --,BBll
Brl we have:
B, B
2
B P(c.)) @ ... B ) gf
+
Σ
B M(ka)
8... ...&B B B(n)) 8P
8
2
+
B p/
||gr|11
sml
k
8
2
(k)
1l |18/|12
8=1
4
(76)
7
where (j) follows from the triangle inequality, and (k) follows from the fact. that B Bpoc ll =
B B
and that ||g = llg l2
Proof of Lemma 6, Consider the set By from Lemma 5, where the bounds in (75) and (76) hold. For the LSR-TGIM model in (13), consider n i.i.d samples, with covariate tensors X; E Rmmxx-xmxx Vi E [n]. where vec(X,) ~ (O, Eg). According to (13), observations y follow an exponential family distribution when conditioned on X, Vi [n]. Consider the vector of observations, y, and the tensor of n samples, X. Define also T(y;; 1122) the mutual information between observations y and index I conditioned side-information X. From (Cover & Thomas, 2012; Wainwright, 2009), we have,
1 EEE I(y; S L
(77)
LP
where (RI((i(y) (ylx) is the Kullback-Leibler (KL) divergence of probability distribution fi(y
and (yl))) of given X' for some B, BB E Br. Denote M,, and W. the link functions associated with fi( ( and (r (w(),,)) X;), respectively. Also denote Mr, the expectation of sufficient statistic T(x)
conditioned on X, under model B, (otherwise known as the canonical parameter). We evaluate the KL
divergence which is given as follows (Nielsen, 2022):
Σ (mu d) a(Tr)) + a(nr FE(m)
(78)
Now, we take the expectation of (78) with respect the side-information X. We have
Ex [m,, - -mr,)p a(m;) + a(wx)) = Ex [[ m))mu) due to the fact that [x(a(ma) Ex(a(nr))) We now have:
Ex ExDK2(((yy)))))06CC
Ex
B X.
(Br. X.) E(T()))) X.,
kE[m]
da(m))
Ex (BB -B,, x,))
Oms
iE(n]
aa(n)))
Σ
8x B., X Ex
(79)
Om:
kE(n)
< En 2
B By 1 M,
(80)
where (1) follows from the fact that μq = B[ X da(m.)) We achieve (79) through Cauchy-Schwartz
Jme
inequality. We make some remarks regarding (80): First, replace the summation samples with
41

---

Published in Transactions on Machine Learning Research (08/2023)
since each sample X; is independent. Secondly, Assumption 2 allows to bound
@a(m.) with M. Thirdly,
0m,
the conditions in (68) mean [ Brllr > thus B B 1
By Ir- Plugging in (80) into
(77) gives us
(m) I(y; n | B Bv M I
4
Ez 21 7-
where (m) follows from (76).
A.2 Numerical Results for Poisson Regression
Figure 10 reports the estimation and prediction accuracy for Poisson regression from the experiments on synthetic data discussed Section 5.1. The shaded regions depict one standard deviation of mcan estimation and prodiction accuracy, basod on 50 replications.
PR: m = 32
PR: m. 32
(m.r)
(m. (32. 4)
Nomnalise Estimation Error
TTR: (m.r) = (32,8)
TTR: (m. (32,8)
(m.r)= (32.4)
LSRIR: (m.r) (32.4)
0.5
LSRTR: (m.rr) (32,8)
0.7
LSRTR: (m.r) (32.8)
Normalised Squared
0.6
0.4
0.5
0.39
04
0.3
0.2
2500
5000
7500
10000
12500
2500
5000
7500
10000
12500
Number observations
Number observations
(a)
(b)
10
TTR: (m.r) (64,
TTR: (m. (64,4)
Normllised Esiimation Error
102
TTR: (m,r) (64.8)
Normalised Squared Emar
TTR: (m. (64. 8)
LSRTR: (m.r))(64., (64.4)
LSRTR: (m.r) (64.4)
LSRTR: (m.rr (64.8)
LSRTR: (m.r) (64.8)
ng/
10-1
10
Number of observations
Number of observations (n)
(c)
(d)
42

---

Published in Transactions on Machine Learning Research (08/2023)
PR: 128
PR: 128
TTR: (m..) - (128.4) 4)
TTR: (m, (128,4) 4)
Normalssed Estimation Error
10
TTR: (m..)
Nomalised Squared Error
TTR: (m. (128.8)
LSRTR: (m..) - (128,
LSRTR: (m. - (128, 4)
LSRTH: (m,r) (128.
LSRTR:: (m.r)
10
10
2500
5000
7500
10000
12500
2500
5000
7500
10030
12500
Number observations (n)
Number of observations (n)
(e)
(f)
Figure 10: Comparison of LSRTR with PR and TTR for two-dimmensional synthetic data when {32, 64, 128), (4, 8) and = 2. Normalised estimation error for 32, 64, and 128 is shown in (a), (c), and (e), respectively. Normalised prediction error for m = 32, 64, and 128 is shown in (b), (d), and (f), respectively. Each marker represents the mean normalised ostimation/p /predictionniiionnnn errors, over 50 repetitions. The shaded regions correspond to one standard deviation the mean normalised errors.
43


Reviews:
[{"review_id": "JoWB7VfjQt", "reviewer": "Reviewer_75Jv", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "The overall strength of this paper lies mainly in its execution and its presentation of experimental data, although there are significant concerns regarding the novelty of the ideas and methods presented.\n\nStrengths:\n\n1. The overall quality of the writing in the paper is commendable. However, I suggest that some sentences could benefit from careful rewording for enhanced clarity.\n2. The minimax lower bound for the proposed structured GLMs is proven with rigor. Even though my expertise does not extend into the realm of statistical learning, I found the proof to be understandable and illuminating.\n3. The authors have conducted extensive and detailed experiments using real-world data, enhancing the value of their findings.\n\nWeaknesses:\n\nThe key area where this paper falls short is in its novelty, or lack thereof.\n\n1. The proposed LSR model does not seem to offer anything new. It appears to be equivalent to the block term decomposition suggested by L.D. Lathauwer in 2006, albeit with what might be a minor variation in the constraint imposed.\n2. The algorithm used for estimation seems rather mundane and doesn't present anything innovative. Moreover, the analysis of its convergence is not adequately addressed, which could limit its applicability or validity.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "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": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This study explores a new low-rank tensor model, named Low Separation Rank (LSR), within the context of Generalized Linear Model (GLM) problems. The LSR model expands upon established tensor models, Tucker and CANDECOMP/PARAFAC (CP), by imposing itself on the coefficient tensor in the GLM model. The researchers propose a block coordinate descent algorithm for parameter estimation in LSR-structured tensor GLMs. Additionally, a minimax lower bound is derived for the error threshold in estimating the coefficient tensor in LSR tensor GLM problems, indicating that its sample complexity could be significantly lower than vector-based GLMs. These results also apply to lower bounding the estimation error in CP and Tucker-structured GLMs. The derived bounds align with the strictest bounds in the literature for Tucker linear regression. The effectiveness of the LSR tensor model is demonstrated through numerical experiments on synthetic datasets for three regression types: linear, logistic, and Poisson. Tests on a selection of medical imaging datasets indicate the utility of the LSR model over other tensor models (Tucker and CP) with real, unbalanced data with limited samples.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "While the paper is generally well-presented, some statements could be better articulated to improve clarity and comprehension.\n\n1. On page 5, “… the factor matrices $\\mathbf{B}_k\\in\\mathbb{R}^{m\\times{}r}\\,\\forall{}k\\in[K]$ are low-rank matrices”. I understand the authors would like to say that $r<<m$, but it does not imply the low rankness of a B. A similar assertion also appears on page 7 and might need to be revisited.\n2. On page 6, the phrase \"- which implies that the coefficient vector b is composed of separable sub-matrices weighted by some vector g.. — is also quite restrictive\" calls for additional explanation. It's not immediately evident that the Tucker model is more restrictive than the proposed LSR. The definition of LSR (as per Eq. (9) or Fig. 1b) suggests that the LSR could be modeled as a Tucker model, where the core is restricted to super-diagonal blocks. From this perspective, the relative restrictiveness of the Tucker and LSR models remains ambiguous.\n3. The claim on page 8 stating, \"....that since find the rank r of a tensor NP-hard, then finding the separation rank of a Kronecker structured matrix is also NP-hard\" appears non-intuitive. A more detailed discussion or a more comprehensive justification would bolster this argument.\n4. In the experimental section (page 26), the paper emphasizes the unbalance issues prevalent in real-world medical imaging data. However, the rationale as to why the proposed GLM should outperform other methods under these circumstances is not clear. Revisiting and elaborating on this point could strengthen the argument.\n\nQuestions: \n\n1. On page 9, the paper notes, \"the link function in GLMs makes part of our analysis non-trivial and fundamentally different to such works.\" Could the authors clarify which part of the analysis required non-trivial derivation, and what fundamental difference this brought to their work?\n\nTypographical Errors:\n1. In Equation (41), there is an extraneous closing bracket in \"$l'\\in[L]]$.\" This should be rectified for better readability.", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "h8S5grd8kP", "reviewer": "Reviewer_sup2", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "[Strengths]\n\nS1: A new tensor GLM model based on the LSR structure is proposed. \n\nS2: The derived minimax lower bound is new. \n\nS3: The performance of the proposed model is promising. \n\n[Weaknesses]\n\nW1: (Conceptual novelty). The notion of LSR tensor is not new and the combination of GLM and low-rank tensor/matrix structure has already been well studied.  Thus, the prosed model seems of little conceptual novelty. \n\nW2: (Technical novelty). The derived minimax lower bound seems a somewhat straightforward extension of the analysis for CP and Tucker models. It is unclear what are the technical difficulty in deriving this bound provided the well developed analysis for other low-rank tensor models. \n\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "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": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper considers the Low Separation Rank (LSR) tensor models in Generalized Linear Models (GLM) by proposing a block coordinate descent algorithm for parameter estimation in LSR-structured tensor GLMs. It also derives a minimax lower bound on the error threshold on estimating the coefficient tensor in LSR tensor GLM problems. Numerical experiments on both synthetic and real data show the effectiveness of the proposed algorithm.  ", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "Q1: It is suggested to highlight the conceptual novelty of the proposed LSR tensor models.  \n\nQ2: It is also suggested to discuss the detailed differences and difficulties in deriving the minimax lower bound in comparison with results for other tensor low-rank models. \n\nQ3: I am more interested in the estimation error of an iterative solution in the algorithm, i.e., the optimization error. Can you give some ideas about how to control/bound it?\n\nQ4: This is a general issue of tensor-based ML models. The low-rank tensor models may be interesting in statistical analysis. However, the empirical performance of general tensor-based ML models is not comparable to the DNNs in most settings. Although this paper gives some hopeful results on real, imbalanced and small-size data, I still appreciate it if the authors give more honest and calm discussions on the effectiveness of tensor-based ML models especially on the following points:\n\nQ4.1) Can tensor-based models be comparable or better than SOTA DNNs in popular classification tasks? \n\nQ4.2) If not, are there any settings where tensor-based learning is a must or very beneficial to the empirical performance? ", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "Xrl6o5HvYi", "reviewer": "Reviewer_zqEm", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strength:  The analysis of generalized linear model problem is comprehensive.\n\nWeakness: There should be a brief analysis of the proposed algorithm in the paper. In addition, the experimental section in the paper is excessively lengthy.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "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": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "In this paper, the authors impose the Low Separation Rank (LSR) model onto the coefficient tensor in the Generalized Linear Model (GLM) problem. They further apply a block coordinate descent algorithm to estimate parameters. In addition, they derive a minimax lower bound on the error threshold on estimating the coefficient tensor in LSR tensor GLM problems. Extensive experiments show the superiority of the LSR model.\n\nOverall, the LSR model is shown to be one effective and efficient representation in the GLM problem.\n", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "1.In the LSR model described by equation (9), the use of different matrices $\\underline{G}_{s}$ is not implemented. One might wonder if including such matrices would enhance the model's generalization capability.\n\n2.When using a block coordinate descent algorithm for the tensor model, the initialization should be close to the ground truth. What method do you use to choose initialization? Furthermore, how do you ensure, theoretically, that the initialization is close to the true values?\n\n3.To ensure the orthogonal structure of each matrix $\\bf{B}$, the authors employ QR decomposition. It is advisable for the authors to demonstrate that the QR step in the algorithm is non-expansive. Otherwise, concerns may arise regarding potential error amplification following the QR decomposition.\n\n4.The authors mentioned that the minimax risk increases with the variance of $y$. Intuitively, for the linear regression, when we increase the variance of input $x$, the variance of output $y$ will also be improved, further increasing the minimax risk. However, there exist $||\\Sigma_x ||_2$ in the denominator of (38). This means that the minimax risk is decreasing with increasing variance of $x$.\n\n5.Why is the normalized squared error corresponding to (m,r) = (256,8) in Figure 2(e) and (f) greater than that of 1?\n\n6.The reviewer suggests that the authors should give the meaning of $M$ in the theorem 6.\n\n7.The reviewer suggests that the authors should simplify the introduction of the experiment and keep necessary experimental results to increase readability. In addition, sections 5 and 7 should be combined.\n\n8.There are some typos: (1) in the beginning of Problem Statement, the dimension of $\\underline{X}$ and $\\underline{B}$ should be $m_1\\times \\cdots \\times m_k$; (2) in the description of equations (24)-(26), (24) should be (23); (3) in the Algorithm 2, there are extra square brackets in $\\hat{\\mu}$ and $\\underline{\\mathcal{X}}$; (4) in Theorem 8, \"destributed\" should be \"distributed\".; (5) The equation (79) should be aligned.\n", "broader_impact_concerns": "No", "additional_comments": "", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: This work considers a low-rank tensor factorization model called low separation rank (LSR) for generalized linear models. The authors provide a comprehensive analysis of the resulting nonconvex problems involved in estimating the tensor factors. Specifically, they derive a tight minimax lower bound for estimating the coefficient tensor, which captures the intrinsic degrees of freedom in the LSR tensor. The work also proposes a block coordinate descent algorithm to solve the nonconvex problem and estimate the tensor factors. Numerical experiments on both synthetic and real datasets demonstrate the usefulness of the LSR model for various generalized linear models.

The main concern raised by the reviewers is the clarification of the novelty, as the LSR tensor model is not new. To address this comment, the authors have revised the paper and clearly highlighted the major contribution, which is the introduction of the LSR tensor model for generalized linear models rather than the LSR itself, and the comprehensive analysis of the resulting nonconvex problem. All the reviewers acknowledge the revisions and the contribution of the paper, and recommend acceptance or leaning towards acceptance. I agree with the reviewers and recommend acceptance. 

Very minor comment: both $\Sigma_x$ and ${\mathbf\Sigma}_x$ (bolded) are used in Section 6. The authors may unify them. 
Decision: Accept (poster)