In this work, we present an efficient method for estimating active subspaces using only random observations of gradient vectors. Our method is based on the bi-linear representation of low-rank gradient matrices with a novel initialization step for alternating minimization

Becker, Stephen

Pourkamali-Anaraki, Farhad

Undetermined

CU Scholar Institutional Repository

ESTIMATING ACTIVE SUBSPACES WITH RANDOMIZED GRADIENT SAMPLINGFarhad Pourkamali-Anaraki, Stephen BeckerSIAM Workshop on Dimension Reduction 2017July 9-10 · PittsburgSummaryIn this work, we present an e cient method for estimat-ing active subspaces using only random observations ofgradient vectors. Our method is based on the bi-linearrepresentation of low-rank gradient matrices with a novelinitialization step for alternating minimization.Active SubspacesIn modern computer simulations, scientists and engineersseek to study the relationships between high-dimensionalspaces of input parameters and quantities of interest. Dueto the large number of input parameters and high costof simulations, many methods have been proposed to re-duce the dimension of the input parameter space. Thesemethods often find small subsets or linear combinations ofthe input parameters that approximately preserve input-output relationships. This low-dimensional characteriza-tion of complex problems with hundreds or thousands ofinput parameters is a crucial tool for modern computersimulations.Active subspaces are powerful tools for identifying im-portant directions in the high-dimensional space of inputparameters [1, 2]. Let x 2 Rm be a vector of simulationinputs and assume that f(x) : Rm 7! R is the mappingbetween x and a quantity of interest. The active subspaceis defined by the top n < m eigenvectors of the followingm⇥m symmetric positive semidefinite matrixC =Zrf(x)rf(x)T ⇢(x)dx, (1)where rf(x) 2 Rm is the gradient vector and ⇢ is auser-specified probability density function. Consider theeigenvalue decomposition of C = W⇤WT , where ⇤ 2Rm⇥m is the diagonal matrix of eigenvalues, listed indecreasing order, andW 2 Rm⇥m containsm orthonormaleigenvectors. The matrix W can then be partitioned:W = [W1,W2]. The column space of W1 2 Rm⇥n isthe n-dimensional active subspace, where n is usuallychosen so the first n eigenvalues are much larger than theremaining m  n eigenvalues.In high-dimensional settings, computing the integral in(1) for constructing the matrix C is impractical. Moreover,in some applications, the gradient vector rf(x) may nothave a closed-form expression. In such cases, gradientscan be approximated by the first-order finite di↵erencewith m+ 1 function evaluationseTj rf(x) ⇡ (f(x+ hej)  f(x)) /h, j = 1, . . . ,m, (2)where ej is the j-th canonical basis vector in Rm andh > 0 is the finite di↵erence parameter.In [4], it is shown that the matrix C and its eigenpairscan be approximated using the following Monte Carlomethod. First, M samples x1, . . . ,xM 2 Rm are drawni.i.d. according to ⇢. Then, rfi := rf(xi) are estimatedvia (2) using M(m+ 1) function evaluations to formbC = 1MMXi=1rfirfTi = cWb⇤cWT . (3)The leading n eigenvectors of bC provide an accurate es-timate of W1 when M is su ciently large [4]. Theseeigenvectors are equivalent to computing the top left sin-gular vectors of the gradient matrixG := [rf1, . . . ,rfM ] 2 Rm⇥M . (4)In this work, we show that the low-rank structure ofthe gradient matrix G allows us to find accurate estimatesof active subspace using fewer function evaluations. Inparticular, we consider a scheme where only k entries ofeach gradient vector rfi 2 Rm are computed uniformlyat random, thus the total number of function evaluationsrequired would be M(k + 1). To estimate the activesubspace, G is written in a bi-linear form G = AB andthen alternating minimization [5] is used to find A and Bthat best fit the observed entries of gradient matrix. Tofurther improve the performance for small values of k, weuse an unbiased estimate of the left singular vectors of Gas the initial point for alternating minimization [6].Estimating Active Subspaces with Gradient SamplingLet us define the linear measurement operator L(·) asL(G) :=⇥RT1rf1, . . . ,RTMrfM⇤2 Rk⇥M , (5)1where each sampling matrix Ri 2 Rm⇥k contains k canon-ical basis vectors in Rm chosen uniformly at random with-out replacement. Given the incomplete observations andthe bi-linear parameterization ofG, our goal is to minimizeminA2Rm⇥n,B2Rn⇥MkL(G)  L(AB)kF . (6)This problem can be reformulated by using the vectoroperator and kronecker productkL(G)  L(AB)kF= ky  R(IM⇥M ⌦A)vec(B)k2= ky  R(BT ⌦ Im⇥m)vec(A)k2, (7)where y = Rvec(G) and R 2 RkM⇥mM contains allsampling matrices RTi , i = 1, . . . ,M , on its main diagonal.Thus, our method iteratively keep one of A,B fixed andoptimize over the other. Each subproblem is convex andcan be solved e cientlyvec(B) [R(IM⇥M ⌦A)]† y,vec(A) ⇥R(BT ⌦ Im⇥m)⇤†y. (8)In [3], estimation of active subspaces is considered in asimilar framework where eachRi 2 Rm⇥k has independentstandard Gaussian entries. As we see from the update rulein (8), our method is more e cient in terms of computationand memory. The proposed method in [3] must storeMmknonzero entries of matrix R, whereas our method requiresonly Mk nonzero entries to be stored. Similarly, ourmethod reduces the cost of matrix-matrix multiplications.We expect our proposed sampling method to perform aswell as Gaussian samples when the active subspaces areincoherent with the standard basis.Another contribution of our work is the novel initializa-tion step based on [6]. The initial iterate for alternatingminimization is preferred to be chosen based on a goodestimate of A, rather than random initializations, to guar-antee the convergence [5]. The recent work [6] presentsan unbiased estimator for the matrix bCe⌃ := eC  ⌘ diag(eC), ⌘ = m  km  1 , (9)where diag(eC) represents the matrix formed by zeroingall but the diagonal elements of eC, which is defined aseC := m(m  1)k(k   1)1MMXi=1(RiRTi rfi)(RiRTi rfi)T . (10)Numerical ExperimentsLet H 2 Rm⇥m be symmetric positive semidefinite andf(x) = 12xTHx, defined on the domain x 2 [ 1, 1]m witha uniform density ⇢. Thus, the gradient is rf(x) = Hx.The eigenvalues of C in are the eigenvalues of H, squaredand divided by 3. Moreover, the eigenvectors of C andH are identical. The matrix H is constructed so that itseigenvalues decay at a slow rate, except for a large gapbetween the fifth and sixth eigenvalues. We set parametersm = 100, M = 2000, and n = 5. The subspace estimationerror is defined as E := kcW1cWT1   fW1fWT1 k2, where fW1is the active subspace estimate using incomplete gradients.In Fig. 1, the mean estimation error over 100 trials isreported for various values of measurements k. Alternatingminimization with 20 iterations is used in two cases: (1)our proposed initial point based on (9), and (2) randominitialization based on a Gaussian matrix.10 20 30 40 50 60 70 80Number of measurements k10-210-1100ActivesubspaceerrorAltMin (random initial)AltMin (our initial)Figure 1: Active subspace estimation error for varyingnumber of measurements k and fixed dimension m = 100.A good initialization point becomes more important forsmall values of k, which are crucial for large-scale problems.For example, at k = 20, the mean estimation errors for ourproposed initial point and random initialization are 0.08and 0.18, respectively. Thus, our initialization procedurereduces the error by almost a factor of 2 in this case.References[1] P. Constantine. Active subspaces: Emerging ideas for dimensionreduction in parameter studies. SIAM, 2015.[2] P. Constantine, E. Dow, and Q. Wang. Active subspace methodsin theory and practice: applications to kriging surfaces. SIAMJournal on Scientific Computing, 36(4):A1500–A1524, 2014.[3] P. Constantine, A. Eftekhari, and M. Wakin. Computing activesubspaces e ciently with gradient sketching. In CAMSAP,pages 353–356, 2015.[4] P. Constantine and D. Gleich. Computing active subspaces withMonte Carlo. arXiv preprint arXiv:1408.0545, 2014.[5] P. Jain, P. Netrapalli, and S. Sanghavi. Low-rank matrix com-pletion using alternating minimization. In ACM symposium onTheory of computing, pages 665–674, 2013.[6] F. Pourkamali-Anaraki and S. Becker. Preconditioned datasparsification for big data with applications to PCA and K-means. IEEE Transactions on Information Theory, 2017.2

ESTIMATING ACTIVE SUBSPACES WITH RANDOMIZED GRADIENT SAMPLING

https://scholar.colorado.edu/downloads/pg15bf59j

ESTIMATING ACTIVE SUBSPACES WITH RANDOMIZED GRADIENT SAMPLING

Abstract

Similar works

Full text

Available Versions

CU Scholar Institutional Repository