A Highly-Efficient Group Elastic Net Algorithm with an Application to Function-On-Scalar Regression
Authors: Tobia Boschi
Matthew Reimherr
Francesca Chiaromonte
Venue: NeurIPS 2021 Poster (2021)
Keywords: Feature selection
Group Elastic Net
Function-on-Scalar regression
Functional Principal Components
Augmented Lagrangian
GWAS

Abstract: Feature Selection and Functional Data Analysis are two dynamic areas of research, with important applications in the analysis of large and complex data sets. Straddling these two areas, we propose a new highly efficient algorithm to perform Group Elastic Net with application to function-on-scalar feature selection, where a functional response is modeled against a very large number of potential scalar predictors. First, we introduce a new algorithm to solve Group Elastic Net in ultra-high dimensional settings, which exploits the sparsity structure of the Augmented Lagrangian to greatly reduce the computational burden. Next, taking advantage of the properties of Functional Principal Components, we extend our algorithm to the function-on-scalar regression framework. We use simulations to demonstrate the CPU time gains afforded by our approach compared to its best existing competitors, and present an application to data from a Genome-Wide Association Study on childhood obesity.

A Highly-Efficient Group Elastic Net Algorithm with an Application to Function-On-Scalar Regression
Tobia Boschi Department of Statistics Penn State University tub37@psu. edu
Matthew Reimherr Department of Statistics Penn State University m1r360psu. edu
Francesca Chiaromonte Department of Statistics Penn State University,
EMbeDS
Sant' Anna School of
Advanced Studies fxc110psu. edu
Abstract
Feature Selection and Functional Data Analysis are two dynamic areas of research, with important applications in the analysis of large and complex data sets. Strad- dling these two areas, we propose new highly efficient algorithm to perform Group Elastic Net with application to function-on-see on-scalar feature selection, where functional response is modeled against a very large number of potential scalar predictors. First, we introduce a new algorithm to solve Group Elastic Net in ultra- high dimensional settings, which exploits the sparsity structure of the Augmented Lagrangian to greatly reduce computational burden. Next, taking advantage of the properties Functional Principal Components, we extend our algorithm to the function-on-scall regression framework. We use simulations to demonstrate the CPU time gains afforded by our approach compared its best existing competitors, and present application data from Genome Wide Association Study on childhood obesity.
Introduction
As problems involving very large and potentially structured data become ever more ubiquitous, attention is being devoted to the integration of approaches and techniques from the areas of Feature Selection and Functional Data Analysis (FDA). Indeed, more and more regression applications comprise large number number of variables some of which scalar and some of which suitable for a functional representation, such as longitudinal measurements or biomedical images (Sorensen et al. [2033] (llah and Finch 2013] Cremona et al. [2019) A great deal of recent work has been concerned with feature selection in these applications. Matsui and Konishi (201 Gertheiss et al. (2013); Fan al. (2015) study the case where the response scalar and the features are functional. Chen et al. (2016) 6 Fan and Reimherr (2016); [Babber et al. )); Parodi et al. (2018); Mirshani and Reimherr (2019) tackle the called function-on-scalar case, where the response is functional and the features are scalar focusing on settings in which the number of features is bigger than the number of observations. However, recent developments in optimization have demonstrated that substantial computational gains can still be made when the number of features is massive (e.g ~16)). In this article, we present Functional Group Elastic Net (fgen), novel and highly efficient method settings where
to solve the function-on-scan feature selection problem in Altra-hig-ddimension the number of features indeed massive and much larger than the number of observations. The ability to solve these problems with lower computational burden is increasingly critical. Given the complex, noisy nature much contemporary data, changing some aspects of their pre-processing or some of the tuning parameters involved in the analysis can lead to completely different results (Krawezyk and Cano 2018 Murdoch et al. 2019). For this reason, repeating an analysis multiple
35th Conference on Neural Information Processing Systems (NcurIPS 2021).

---

times (e.g., with different choices of data preprocessing pipelines, or to tune certain meta- parameters) paramount to capture significant signals and ensure the stability of outcomes (Yu and Kumbier 2020). Substantial reductions in computational burden enable such repetition, allowing scientists and practitioners to conduct truly meaningful and reproducible analyses.
Group Elastic Net incorporates the group structure (Yuan and Lin 2006) and the Elastic Net penalty (Zou and Hastie, 2005) into a penalized regression framework. The former allows one represent each feature (or component) by group of variables. The latter induces sparsity and regularizes the estimates. We consider the case where all groups have the same size k. The minimization problem is formulated as follows:
P min (1/2) (XB + A1 E1BBBll (A22222(22/2) (22/2) Σ B.ll
(1)
B
Let p be the number of features, the number of statistical units, and ll the 12 norm for matrices, i.e the Frobenius norm, and vectors. Then, X E R"XP is the design matrix (that we assume to have standardized columns), Y E T the response matrix, and B E TR the coefficient matrix. In other words, (I) describes a sparse multi-task model where the response and each of the features are represented by a group of & coefficients (Zhang and Yang 2018). Throughout this article, we follow the notation in Johnson et al. (2014) and we use the subscripts and (i) to indicate the i-th row and the i-th column of a matrix, respectively. Thus, B; E IR are the coefficient values associated with the i-th group, and X()) E 12d are the observed values relative the i-th feature. Before proceeding. note that ① can be expressed as
min (h(XB) + +B))) B
where = = - Y12 is the least-squares loss function and TTB) - Ei (B,)
_ is the Group Elastic Net penalty function. The first term in π not differentiable and creates sparsity at the group level, i.e, if component selected, then all its coefficients are selected and vice-versa. The second term is Ridge-type penalty which reduces model complexity and tries to control variance inflation due to multicollinearity feature selection models are indeed known to be less effective and not reliable in scenarios characterized by very high collinearity among features (Katrutsa and Strijov 2015). a1 and X2 are penalty parameters > 0 and control the weight of the two penalties with respect to the least square loss.
A1 E B (>e/2 [
To solve M.. we develop a new Semi-smooth Newton Augmented Lagrangian (SsNAL) algorithm. We then extend it to the function-on -oniiiol----ccllrr regression framework by means of Functional Principal Components (FPC) (James et al.] 22000] [Chiou al.) [004]] Hall and Hosseini-Nasab) (006). S&NAL exploits the sparsity induced by the augmented Lagrangian second order information to guarantee a super-linear convergence and greatly reduce the computational cost. This methodology, first introduced by Tomioka and Sugiyama (2009) and Tomioka et al. (201T), has been recently used in several applications, e.g., to regular Lasso (Li et al.) 2018), constrained Lasso (Deng and So 2019), and Elastic-Net (Boschi et al. 2020). However, incorporating the group structure significantly increases the dimension of the problem. Indeed, (I) is not separable and the optimization must be carried out jointly across the coordinates of the outcome. Therefore, considering new group penalty while preserving the efficiency of the method requires that we carefully redefine a set of all-new mathematical operators and the theory behind them.
We implemented an efficient version of fgen in python and benchmarked against the two best Group Elastic Net solvers we found in the literature: the python package sklearn (Pedregosa) jet al. T) and the R package glmnet (Friedman et al.) 2010), which is written in fortran Both of these solvers implement a highly optimized coordinate descent algorithm (Friedman et al.
2010
Breheny and Huang 2015) and outperform competitors such us FISTA (Beck and Teboulle
2009 Bonnefoy et al.) 2015), ADMM (Deng et al. 2013 Zhu 2017), and proximal gradient (Chen ct al. 2010) by at least one order of magnitude in terms of CPU time. Our simulation results
demonstrate that in sparse scenarios fgen is at least 3 times faster than glmnet and more than 10 times faster than skl earn. We also applied fgen to the Intervention Nurses Start Infants Growing on Healthy Trajectories (INSIGHT) study (Paul et al.) 2014),, which investigates risk factors for childhood obesity. Specifically, we examined the association between hundreds of thousand of Single Nucleotide Polymorphisms (SNPs) and growth eurves, which represent functional outcome.
The remainder of the article organized as follows. In Section we describe the Group Elastic
Net problem and introduce some preliminary results. In Section 3 we present our new methodology
2

---

and illustrate how to extend it to a function-on-se -n-alr feature selection problem. In Section 4 we investigate the performance of our method on simulated data and apply data from INSIGHT. In Section 5 we provide final remarks and discuss future developments. Proofs of theoretical results and additional simulations are included in the Supplement. The fgen code available at https://gg gitggtttuu..... www//oobbbabschh/eoeeeee.e..
2
Preliminaries
In this section we define the Group Elastic Net problem and we introduce some results related to Fenchel conjugate functions and proximal operators, which are essential tools in our developments.
2.1 Fenchel conjugate function and proximal operator of TTB)
Fenchel conjugate functions (Fenceel 1949) allow one to more readily define the dual problem
(Boyd and Vandenberghe 2004) of (P) which is called the primal problem. Let X' S R² be convex set and f X R. Then, the conjugate function of f X - E defined as
(=) suP=e> ((2, the inner product. i.e. the dot product. If X' C P i.e. if is matrix, the definition is still valid but ( is the Frobenius inner product. In our first proposition we provide closed form solution for the Group Elastic Net penalty conjugate function (see Supplemental Section A.1 for a proof).
- where X = {2 suPzex ((2,z) - f(z)) < <0 (..)) indicates
Proposition 1. Given Z E RRxx,, the conjugate function of π has the form
x (2)
- =2x))
Σ (110112~1144)
where [ ]+ is the positive part operator; [s]+ = 8 if > 0 and 0 otherwise.
Note that x² (Z) is continuous differentiable function. This is more general results than the one presented in Li et al. (2018); Boschi et al. 2020). because we extend the definition of T to the case where 2 is matrix and not just vector. In the simple scenario where k 1, i.e. when every group consists of just one variable, we obtain again the conjugate function of the standard Elastic Net penalty. Notably, starting from non-separable objective function, derive x² which does separate. As we will see in Section 3.1 this allows one to induce a new level of sparsity in the Lagrangian problem which actually key for the massive computational advantage offered by fgen.
Proximal operators (Rockafellar 1976a b) are fundamental in many optimization algorithms. Given a lower semi-continune convex function f IP + R. the proximal operator of f at x with parameter >> 0 is denoted prox,, R - R and defined prox, (2) arg min, (f( (2a)-1 =112). If f Rpxk + IR, i.e. if i is a matrix, then prox,; IR + ggxk and ||-11a is the Frobenius norm. Parikh al. (2014) and [eck (2017) (Chapter 6) provide numerous examples and properties. Combining their results, one can easily find the form of the proximal operator of (B) provided in our second proposition (see Supplemental Section A.2 for a proof).
Proposition 2. The proximal operator of T (BB) is: prox, n (B) (prox,, (Bi).... (Bro),((Bi.... where
prox, (Bp))
prox (B,) = (1 (12 1 |B,11 ad1 B,
(3)
Note that prox, (B) : Rxpx + Rpxk, To implement fgen, one also needs the proximal operator of T²" which can be obtained through the Moreau decomposition:
propppoxxxxo (x) + (=0/=+/o(//0) 0> 0.
(4)
2.2 Dual formulation and Augmented Lagrangian
Here we introduce the dual Group Elastic Net problem and its augmented Lagrangian. From Boyd
and Vandenberghe (2004), a possible dual formulation of (P is
min (v)) s.t. x =
(D)
V.2
3

---

Algorithm Semi-smooth Augmented Lagrangian (SsNAL) method
Augmented Lagrangian method
Semi-smooth Newton method for 7
Start from the initial values V - B 0
solve 7 and find (V 2
while not converged do
while not converged do
(I) Given B", find V*+1 and = which ap- proximately solve the inner subproblem
(1) Find the descent direction solving exactly
by conjugate grudient the linear system
(Vk+1, z
ss min Lo (V, z B²) (7) V,2
8 vec")) = vee (Vé(V*))
(9)
(2) Line search (Li et all 2018); choose P
(2) Update the Lagrangian multiplier B and the parameter o:
(0, 1/2) and reduce the step size until
1 (V* + + D') S (VV) + us (VN(V>), D))
Bh+ B² -a -x
z2
(3) Update V: V:1+ = V D'D
(8)
t 00
(4) Update z prox ( - x
end while
end while
where V E R and Z E Rp*k the dual variables matrices. In particular V,, Zj R are the dual variables associated with the i-th group. h" and π are the Fenchel conjugate functions of h and t, respectively. Specifically, A (VV (1/2) | + (Y, V) (Dünner et al. 2016) and x" xZ) is given in Proposition We can now define the augmented Lagrangian function and the Ktarush-Kuhin Kahn-Tucker (KKT) system associated with (). The augmented Lagrangian is given by
C. (V, +x)) Σ (B. VT X(0 + Z.) + (a/2) Σ C Zill (5)
where a > 0. both the primal variable and the Lagrangian multiplier which penalizes the constraints' violations. The KKT system is given by the following three equations:
V XB= 0, 0, 0= 0= (Z) B=0,
X 7
(6)
Note that (V) V + Y. A closed form of (Z) is not essential for SsNAL method. The KKT equations will be useful to determine the convergence of our algorithm, since the set (** B*) solves the KKT (6) if and only if (V*, 2²) and 8 are the optimal solutions of (D) and
respectively (Boyd and Vandenberghe 2004).
3
Methodology
In this section we present our new methodology. First, introduce SsNAL algorithm to solve the Group Elastic Net problem. Next, we illustrate how to extend it to the function cecction-onscaal on-scalar regression framework. Finally, we describe how implement solution path over different values of Å₁.
3.1 SsNAL method
The S&NAL method is summarized in Algorithm I consists an Augmented Lagrangian method characterized by inner subproblem. The subproblem is solved with Semi-smooth Newton method which exploits the sparsity of the augmented Lagrangian second order information and greatly reduces computational costs. We now provide the details of its implementation and some important theoretical results. From Rookafellar ( (1976)). one can find the optimal solution of (D) by solving the Augmented Lagrangian method described in Algorithm[] [ The essential part of the algorithm is the subproblem
7 As described in Li et al. (018), an approximate solution (V, 2) for a given B can be found as
= arg min L. (V 2, B), = arg min C. V2 2
(10)
V
With slight abuse of notation, we indicate by L. (VIZ, B) the function L, (V, B) where the parameter Z and B fixed. Similarly for L, V, B). Our third proposition provides explicit forms for C. (V Z, B) and Z (see Supplemental Section A.3 for proof).

---

Proposition 3. Define (VV = Lo (V Z, B). Then, for the Group Elastic Net problem we have
(a) (v) =((() 1 σ2
Σ
p
prox,x
oV (Bi-oV"X; (1))
1
p
|B:112
20
20 i=1
(11)
(b) 2 - ppoxx+ /a (B/o - V)),
where po (B/o X = (prox,- (Bi/o - V*)a)) po (Bp/a V
Z has closed form. To find V one has to minimize to or, equivalently, find the solution of V 0. Note that v is continuous and differentiable, and thus Vi is well defined.
To solve the subproblem
7 we propose the Semi-smooth Newton method in Algorithm | and
according to the rule in Proposition 3 and V by minimizing & through one Newton step. The main computational cost is solving the linear system ( This leads to our next crucial result (see Supplemental Section A.4 for a proof).
Z updated iteratively
Theorem 1. Let T = B - oXTV, X = X @ ly (the nk x pk Kronecker product between X and the k x k identity matrix). O2rb be the generalized Hessian of v, and 0 prox, be the Clarke sub-differential of prox,= (Clarke) 1990). Then we have
(V)(()) = =V VY X ox (T) (())&(V))
Ink prox (T)
(12)
Pi
Moreover, let Q E Rokkpp be the block- diagonal matrix Q =
where each P, is
Pp.
squared is x k matrix defined
(1+0 (1- ||T. aA1) - 1 oATTTT))
1T > 0d1
(ii)
(13)
0
o.w.
Then Q E ô prox,- (T) and 8 (V) vec(D) = (Ink + =xQXT) vec(D) for every D Bx in the domain of V - where vec( D) E Rin obtained by stacking all the columns of D.
Note that, while in Li et al. 2018): Deng and So (2019); Boschi et al. (2020) V and 8² are vector and a matrix, respectively. here the dimensions of these operators increase due to the group nature of the problem. In particular, Vei becomes a matrix and 82 higher order tensor which we express as an nk X nık matrix by stacking its dimensions. Moreover, Q is not simple diagonal matrix as in the previous SsNAL algorithms, but is now characterized by blocks associated to the different groups of variables.
Theorem is critical for preserving the efficiency of fgen, while integrating groups into the problem. First, it states that solving (9) is equivalent solving (Ink + xxQX") vec( vec (V)( Second, the form of Q still allows one to induce sparsity in the linear system and drastically reduce the computational cost. Indeed, let J {i : [22 2 σ and let [J] be the cardinality of J. Then the linear system (9) is equivalent to
+ Tok+oxxQu0
Qx vec(D) vec (VV(V)).
(14)
Ls. with Xy gaxr being the sub-matrix of X restricted to the columns in J. In addition, Q3 E Brrxrk is the block-diagonal matrix formed by all the P. such that E J. Using the Cholesky factorization (the generalized Hessian is positive
X, Here, X- Bknk is defined as x==
semidefinite) the total cost of solving the linear system reduces from O (nk (n² + p² + ni)) to O (nk33 + + nr)). This includes computing X which is O (nrk3 (n + r)), and the Cholesky factorization, which ( Because the sparsity induced by the Group Elastic Net penalty, is usually much smaller than p implying substantial computational gain. Even when p is very large (~ 10%) one can still solve the linear system efficiently, as long as the dimension k of each group is relatively small (< 102). Furthermore, if r < n, which is often the case when
5

---

Algorithm 2 Functional Group Elastic Net method
(1) Perform FPC of B and find the first & basis components (76, 7w)) with their eigenvalues (p pk)
(2) Find the first & FPC scores for each response function Vy: Y: = (i)i,,)........ (Di, 7h)2)
(3) Using Y as response matrix, apply SsNAL to solve ① and find the coefficient scores estimates B.
(4) Project B into the FPC basis to find the coefficient curve estimates: B, Σラ= B(,,5)
the solution of the Group Elastic Net problem is sparse, one can factorize an rk x rk (instead of nk x nk) matrix using the Sherman-Morrison Woodddbryy formula (Van Loan and Golub 1983):
(Ii ( + +xx3xx3333
Inn -Xa ((095)'+x338))
X3.
(15)
The total cost is further reduced from O (n1=3(n² 22 nr)) O (rk- + r2 + + 1)), in- cluding the computation of Q77 which can be done with a cost of ((rk2) by inverting each one of the P. blocks independently. Finally, if in the first iterations of the algorithm and r are both larger than 10², one can solve ( approximately using the conjugate grudient method (Polyak) 1969)
To determine the convergence of the Augmented Lagrangian and the Semi-smooth Newton methods, we check the residuals of the third and first KKT in (6) respectively, i.c.
E E+
VT (1) + 2:112
res(kk)) Di=1 1 + Y. - X,BI1, Σ |2
res(kkts) 1+ 11
(16)
VI 1
Taking the 12-norm of the KKT residuals, normalizing them and using them to assess convergence a common procedure the literature (Li et al., 2018 Deng and So. 2019). Both methods have super-linear convergence rate. Accordingly, the convergence rate of the entire algorithm, which the sum of the convergence rate of the two sub-problems (Tomioka and Sugiyama 2009), still super-linear. Thus, we show in Section 4 fgen typically converges in very few iterations. The convergence analysis here follows directly from that in Tomioka et al. (20TT) and Boschi et al. (2020), where SsNAL convergence is proved for the standard Elastic Net. The proof leverages results in Rockafellar 1976a b): Luque [984);; [L al. (2018) and the fact that π is a continuous differentiable function (wiich is also true for the Group Elastic Net).
3.2 Extension to function-on-scalar regression
We now extend fgen to the function-on-scalas features selection problem. In function-on scalar regression, a functional response is regressed on a set of scalar predictors. Assuming the response belongs the Hilbert Space [[a b)), the optimization problem ①
becomes
min
(1/2) VI + A1 Σ |B + (A2/2) Σ | B
(17)
1=1
where y and B are functional objects with and p rows, respectively. Each row )1 is response function and each row B; is coefficient function. The squared L of function f is II
(f.,,,, where the inner product between two functions ff and g is (f.g) = So fg.
Applying SsNAL directly (17) is not straightforward and would substantially hinder its efficiency. First, the definition of conjugate functions and proximal operators in functional spaces would require theoretical background. Second, and perhaps most important from practical standpoint, computing integrals is much more expensive than computing euclidean norms. For these reasons, in Algorithm 2 we take advantage of Functional Principal Components (FPC) (Horváth and Kokoszka (0017) to solve an optimization problem of the same type as (I) which (7). In particular, we build response matrix Y. where each group is formed by the first l FPC scores of the the function V, (Fan and Reimherr (016)). The level of approximation of fgen thus depends on the number of FPC scores k, i.e the dimension of each group. Indeed, given function f and its FPC basis (%a 10c 1. we have
2012 Kokoszka and Reimherr
fact a very close approximation
11 Σ IIF rheell
(18)

---

This property- which is true for every orthonormal basis system -plays crucial role in the extension of our SsNAL approach to the function-on catioon-o--ccllat regression, since it allows one to approximate the 1.2 function norm with the standard I2 matrix norm, Consequently, one can use the FPC scores to construct the response matrix Y and the coefficients matrix B in (T starting from the response functions Y' and the coefficient functions B. The number of selected FPC scores determines the dimension &or of each group in the Group Elastic Net problem. In many applications just few FPC scores allow one to obtain a very close approximation of the original functions. Indeed, among the many orthonormal bases one could envision, FPC has the advantage of being the most parsimonious allowing one to reconstruct the response curves using fewer coefficients than any other orthonormal basis. In scenarios investigated by simulation in Section 4 k = 5 is sufficient to capture more than the 99% of the -oono This produces an almost perfect approximation of (17) while fully preserving f fgen efficiency.
3.3 Solution path implementation
To evaluate different values of the penalty parameter '1 we implement an efficient solution path search. We compute the solution for decreasing sequence sequence of λ starting from Xmar (xx; TYll which selects O active features. When we move to the next Y value, we use the solution obtained the previous value for initialization (warm start). The two consecutive solutions tend to be close. and fgen converges in very few iterations usually just one. We also allow the user to specify a maximum number of selected features; when this number is reached the path search is stopped, further reducing computation.
To guide the choice of (2) we propose two quantitative criteria: k-fold Cross Validation (cv) and an Extended Bayesian Information Criterion (E-bic) (Chen and Chen 2012). which modifies the standard BIC also include the number of features p. In symbols, we have
e-bic B) = k log (rss(B)/(mkk)) + (kv) (log(nk) + log /n
(19)
where rss rss(B) is the residual sum of squares associated with the solution B, and v are the Group Elastic Net degrees of freedom. From Tibshirani al. 2012))
tr(X, (X3X, >2I,)) -xx))))x3) Note that cv can be very computationally expensive because it requires to run fgen multiple times for each value of Åι and x2 under consideration. In contrast, e-bic can be computed directly from the original solution. Before evaluating both criteria, we de-bias the fgen estimates following the approach suggested by Belloni et al. (2014); Zhao al. (2077). First, we fgen, then, we fit a standard least squares on the selected features. In the next section, following standard practice in the literature e.g.. Friedman al. (2010); Pedregoss et al. (201T)-
rewrite \₁ and X2 as A = Cy X and x2- = (1-a)cxAthar Amaa with cx E (0, 1] and a E (0, ex determines the reduction with respect of Xmaz or controls the relative weight of the two penalties.
4 Simulation study and INSIGHT data
In this section we use synthetic data to illustrate the computational efficiency of fgen, and apply our new method to Genome Wide Association Study (GWAS) on childhood obesity. In the simulations, we benchmark fgen against the two best Group Elastic Net solvers we found in the literature: the python package sklearn and R package glmnet, which is written in fortran. Other functional-od -onsclll feature selection methods, such as the ones proposed by Barber et al. (0017); Parodi et al. (2018); MMrshani and Reimherr (2019), have computational burden more than two order of magnitudes larger than fgen and could not complete instances with p > 10
4.1 Simulation results
We generate synthetic data as follows. The entries of the design matrix X E Raxp are each drawn independently from a standard normal distribution. The response curves are created as y XB + e. B contains Po non-zero curves. These and the errors c are generated from a O mean Gaussian process with a Matérn covariance function (Cressie and Huang] 1999) of the form
C(t, s) = w2 (f(m)2--1)~ -1 (10- (20) s1) " K, (2)
(20)
where K2 is a modified Bessel function. In particular, we set the point-wise variance 2² and the range | = 0.25 (this determines how fast the curves dependency decays). The smooth parameter
7

---

Table I: b and report CPU time in seconds for fgen, sklearn and glmnet, respectively. For fgen we also report the number of iterations parenthesis. is the number of selected features, I is the range parameter of the Matern process used to generate the coefficients.
a = 0.8,
7
0.25
n=500
n= 1000
n=50000
pi Po
Cx
a
r
b
a
C
0.8
2
0.1(3)
0.5
0.3
(.1
0.9
0.6
0.5
9.1
3.2
5
0.4
10
0.2(4)
0.5
0.3
4
0.3(3)
1.1
0.6
1.1(3)
10.9
3.1
2(10"); 10
21
0.4(4)
0.3(4)
1.2(3)
10.4
3,4
0.8
2
0.2
0.9
0.5
0.2(2)
1.9
1.2
1.3(2)
18.8
5.4
10
0.4
10
(43
0.8
0.5
0.4(3)
2.1
1.0
6
2.0(3)
22.3
4.9
0.2
21
(.9(4)
0.9
0.5
8
0.5(4)
2.2
1.2
10
1.9(3)
23.2
5.3
0.8
(23(2)
4.5
1.5
0.7(2)
9.9
3.3
2.6(2)
97.9
17.2
5
0.6
29
0.6(2)
4.7
1.6
17
0.7(2)
9.6
3.3
12
3.0(2)
85.7
19.8
10"; 10
0.4
486
(3.7(3)
4.8
0.2
53
1.4(2)
10.2
3.1
42
4.1(2)
99.7
17.8
0.8
0.6(2)
10.2
2.6
1.4(2)
22.6
5.2
3.5
175.2
27.9
10
29
(29(2)
10
2.4
17
(.5
19.6
4.8
3.8(2)
183.7
27.6
0.4
486
79.2(3)
9.1
2.4
53
1.8(2)
20.2
5.3
42
6.7(2)
190.7
28.0
0.5,
= 0.25
n=500
n= 1000
n=500
p. po
k
Cx
a
b
c
r
a
b
c
a
b
C
0.8
2
0.1
1.0
0.6
0.2(2)
1.8
1.3
1
0.9
18.3
5.2
2(10");
(32(3)
(.2(2)
0.9
7
(.9
20.6
0.2
73
1.1(3)
0.9
0.5
0.4(3)
2.0
1.3
(39
22.7
5.2
0.8,
= 0.10
n=500
n= 1000
n= 5000
p: Po
CA
r
b
r
a
C
0.8
3
0.2(3)
0.9
0.6
0.2
1.8
1.3
0.9
19.7
5.1
2(10*);
10
0.4
(4)
0.6
0.4
2.6
10
1.9
19.8
5.2
0.2
1.5(4)
1.2
0.5(4)
2.4
1.2
2.0(3)
22.2
4.0
responses
errors
coefficients
Lup estimated
T
49
?
Figure 1: Curves related to a simulation scenario described in Table ( row, 500). The left and center panels display sample of 10 response curves and IO linear model errors terms, respectively. The right panel depicts the true non-zeros coefficients curves (dashed lines) and the fgen de-biased estimates (solid lines).
equal 3.5 for and to 1.5 for c, i.e the errors are rougher than the coefficients. Each curve is sampled at 1000 evenly spaced points between 0 and 1. Figure shows instances of the response and error curves, and e, and the true non-zero coefficient curves in along with their de-bia estimates produced by fgen (the underlying simulation parameters are those in Table D second
row, n 500). all scenarios, fgen is rum with both the tolerances in (T6) set to 10 (we set the same tolerance for sklearn and glmnet) and and increase by a factor of every iteration. we start from smaller values of o, the algorithm needs more iterations to converge, while if o is too large, fgen does not converge to the optimal solution. We set D cyXmez and X2 = (1 a)cxA>m where cx E (0. 1), a E (0, 1). and mar max, [Xx, (X,) TYI. Note that for glmnet and sklearn need divide A mar by since both solvers divide the least squares loss in ( by the number of observations.
(2) set to 0.2. We start from の0 po/p
Table I reports CPU times for fgen, sklearn and glmnet under different simulation settings. fgen the fastest solver in almost every instance. When both n and p are large and the solution is sparse, fgen is approximately 6 times faster than glmnet and more than 30 times faster than sklearn. Note

---

response
10-fold
e-bic
0.18
==
9
m0
0.1.4
0.00608
D =0
7
0.10
10
0.00602
0.06
10
20
30
0.99
0.9
0.8
0.7 0.99 0.9
0.8
0.7
weeks
C)
Cy
Figure 2: Plots related to the INSIGHT study. The left panel displays the growth curves. The center and right panels depict values of the 10-fold co Mean Squared Error and the e-bic, respectively, against CA These are obtained from fgen run with 3 different values of a; 0.8 (green line), 0.6 (red line), 0.4 (blue line). fgen estimates are de-biased prior computing both criteria.
that the super-linear convergence rate allows fgen to converge in very few iterations (no more than 4 in all cases). The CPU time increases with k for all solvers. However, A: = 5 already captures more than the 99% of the -orm in all the scenarios considered. we decrease a from 0.8 to 0.5, fgen need even fewer iterations to converge, increasing its computational gain with respect to the competitors. Considering rougher coefficients (created with a Matern process with range parameter I 1= 0.1) does not affect the relative performance of the algorithms. The instance with an active set of 486 features is the only one where fgen performs worse than its competitors. As expected, in the presence of non-sparse solutions fgen looses some of its efficiency. However, tune the penalty parameters in practice, one evaluates a sequence of сΛ values starting from very sparse solutions. In the first steps of the solution path, fgen exploits sparsity and is very efficient. In the following steps, still converges very quickly thanks to the warm-slarl approach described in Section 33. In Supplemental Table B.2 compare the solution path computing time. fgen outperforms the other solvers in every scenario, being approximately 2 times faster than glmnet and from 10 to more than 30 times faster than sklearn. Finally, gauge uncertainty in our CPU time evaluations, we replicated a subset of the instances explored in Table independent times. Mean CPU times and standard errors over such replicates are reported in Supplemental Table B.I Results agree with those obtained considering just one replication. Furthermore, one can notice that fgen has also smaller variability CPU time when n 5000.
Taking into account that glmnet (written in fortran) and sklearn are highly optimized packages, the results above provide strong evidence in support of method. We also tracked prediction performance for all methods and in all simulation settings considered, but we do not report them here since all three solvers solve the same convex minimization problem and therefore converged the same solution in all settings.
4.2 INSIGHT study
Here, apply fgen to data from the Intervention Nurses Start Infants Growing Healthy Trajec- tories (INSIGHT) study (Paul et al.) 2014).. In particular, we focus on data collected to investigate genetic variants that may affect the risk of childhood obesity. As the prevalence of obesity increases also among children, examining possible causes and risk factors has become an essential public health concern. INSIGHT provides genome-wide Single Nucleotide Polymorphisms (SNPs) information for a cohort of very young children, along with longitudinal information on their growth. Selecting SNPs that may affect growth is thus GWAS Genome-Win Association Study) where the outcome is a growth curve. In recent years, many GWASs have identified SNPs strongly associated with obesity phenotypes (Locke et al. 2015). Before proceeding with the analysis, we point out that due high feature collinearity, low -ggallto-noise ratio, and ultra-high dimensionality. GWAS data are very hard examine and users should be very careful in interpreting results; e.g., selected SNPs may just be proxies for other causal SNPs in their vicinity. While being well aware of all its complexities and potential pitfalls, our main aim presenting GWAS analysis is to show the efficiency and broad applicability of fgen. Our functional outcome captures the evolution of weight/heighl ratios (Daniels et al.] 2015) measured at birth and at 4, 16, 28, and weeks for total of n 210 children. The growth curves shown in the left panel of Figure 2 are fitted as in Craig et al. (2019) using

---

Principal Analysis by Conditional Estimation (Chen et al., 2017). Building smooth curve for each child allows one to capture information along the entire time domain and at the same time de-noise and mitigate the effect of outlyng//aoomalon are available upon request at affibGaaP using the access number (as000149.vvl The growth curves are based on privacy protected data and cannot be made publicly available.
raw measurements. The SNPs collected in INSIGHT
Craig al. 2019) used flame (Parodi et al. 2018) solve the function-on-seit feature selection problem. To do so, they had to reduce the analysis from p 342325 to 10000 SNPs with various preliminary screening steps. The computational efficiency of fgen allows us to inspect all the 342325 SNPs simultaneously. The center and left panels of Figure 2 display 10-fold cV and e-bic for different values of a and ex. Both criteria identify just one (cx- 0.99) dominant SNP, rs79187646. Without drawing any strong domain conclusion, we remark that the selected SNP appears to be relevant in the literature. Notably, this the same SNP selected by [Boschi et al. (2020), where the same SNP data were associated to BMI at age 3 a scalar response. Also notably, this dominant SNP was not among the SNPs identified by Craig et al. (2019). Based on the U.S. National Library of Medicine rs79187646 is located in NTM, a well known gene. According to the NHGRI-EBI GWAS Catalog many GWASs, including two recent studies (Kichaev et al. 2019 Pulit et al. 2019), have connected NTM to body mass, food addiction, intake of sweet substances and other obesity-related traits.
5 Conclusions
In this article we proposed new Function Group Elastic Net method (fgen) to solve the function-on- scalar feature selection problem. Our proposal starts with the development of a novel, highly-efficient SsNAL algorithm solve the Group Elastic Net which is then extended to the function-on-scall regression framework using Functional Principal Components representation. Though we could rely on critical prior results (Tomioka and Sugiyama 2009: Li al., 2018 Boschi et al. 2020), in order to integrate the group structure into SsNAL, we had to tackle more complex mathematical operators and redefine the theoretical foundation. Our simulations show a substantial reduction in CPU time with respect to the best existing Group Elastic Net solvers. Finally, we applied fgen to a GWAS study detecting SNP that may affect obesity risk in children.
The current version of fgen is limited to the case where each group has the same size k and to the funetional-on-seallr feature selection problem. In the future, we plan to further extend work investigating more complex optimization problems (e.g., allowing each group to have a different size) and functional regression frameworks. In particular, we aim to adapt our methodology the function-on-functee feature selection scenario, where both the response and the predictors can be represented functional curves.
Acknowledgments and Disclosure of Funding
We thank Kateryna Makova and her laboratory Penn State for access to the INSIGHT data, Ana Kenney for help with the data and useful discussions. The work of Tobia Boschi and Francesca Chiaromonte was partially supported by the Huck Institutes of the Life Sciences at Penn State, the work of Matthew Reimher was partially supported by the Grant NSF SES-1853209.
References
Barber, R. F., M. Reimherr, T. Schill, et al. (2017). The function-on-ssallrt lasso with applications to longitudinal gwas. Electronic Journal of Statistics 11(1), 1351-1389,
Beck, A. (2017). First-order methods in optimization, Volume 25. SIAM.
Beck, A. and M. Teboulle (2009). A fast iterative serinkage-thresooldin algorithm for linear inverse problems. SIAM journal on imaging sciences 2(1), 183-202.
Belloni, A., V. Chernozhukov, and C. Hansen (2014). Inference on treatment effects after selection among nigh-dimensional controls. The Review of Economic Studies 81(,), 608-650
Bonnefoy, A., V. Emiya, L. Ralaivola, and R. Gribonval (2015). Dynamic screening: Accelerating first-order algorithms for the lasso and gassp- IEEE Transactions on Signal Processing 63(19), 5121-5132.
10

---

Boschi, T., M. Reimherr, and F. Chiaromonte (2020). An efficient semi-smooth newton augmented lagrangian method for elastic net. arXiv preprint 20000.0000..3970.
Boyd, S. and L. Vandenberghe (2004). Convex optimization. Cambridge university press.
Breheny, P. and J. Huang (2015). Group descent algorithms for nonconvex penalized linear and logistic regression models with grouped predictors. Statisties and computing 25(2), 173-187.
Chen, and Z. Chen (2012). Extended bic for small-n-large-p sparse glm. Statistica Sinica, 555-574.
Chen, K., X. Zhang, A. Petersen, and H.-G. Müller (2017). Quantifying infiite-diimensione data: Functional data analysis in action. Statistics in Biosciences 9(2), 582-604.
Chen, X., Q. Lin, S. Kim, J. G. Carbonell, and E. P. Xing (2010). An efficient proximal gradient method for general structured sparse learning. stat 1050.
Chen, J. Goldsmith, and R. T. Ogden (2016). Variable selection function-on-scalar regression. Stat 5(1), 88-101.
Chiou, J.-M.. H.-G. Müller, and J.-L. Wang (2004). Functional response models. Statistica Sinica, 675-693.
Clarke, F. H. (1990). Optimization and nonsmooth analyvis, Volume 5. Siam.
Craig, S. A. M. Kenney, J. Lin, M. Paul, L. L. Birch, J. Savage, M. E. Marini, F. Chiaromonte, M. L. Reimherr, and K. D. Makova (2019). Polygenic risk score based on weight gain trajectories is strong predictor of childhood obesity. bioRxiv, 606277.
Cremona, M. A., H. Xu, K. D. Makova, M. Reimherr, Chiaromonte, and P. Madrigal (2019). Functional data analysis for computational biology. Bioinformatics (Oxford, England) 35(17), 3211.
Cressie, N. and H.-C. Huang (1999). Classes of nonseparable, spatio-temporal stationary covariance functions. Journal of the American Statistical Association 94(448), 1330-1339,
Daniels. S. R., S. G. Hassink, et al. (2015). The role of the pediatrician in primary prevention of obesity. Pediatrics 136(1), e275-e292.
Deng, W., W. Yin, and Y. Zhang (2013). Group sparse optimization by alternating direction method. In Wavelets and Sparsity XV, Volume 8858, pp. 88580R. International Society for Optics and Photonics.
Deng, Z. and A. M.-C. So (2019), An efficient augmented lagrangian based method for constrained lasso. urXiv preprint 1rrXiv1903..00
Dünner, C., S. Forte, M. Takáč, and M. Jaggi (2016). Primal-dual rates and certificates. arXiv preprint
Fan, Y., G. M. James, P. Radchenko, et al. (2015). Functional additive regression. The Annals of Statistics 43(5), 2296-2325.
Fan. Z. and M. Reimherr (2016). High- digh-dimmnsional function-on function-on-scale regression via the alternating direction method of multipliers. 2016 3rd International Conference on Information Science and Control Engineering (ICISCE), pp. 397-399. IEFF.
Fenchel, W. (1949). On conjugate convex functions. Canadian Journal of Mathematics I(1), 73-77.
Friedman, J., T. Hastie, and R. Tibshirani (2010). Regularization paths for generalized linear models via coordinate descent. Journal of statistical software 33(1), 1.
Gertheiss, J., A. Maity, and A.-M. Staicu (2013). Variable selection in generalized functional linear models. Stat 2(1), 86-101
Hall, P. and M. Hosseini-Nasab (2006). On properties of functional principal components analysis.
Journal of the Royal Statistical Society: Series B (Statistical Methodology) 68(1), I), 109-126.

---

Horváth, L. and P. Kokoszka (2012). Inference for functional data with applications, Volume 200. Springer Science Business Media.
James, G. M., T. J. Hastie, and C. A. Sugar (2000). Principal component models for sparse functional data. Biometrika 87(3), 587-602.
Johnson, R. A., D. W. Wichern, et al. (2014). Applied multivariate statistical analysis, Volume 6. Pearson London, UK:
Katrutsa, A. and V. Strijov (2015). Stress test procedure for feature selection algorithms. Chemomet- rics and Intelligent Laboratory Systems 142, 172-183.
Kichaev, G., G. Bhatia, P.-R. Loh, S. Gazal, K. Burch, M. K. Freund, A. Schoech, B. Pasaniuc, and A. L. Price (2019). Leveraging polygenic functional enrichment to improve gwas power. The American Journal of Human Genetics 104(1), 65-75.
Kokoszka, P. and M. Reimherr (2017). Introduction functional data analysis. CRC Press.
Krawezyk, B. and A. Cano (2018). Online ensemble learning with abstaining classifiers for drifting and noisy data streams. Applied Soft Computing 68, 677-692.
Li, X., D. Sun, and K.-C. Toh (2018). A highly efficient semismooth newton augmented lagrangian method for solving lasso problems. SIAM Journal Optimization 28(1), 433-458.
Locke, A. E., B. Kahali, S. Berndt, A. E. Justice, T. H. Pers, F. Day, C. Powell, S. Vedantam, M. L. Buchkovich, J. Yang, et al. (2015). Genetic studies of body mass index yield new insights for obesity biology. Nature 518(7538), 197.
Luque, F. J. (1984). Asymptotic convergence analysis of the proximal point algorithm. SIAM Journal on Control and Optimization 22(2), 277-293
Matsui, H. and S. Konishi (2011). Variable selection for functional regression models via the I1 regularization. Computational Statistics & Data Analysis 55(12), 3304-3310.
Mirshani, A. and M. Reimherr (2019). Adaptive function-on-scalar regression with a smoothing elastic net. arXiv preprint arXiv: 1905.09881
Murdoch, W. J., C. Singh, K. Kumbier, R. Abbasi-Asl, and B. Yu (2019). Definitions, methods, and applications in interpretable machine learning. Proceedings of the National Academy of Sciences 16(44), 22071-22080.
Parikh, N., S. Boyd, al. (2014). Proximal algorithms. Foundations and Trends in Optimiza- tion 1(3), 127-239.
Parodi, A., M. Reimher,, et al. (2018). Simultancous variable selection and smoothing for high- dimensional function-on-scalar regression. Electronic Journal of Statistics /2(2), 4602-4639.
Paul, I. M., S. Williams, S. Anzman-Frasca, S. Beiler, K. D. Makova, M. E. Marini, L. B. Hess, S. E. Rzucidlo, N. Verdiglione, J. A. Mindell, et al. (2014). The intervention nurses start infants growing healthy trajectories (insight) study. BMC pediatrics 14(1), 184.
Pedregosa, F., G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Pretten- hofer, R. Weiss, V. Dubourg, et al. (201 1). Scikit-learn: Machine learning in python. the Journal of machine Learning research 12, 2825-2823
Polyak, B. T. (1969). The conjugate gradient method in extremal problems. USSR Computational Mathematics and Mathematical Physics 9(4,, 94-
Pulit, S. L.. C. Stoneman, A. P. Morris, A. R. Wood, C. A. Glastonbury, J. Tyrrell, L. Yengo, T. Ferreira, F. Marouli, Ji, et al. (2019). Meta-analysis genome-wide association studies for body fat distribution in 694 649 individuals of european ancestry. Human molecular genetics 28(1), 166-174.
Rockafellar, R. T. (1976a). Augmented lagrangians and applications of the proximal point algorithm in convex programming. Mathematics of operations research I(2), 97-116.
12

---

Rockafellar, R. T. (1976b). Monotone operators and the proximal point algorithm. SIAM journal on control and optimization 14(5), 877-898.
Sørensen, H., Goldsmith, and L. M. Sangalli (2013). An introduction with medical applications to functional data analysis. Statistics in medicine 32(30), 1222-52400
Tibshirani, R. J., J. Taylor, et al. (2012). Degrees of freedom in lasso problems. The Annals of Statistics 40(2), 1198-1232.
Tomioka, and M. Sugiyama (2009). Dual-augmented lagrangian method for efficient sparse reconstruction. IEEE Signal Processing Letters 16(12), 1067-1070.
Tomioka, R., T. Suzuki, and M. Sugiyama (2011). Super-linear convergence of dual augmented lagrangian algorithm for sparsity regularized estimation. Journal of Machine Learning Re- search 12(5).
Ullah, S. and C. F. Finch (2013). Applications of functional data analysis: systematic review. BMC medical research methodology I3(1), 1-12.
Van Loan, C. F. and G. H. Golub (1983). Matrix computations. Johns Hopkins University Press.
Yu, B. and K. Kumbier (2020). Sciences 117(8). 3920-3929.
Veridical data science. Proceedings of the National Academy of
Yuan, M. and Y. Lin (2006). Model selection and estimation in regression with grouped variables. Journal of the Royal Statistical Society: Series B (Statistical Methodology) 68 1), 49-67.
Zhang, Y. and Q. Yang (2018). An overview of multi-tast learning. National Science Review 5(1), 40-43.
Zhao, S., A. Shojaie, and D. Witten (2017). In defense of the indefensible: very naive approach to high-dimensional inference. arXiv preprint arXiv:1705.0 1755.0553.
7hu, Y. (2017). An augmented admm algorithm with application to the generalized lasso problem. Journal of Computational and Graphical Statistics 26(1), 195-204.
Zou, and Hastie (2005). Regularization and variable selection via the elastic net. Journal of the royal statistical society: series B (statistical methodology) 67(2), 301-320.
NeurIPS paper checklist
1. For all authors...
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope? [Yes]
(b) Did you describe the limitations of your work? [Ye]] We briefly describe some
limitations our methodology at the end of Section 5
(c) Did you discuss any potential negative societal impacts of your work? [No] We do not the potential for adverse societal impacts
(d) Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yes]
2. If you are including theoretical results...
(a) Did you state the full of assumptions all theoretical results? [Yes]
(b) Did you include complete proofs of all theoretical results? [Yes] All proofs are included
in Supplemental Section A
3. If you ran experiments..
(a) Did you include the code, data, and instructions needed reproduce the main ex- perimental results (either in the supplemental material or URL)? [Yes] Code to reproduce the simulation results included the Supplemental Material.
(b) Did you specify all the training details (e.g., data splits, hyperparameters, how they were chosen)? [Yes] See Section 4
13

---

(c) Did you report error bars (e.g., with respect to the random seed after running experi- ments multiple times)? [Yes] We rum subset of experiments multiple times computing mean and standard errors. These results are included in Supplemental Section
(d) Did you include the total amount of compute and the type of resources used (e.g.. type of GPUs, internal cluster, or cloud provider)? [Yes] See Supplemental Section
4. If you are using existing assets (e.g., code, data, models) or curating/releasing new assets...
(a) If your work uses existing assets, did you cite the creators? [Yes] We cited the creators of data from INSIGHT in Section 4.2
(b) Did you mention the license of the assets? [N/A]
(c) Did you include any assets either in the supplemental material URL? [Yes] The code to run fgen is included in the Supplementary Material and will be made available GitHub
(d) Did you discuss whether and how consent was obtained from people whose data you're using/curating? [No] The INSIGHT study operated under an appropriate IRB. We did not reference it or provide details here. We do thank the individuals who provided access to the data in the Acknowledgments Section.
(e) Did you discuss whether the data you using/curating contains personally identifiable information or offensive content? [Yes] The INSIGHT data we used for our application is privacy potected. We discuss how parts of it are available upon request through dbGaP, and parts of are not available, in Section 4.2
5. If you used crowdsourcing or conducted research with human subjects...
(a) Did you include the full text of instructions given participants and screenshots, if applicable? [N/A]
(b) Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals, applicable? [N/A]
(c) Did you include the estimated hourly wage paid to participants and the total amount spent participant compensation? [N/A]
14


Reviews:
[{"review_id": "y6XLTcWOu1", "reviewer": "Reviewer_2b9Y", "summary": "This paper adds to the various existing solvers for group elastic net regression. Here, a grouping of input variables is jointly penalized with both L1 and L2 penalties in order to yield an estimate of sparse penalized regression coefficients. The newly developed optimization algorithm — semi-smooth newton augmented langrangian — incorporates group structure and leverages sparsity to solve the optimization problem an order of magnitude faster than off-the-shelf methods from python sklearn and R glmnet. The method is further extended to function-on-scalar regression by leveraging functional principal components analysis. This extension allows for analyzing real data where the target is functional curves, in this case growth curves of childhood obesity development. In the INSIGHT study, one single nucleotide polymorphism was found to be significantly associated with childhood obesity, a locus occurring in a gene that has been previously associated with anthropometric-related phenotypes.", "questions": "", "limitations": "", "rating": 6, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The authors present a novel algorithm for fast optimization of group elastic net regression. This is further extended to function-on-scalar regression to increase the relevance of the method. The motivation is clear and the paper is well-organized and clearly written. The experiments are clearly laid out and easy to follow. The empirical results comparing computational efficiency are strong, and it seems the new methodology here could become a new baseline, given its comparison to already highly adopted packages. Theoretical results further confirm the efficiency of the proposed method. The applications to a childhood obesity GWAS provide a novel perspective on traditional GWAS framework, and while the result of confirming a previously discovered significant association is encouraging, there are some areas that could be improved to increase impact.\n\nFirst, GWAS is not typically considered with penalized regression due to identifiability issues. For this reason, it may be more applicable to consider this methodology for prediction rather than hypothesis testing, which is difficult to interpret in this setting. The signal-to-noise ratio in GWAS is often very low (as defined by the overall variance explained or variance explained per SNP), and it appears the current methodology might struggle in this setting as well. In particular, there is only a single SNP that’s identified by the proposed method, and the $c_\\lambda$ parameter is effectively at its maximum, so the model is being constrained as heavily as possible. It might be more fruitful to consider applications such as polygenic risk scores (where prediction is considered), or to look into gene expression experiments with measurements taken over time in order to increase the sample size. The sample size here (n=210) is bound to be underpowered to detect associations.\n", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "The authors have addressed limitations of their method. Societal impacts of the GWAS application should be discussed a bit more thoroughly considering this methodology paves the way for more function-on-scalar type GWAS to be investigated, and the current assumptions of the method should be made clear for those interested in hypothesis testing.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "4"}, "extra_text": {}}, {"review_id": "nFW5f0qqz8a", "reviewer": "Reviewer_FbcU", "summary": "The main contribution of the paper is two-fold.\n- First, the paper proposed a fast algorithm to learn the regression coefficients with the group elastic net regularization. Suppose that the dataset has k-dimensional response variables, and therefore the regression coefficients B has the form $B\\in R^{p\\times k}$ (p: number of features). Then the group elastic net regularization enhances the sparsity of features for all k response variables (i.e., some rows in B are to be zero).\n- Secondly, the paper proved that it can be extended to function-on-scalar setup, that is, we assume a regression problem where the response variable is a functional object of summing $k$ basis functions. (Therefore each element in $B\\in R^{p\\times k}$ is a functional object.)\n", "questions": "", "limitations": "", "rating": 6, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The proposed method is to efficiently learn a sparse function-on-scalar model, and (although it is not specifically stated in the paper) a sparse multitask model. It is theoretically proved that, even if the number of features is huge, the main computational cost depends on the number of only a part of features ($|\\mathcal{J}|$ in the paper), which is of value.\n\nSuggestions:\n\n-   The improvements from the conventional method is better to be stated in Section 1. For example, the paper stated three conventional methods in line 152 and their limitation for our problem setup, but the overview of such discussions is better to be appeared in Section 1 to emphasize how the proposed method achieved the improvements.\n-   The paper mainly consider a fast method to learn a sparse function-on-scalar model, but it solves a sparse multitask model (multiple response variables) as a building block. So, in order to clarify the advantage of the proposed method, it is better to discuss the relationships to conventional learning methods for sparse multitask model. For example, the following paper discusses the group lasso (instead of group elastic net) regularization for multitask learning:\n    -   Karim Lounici, Massimiliano Pontil, Alexandre B. Tsybakov and Sara A. van de Geer, \"Taking Advantage of Sparsity in Multi-Task Learning\". COLT2009.\n\nQuestions:\n\n-   Suppose that we consider a model with ordinary response variable (i.e., non-multitask, non-function-on-scalar). Then, does the proposed method work efficiently compared to conventional elastic net learning methods? In other words, how the multitask setup and/or function-on-scalar setup is employed for efficient computation?\n-   In Section 3.1, the definition of $\\hat{X}_\\mathcal{J}$ is not so obvious: since $T\\in R^{p\\times k}$, $\\mathcal{J}$ must be $\\mathcal{J}\\subseteq\\{1, 2, \\dots, p\\}$ and therefore inappropriate as a subscript for $\\hat{X}\\in R^{nk\\times pk}$. Does it mean that $\\hat{X}_\\mathcal{J} = \\hat{\\mathcal{X}}$ where $\\mathcal{X} = X_\\mathcal{J}$?\n", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "As stated in the main review, the proposed method can largely reduce the learning cost both theoretically and experimentally. So it seems to be of use.\n\nThe limitation of the method is discussed in Section 5. As possible extensions of the method, the authors consider (1) different number of coefficients in a group, and (2) use for other functional models than function-on-scalar model.\n", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "3"}, "extra_text": {}}, {"review_id": "7ZGUFloZ_6l", "reviewer": "Reviewer_ZQkA", "summary": "This paper proposed a novel optimization algorithm for the Group Elastic Net problem, with additional application to the function-on-scalar regression problem. By exploiting the sparsity structure of the Augmented Lagrangian, the new algorithm significantly reduces the computational load. Experiments also showed supreme results, comparing with multi-task elastic net procedures in sklearn and glmnet. The paper is well written and easy to follow.", "questions": "", "limitations": "", "rating": 7, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "The following contributions from the paper are quite novel:\n* Exploiting the analytical form of the Fenchel conjugate of the group elastic net regularization.\n* Using semi-smooth Newton method to solve the dual variable in the subproblem.\n* Applying the algorithm to function-on-scalar regression problem.\n\nThis novel optimization algorithm proposed in this paper can be very valuable to the machine learning practitioners. ", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "Here are a few minor concerns I have for the paper:\n\n1. Table I only showed the speedup compared with sklearn and glmnet, but no comparison in solution accuracy is provided along with it. The speedup is meaningful only if the solution accuracy is comparable to other methods. A better way to compare different algorithms is to show the accuracy (or error) along with the time spent in the computation. \n\n2. The discussion on the convergence rate could be more rigorous. For example, the super-linear convergence rate is for the semi-smooth Newton subproblem, and I did not see any strict theorem or proof of convergence rate for the whole algorithm.\n\n3. The augmented Lagrange algorithm converges as the \\sigma variable goes to infinity. How is this increase in \\sigma scheduled in the algorithm? Will different scheduling scheme impact the convergence property?", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "4 hours"}, "extra_text": {}}, {"review_id": "3e0HmI4OtX3", "reviewer": "Reviewer_nVk8", "summary": "The paper proposes a new group elastic net algorithm solver for linear regression problem with multivariate response variable. The paper then claims its application on function-on-scalar.", "questions": "", "limitations": "", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "I am not an expert in optimization, but I feel that this paper is strikingly similar to [1], which aims to solve a similar elastic net problem but with univariate response variable. While the paper extends the univariate case to multivariate case, I think this extension is not very complicated. Most of methods and results have similar corresponding in [1]. At least for me, this paper seems trivial. I think that the authors should give a thorough related work section and explain what their main contributions are.\n\nBesides, I feel that the application on function-on-scalar regression is also a trivial application. The contribution of using FPCA is overly stated, since it seems very natural and is not the only choice. In fact, by using any function basis and doing dimension reduction via projection, we can always easily convert a function-on-scalar problem to a linear regression problem with multivariate response.\n\nThe paper is technically sound. The results are well proved.\n\nThe paper is also clearly written and well-organized.\n\nIn terms of significance, I think the paper has its own value in application perspective. Since the problem solved by the paper can appear in many places, an efficient implementation will be valuable. However, my main concern is about its limited originality as mentioned above. I would suggest this paper to be submitted to some statistics software focused conferences/journals rather than methodology/theory focused ones like NeurIPS. If the authors insist, it is better to explain its novelty in terms of methods or theoretical analysis.\n", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "The authors claim that they do not see any potential for adverse societal impacts. Since this paper is mainly theoretical, I agree with the authors' claim.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "4"}, "extra_text": {}}]

Author rebuttal:

Meta-review: Three reviewers indicated acceptance, only one reviewer had concerns about the novelty with respect to reference [1]. It turned out, however, that [1] is basically an unpublished pre-print of this paper, and in their rebuttal, the authors could convincingly show that even with respect to [1], the paper contains some important extensions. So I recommend acceptance of this paper.   
Decision: Accept (Poster)