Small Superposition Dimension and Active Set Construction for Multivariate Integration Under Modest Error Demand
aa r X i v : . [ m a t h . NA ] M a r Small Superposition Dimensionand Active Set Constructionfor Multivariate IntegrationUnder Modest Error Demand
A. D. Gilbert and G. W. WasilkowskiOctober 30, 2018
Abstract:
Constructing active sets is a key part of the Multivariate DecompositionMethod. An algorithm for constructing optimal or quasi-optimal active sets is proposedin the paper. By numerical experiments, it is shown that the new method can provide setsthat are significantly smaller than the sets constructed by the already existing method.The experiments also show that the superposition dimension could surprisingly be verysmall, at most 3, when the error demand is not smaller than 10 − and the weights decaysufficiently fast. In this short paper, we consider approximating integrals with infinitely many variables.We focus on approximations with a modest error demand, aiming at problems in, e.g.,
Mathematical Finance and
Uncertainty Quantification , where the underlying stochasticprocess is not known and hence only rough approximations are needed. In our tests weuse ε = 10 − n for n = 1 , , { γ u } u ⊂ N + -weighted tensor product Banachspaces F γ which allow for the decomposition f ( x ) = X u ⊂ N + , | u | < ∞ f u ( x ) . Here the summation is with respect to finite subsets u of positive integers and each f u depends only on the variables x j with j ∈ u . We also assume that the weights γ u have aproduct form.Integrals of such functions can be approximated by the Multivariate DecompositionMethod , which is a refined version of the
Changing Dimension Algorithm introduced in[4]. An essential part of those methods is the construction of an active set U ( ε ) of subsets1 such that the integral of P u / ∈ U ( ε ) f u can be neglected since it is bounded from above by ε (cid:13)(cid:13)(cid:13)(cid:13)(cid:13)(cid:13) X u / ∈ U ( ε ) f u (cid:13)(cid:13)(cid:13)(cid:13)(cid:13)(cid:13) F γ for all f ∈ F γ . In other words, it is enough to approximate integrals of the partial sum X u ∈ U ( ε ) f u . We would like to construct possibly small active sets and such that the largest cardinalityamong its elements u , d ( U ( ε )) := max u ∈ U ( ε ) | u | , is also small. This is because the partial sum P u ∈ U ( ε ) f u , that can be considered instead ofthe infinite sum P u ⊂ N + , | u | < ∞ , has a small number | U ( ε ) | of functions f u , each dependingon no more than d ( U ( ε )) variables.A specific construction of such sets (denoted by U PW ( ε )) was proposed in [7] and itwas shown there that the largest cardinality among all u ∈ U PW ( ε ) grows very slowly withdecreasing ε , d (cid:0) U PW ( ε ) (cid:1) = O (cid:18) ln(1 /ε )ln(ln(1 /ε )) (cid:19) as ε → . Moreover the size | U PW ( ε ) | grows polynomially in 1 /ε . However, the asymptotic constantsin the big- O notation were not investigated and, as we shall see, they could be very large.This is why in this paper we consider constructing possibly smallest active sets denotedby U opt ( ε ). As we will show by examples, the difference between the size of U opt ( ε ) and U PW ( ε ) could be very large. We also provide a construction of quasi-optimal sets, denotedby U q − opt ( ε ), which sometimes are only slightly larger than the optimal U opt ( ε ); however,their construction is less expensive.We are also interested in active sets with the smallest d ( U ( ε )). This leads to the follow-ing concept of ε -superposition dimension (or superposition dimension for short) definedby d sup ( ε ) := min { d ( U ( ε )) : U ( ε ) is an active set } . Since the optimal active sets in our experiments have very small d ( U opt ( ε )), this impliesthat the superposition dimension is also small.Note that our concept of the superposition dimension depends on the integrationproblem as well as the error demand ε . Hence it is in the same spirit as the definition of truncation dimension introduced recently in [3]. They are different from the definitions instatistical literature, see, e.g., [1, 5, 6, 9], where superposition and truncation dimensionsare defined based on ANOVA decompositions and without any relation to the integrationproblem or the error demand ε . Moreover, the dimensions from [1, 5, 9] depend on specificfunctions, whereas the dimensions in [3] and in this paper are defined in the worst casesense, i.e., are relevant to all functions from the space F γ .2lthough the algorithms for constructing U q − opt and U opt work for rather generalproblems and spaces, we applied them to the integration problem and for weighted spacesof functions with mixed first order partial derivatives bounded in L p norms for p ∈ [1 , ∞ ].Such spaces have often been considered (mostly for p = 2) when dealing with quasi-MonteCarlo methods .The results depend on how fast the weights converge to zero. In the experiments, weconsidered γ u = Y j ∈ u j − a for a = 2 , , . For p = 1, the construction from [7] is optimal, and it yields the following results: d sup (10 − ) = a = 4 , a = 3 , a = 2 , and | U PW (10 − ) | = a = 4 , a = 3 , a = 2 .d sup (10 − ) = a = 4 , a = 3 , a = 2 and | U PW (10 − ) | = a = 4 , a = 3 ,
22 for a = 2 .d sup (10 − ) = a = 4 , a = 3 , a = 2 and | U PW (10 − ) | =
10 for a = 4 ,
22 for a = 3 ,
114 for a = 2 . For p > U PW are no longer optimal; however they are not much worse than optimalsets when p is relatively close to 1. Moreover, for all the tests we have performed d ( U PW ( ε ))is very close to the superposition dimension. However the sizes (i.e., cardinalities) of U opt ( ε ) and U PW ( ε ) could be very different, especially for p = ∞ and/or small a .For instance, for p = 2 we have the following results. In the case of ε = 10 − d sup (10 − ) ≤ a = 4 , a = 3 , a = 2 , and | U opt (10 − ) | = a = 4 , a = 3 , a = 2 , whereas | U PW (10 − ) | = a = 4 , a = 3 ,
15 for a = 2 . For ε = 10 − d sup (10 − ) ≤ a = 4 , a = 3 , a = 2and | U opt (10 − ) | = a = 4 , a = 3 ,
30 for a = 2 , whereas | U PW (10 − ) | = a = 4 ,
18 for a = 3 ,
158 for a = 2 . ε = 10 − d sup (10 − ) ≤ a = 4 , a = 3 , a = 2 , and | U opt (10 − ) | = a = 4 ,
24 for a = 3 ,
255 for a = 2 , whereas | U PW (10 − ) | =
20 for a = 4 ,
70 for a = 3 , a = 2 . The results for a = 4 suggest that to achieve an error smaller than 10 − it is enough toapproximate f ( x ) = P u ⊂ N + , | u | < ∞ f u ( x ) by f ∅ + f { } ( x ) + · · · + f { } ( x ) + f { , } ( x , x ) + f { , } ( x , x ) + f { , } ( x , x ) . As for the quasi-optimal sets, they are the same for a = 4 and slightly larger for a = 3 , | U q − opt (10 − ) | = 6 for a = 2 , | U q − opt (10 − ) | = 32 for a = 2 , and | U q − opt (10 − ) | = (cid:26)
26 for a = 3 ,
261 for a = 2 . For p = ∞ we have d sup (10 − ) ≤ a = 4 , a = 3 , a = 2 , and | U opt (10 − ) | = a = 4 , a = 3 ,
33 for a = 2 , whereas | U PW (10 − ) | = a = 4 ,
21 for a = 3 , a = 2 . Now for ε = 10 − d sup (10 − ) ≤ a = 4 , a = 3 , a = 2 , and | U opt (10 − ) | = a = 4 ,
15 for a = 3 , a = 2 , whereas | U PW (10 − ) | =
21 for a = 4 ,
149 for a = 3120 ,
935 for a = 2 . For ε = 10 − d sup (10 − ) ≤ a = 4 , a = 3 , a = 2 , | U opt (10 − ) | =
15 for a = 4 ,
83 for a = 3 , ,
446 for a = 2 , whereas | U PW (10 − ) | = (cid:26)
72 for a = 4 ,
923 for a = 3 . For the tests above, the quasi-optimal active set was different from the correspondingoptimal active set in the following cases only: | U q − opt (10 − ) | = 38 for a = 2 , | U q − opt (10 − ) | = 1904 for a = 2 , and | U q − opt (10 − ) | = (cid:26)
92 for a = 3 , ,
159 for a = 2 . A collection of the active sets constructed above have been listed in full in the Ap-pendix.Our algorithms can also be used to construct the active sets where, instead of thestandard worst case error, the normalized worst case error is used. More precisely, for thenormalized worst case error we would like to have sets U norm ( ε ) such that the integral of P u / ∈ U norm ( ε ) f u is bounded by ε kSk k f k for all f ∈ F γ , where kSk is the norm of the integration operator. Since in our case kSk ≥
1, thecorresponding active sets U Xnorm ( ε ) are subsets of U X ( ε ) for X ∈ { opt , q − opt , PW } andcould be even smaller. We provide in this section basic concepts and definitions for special spaces of functionsthat are very often assumed in the literature, especially in the context of quasi-MonteCarlo methods. The presented algorithms can easily be modified to more general spaces. γ -Weighted Spaces We follow here [2]. For D = [0 , D = D N + be the set of sequences (points) x =[ x , x , . . . ] with x i ∈ D . Here N + is the set of positive integers and we will use u and v to denote finite subsets of N + . We will also use the following notation: For x ∈ D and u ,by [ x u ; u c ] we denote the point in D such that[ x u ; u c ] = [ y , y , . . . ] with y j = (cid:26) x j if j ∈ u , j / ∈ u . Next, f ( u ) = Y j ∈ u ∂∂x j f. p ∈ [1 , ∞ ] , let F γ ,p be the Banach space of functions defined on D with the following norm k f k F γ ,p = X u ⊂ N + , | u | < ∞ γ − p u k f ( u ) ([ · u ; u c ]) k pL p /p . Of course, for p = ∞ , k f k F γ ,p = sup u ⊂ N + , | u | < ∞ k f ( u ) ([ · u ; u c ]) k L ∞ γ u . We assume that the numbers γ u (called weights ) are of product form (see [8]) γ u = Y j ∈ u cj a for positive a and c. (1)In general choosing the weights (in our case choosing a and c ) for a specific integral orapplication is a difficult problem which we do not attempt to address here. We assumethat the parameters a and c are given with the problem.It was shown in [2] that any f ∈ F γ ,p admits a unique decomposition, called the anchored decomposition , f = X u ⊂ N + , | u | < ∞ f u with f u given by f u ( x ) = T u ( h u )( x ) := Z D | u | h u ( t ) Y j ∈ u ( x j − t j ) d t for some h u ∈ L p ( D | u | ) , where ( x j − t j ) is 1 if x j > t j and 0 otherwise. The functions f u belong to the followingBanach spaces F u F u = T u ( L p ) and k f u k F u = k f ( u ) u k L p . Of course, F ∅ is the space of constant functions with the absolute value as its norm. Thespaces F u for u = ∅ are anchored at 0 since f u ( x ) = 0 if there is j ∈ u with x j = 0. Thisis why f ( u ) ([ · u ; u c ]) = f ( u ) u and k f k F γ ,p = X u ⊂ N + , | u | < ∞ γ − p u k f u k pF u /p . The space F γ ,p contains in particular the following class of functions. Example 1
For a smooth function g : R → R and fast decaying numbers a , a , . . . ,consider f ( x ) = g ∞ X j =1 x j a j ! for x j ∈ D. (2)6hen f ( u ) ([ x u ; u c ]) = g ( | u | ) X j ∈ u x j a j ! Y j ∈ u a j . Hence f ∈ F γ ,p if the derivatives of g and the coefficients a j satisfy X u ⊂ N + , | u | < ∞ Q j ∈ u | a j | p γ p u Z D | u | (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) g ( | u | ) X j ∈ u x j a j !(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) p d x u /p < ∞ . Consider the following integration functional S : F γ ,p → R given by S ( f ) = lim s →∞ Z D s f ( x , . . . , x s , , . . . ,
0) d[ x , . . . , x s ] . Let p ∗ denote the conjugate of p , 1 p + 1 p ∗ = 1 . We assume that X u ⊂ N + , | u | < ∞ γ p ∗ u ( p ∗ + 1) −| u | /p ∗ < ∞ (3)since the left hand side of (3) is the norm of S , i.e., (3) is a necessary and sufficientcondition for continuity of S . Indeed, letting S u be the restriction of S to F u , we havethat k S u k F u = sup k f u k F u =1 S u ( f u ) = 1( p ∗ + 1) | u | /p ∗ which, with an application of H¨older’s inequality, yields (3). For product weights of theform (1), we have kSk = X u γ p ∗ u ( p ∗ + 1) −| u | ! /p ∗ = ∞ Y j =1 (cid:18) c p ∗ j a p ∗ ( p ∗ + 1) (cid:19) /p ∗ . Hence for product weights, (3) is equivalent to a > /p ∗ . For the remainder of the paperit is assumed that a > /p ∗ .A very important part of the Multivariate Decomposition Method ( MDM for short) isa construction of active sets U ( ε ), i.e., sets that satisfy (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) S X u / ∈ U ( ε ) f u (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ≤ ε (cid:13)(cid:13)(cid:13)(cid:13)(cid:13)(cid:13) X u / ∈ U ( ε ) f u (cid:13)(cid:13)(cid:13)(cid:13)(cid:13)(cid:13) F γ ,p for all f ∈ F γ ,p . (4)7he essence of (4) is that, when approximating S ( f ), it is enough to restrict the attentionto functions X u ∈ U ( ε ) f u , since any algorithm approximating P u ∈ U ( ε ) S u ( f u ) with the worst case error on L u ∈ U ( ε ) F u bounded by ε has its worst case error on the whole space F γ ,p bounded by2 /p ∗ ε. The factor of 2 /p ∗ is the result of applying H¨older’s inequality, see, e.g., [3]. Clearly, thereare many sets satisfying (4), and we would like to construct possibly small active sets. Definition 2
We say that an active set, denoted by U opt ( ε ) is optimal , if | U opt ( ε ) | = min {| U ( ε ) | : U ( ε ) satisfies (4) } . We also define the ε - superposition dimension as the smallest d ( U ( ε )) among all activesets, d sup ( ε ) := min { d ( U ( ε )) : U ( ε ) satisfies (4) } . U ( ε ) A construction of active sets was first proposed in [7]. The corresponding sets will bedenoted by U PW ( ε ). It was shown there that d ( U PW ( ε )) = O (cid:18) ln(1 /ε )ln(ln(1 /ε )) (cid:19) as ε → U PW ( ε ) is polynomial in 1 /ε . It is easy to see that for p = 1, U PW ( ε ) are optimal. However, as we shall see, their size might be too large for large valuesof p , especially for p = ∞ .Let U be a set of subsets u . Then (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) S X u / ∈ U f u !(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ≤ X u / ∈ U k f u k F u k S u k F u = X u / ∈ U k f u k F u γ u γ u k S u k F u ≤ (cid:13)(cid:13)(cid:13)(cid:13)(cid:13)X u / ∈ U f u (cid:13)(cid:13)(cid:13)(cid:13)(cid:13) F γ ,p X u / ∈ U γ p ∗ u k S u k p ∗ F u ! /p ∗ . Hence we are looking for U ( ε ) such that X u / ∈ U ( ε ) γ p ∗ u k S u k p ∗ F u /p ∗ ≤ ε. (5)Since H¨older’s inequality is sharp, (5) is equivalent to (4).For the sake of completeness, we recall the construction for p = 1, see [7].8 .1 Case p = 1 For p = 1, we have p ∗ = ∞ and k S u k F u = 1 for all u . Hence X u / ∈ U ( ε ) γ p ∗ u k S u k p ∗ F u /p ∗ = sup u / ∈ U ( ε ) γ u which for product weights reduces to sup u / ∈ U ( ε ) Q j ∈ u c j − a . Therefore U PW ( ε ) = ( u : Y j ∈ u cj a > ε ) . (6)It is easy to see that U PW ( ε ) is the smallest set satisfying (4), i.e., it is a subset of any U ( ε ) satisfying (5).The examples of U PW for specific values of a and ε are presented in the Appendix. Forsimplicity we use c = 1 there. p > For p >
1, the conjugate p ∗ is finite and the construction of U ( ε ) is more complicated.We begin by recalling the construction of U PW ( ε ) in [7]. To simplify the notation, let γ u = γ p ∗ u ( p ∗ + 1) | u | = (cid:18) c p ∗ p ∗ + 1 (cid:19) | u | Y j ∈ u j − a p ∗ . For given ε and p , a special threshold is computed and all u with γ u exceeding the thresholdare included in the active set. More precisely, for t ∈ (1 / ( ap ∗ ) ,
1) a threshold is given byThreshold( ε, t ) = ε p ∗ P u ⊂ N + , | u | < ∞ γ t u ! / (1 − t ) . Note that the interval (1 / ( ap ∗ ) ,
1) is non-empty by the assumption that a > /p ∗ in-troduced in Section 2.2. In our numerical experiments we approximated the sum of γ t u for t = i/
40 (39 ≥ i > / ( ap ∗ )) and selected the value which resulted in the largestThreshold( ε, t ). The approximations are calculated in a similar way as the computationof A s explained later (see (7)).Clearly, U PW ( ε ) contains a number of u ’s with the largest γ u ; however, the number ofthem could be much larger than needed.The collection of those u with the largest γ u that are necessary for (4) would resultin the optimal set U opt ( ε ). Since this optimal set is always a subset of U PW ( ε ), the goodproperty d ( U opt ( ε )) = O (cid:18) ln(1 /ε )ln(ln(1 /ε )) (cid:19) as ε → , is preserved. More precisely, let ( u j ) j ∈ N + be a sequence of all subsets u ordered so that γ u j ≥ γ u j +1 j = 1 , , . . . . U opt ( ε ) = { u , . . . , u k } with k = k ( ε ) such that kSk p ∗ − k X j =1 γ u j ≤ ε p ∗ < kSk p ∗ − k − X j =1 γ u j . The problem with this approach is that we do not know a priori the number k = k ( ε ) andordering a large number of γ u might be too expensive. This is why the numbers γ u willbe ordered on-line. Actually, we propose two ways of constructing active sets. The firstand simpler one produces what we call, quasi-optimal sets U q − opt ( ε ) and it uses a partialordering of γ u . The second one, uses ordering of γ u and produces optimal sets U opt ( ε ).However, as we will see the difference between both sets is very small; sometimes thesesets are equal.The numbers γ u have the following properties that are crucial for our construction ofquasi-optimal and optimal sets. Let ℓ be a given cardinality. In what follows we will write u = { u , . . . , u ℓ } , where u < · · · < u ℓ . The first property is: If u = { u , . . . , u ℓ } and v = { v , . . . , v ℓ } with v j ≥ u j for all j then γ u ≥ γ v . The other property is: For ℓ + 1 ≥ c /a , γ { u ,...,u ℓ } ≥ γ { u ,...,u ℓ ,u ℓ +1 } . We are ready to describe the constructions of active sets. First we need to approximate A = X u ⊂ N + , | u | < ∞ γ u from above and with the relative error significantly smaller than ε p ∗ . This can be done asfollows. For a large natural number sA = exp ln ∞ Y j = s +1 (cid:18) c/j a ) p ∗ p ∗ + 1 (cid:19)!! s Y j =1 (cid:18) c/j a ) p ∗ p ∗ + 1 (cid:19) ≤ exp c p ∗ p ∗ + 1 ∞ X j = s +1 j − a p ∗ ! s Y j =1 (cid:18) c/j a ) p ∗ p ∗ + 1 (cid:19) ≤ exp (cid:18) c p ∗ p ∗ + 1 Z ∞ s +1 / x − a p ∗ d x (cid:19) s Y j =1 (cid:18) c/j a ) p ∗ p ∗ + 1 (cid:19)
10 exp (cid:18) c p ∗ ( p ∗ + 1) ( a p ∗ −
1) ( s + 1 / a p ∗ − (cid:19) s Y j =1 (cid:18) c/j a ) p ∗ p ∗ + 1 (cid:19) =: A s . (7)It is easy to see that the relative error between A and its approximation A s is proportionalto 1 /s a p ∗ − with the asymptotic constant c p ∗ / (( p ∗ +1) 2 ap ∗ − ) Q ∞ j =1 (1+( c/j a ) p ∗ / ( p ∗ +1)).A general idea of our construction is to select sets u with large γ u and subtract γ u from A s . This is repeated until A s is reduced to or below ε p ∗ .More specifically, consider a partition of R + into intervals I i such that the numbers in I j are greater than those in I j +1 . For simplicity, we used I = [10 − , ∞ ) , and I j = [10 − j , − j +1 ) for j = 2 , , . . . . in our numerical experiments when constructing quasi-optimal sets. However, we thinkthat a better partition is possible, especially when constructing optimal active sets. Wealso associate with every interval a list L j that contains those u for which γ u has beensubtracted from A s in j th step.In the first j = 1 step, add the empty set to L and subtract γ ∅ = 1 from A s . If thenew A s satisfies A s ≤ ε p ∗ , then terminate. Otherwise consider non-empty sets u in theorder of increasing cardinalities. Hence start with singleton sets u = { i } for i = 1 , . . . , k ,where k is the largest integer such that γ { k } is in I . Place { k + 1 } into list L , and startsubtracting from A s the values γ { i } and store { i } in L until either the difference becomesless than or equal to ε p ∗ , in which case we terminate, or i = k . Next repeat the same forsets of cardinality 2, starting with sets { , i } for i ≤ k , where now k is the largest integersuch that γ { ,k } ∈ I . Store { , k + 1 } in L . Next consider sets { , i } , { , i } , etc., untileither all cardinality 2 sets corresponding to the current interval have been visited or thenew value of A s is ≤ ε p ∗ , in which case we terminate. Continue working through the setsin order of increasing cardinality ℓ until γ { , ,...,ℓ } / ∈ I and ℓ ≥ c (of course, ℓ is always atleast c for c ≤ j = 2. The procedure in this (and later) steps isvery similar except that for a fixed cardinality of u , we check if any such set has alreadybeen placed in L in the 1st step. If it has, we start working with such sets first. Forinstance, for cardinality 1, if { k + 1 } ∈ L , then we begin with sets { i } for i ≥ k + 1.For cardinality 2, if { i i , i } has been placed in L , then we inspect sets { i , i } for i ≥ i ,before any other sets of cardinality two are considered. Once we find γ { i,i +1 } / ∈ I , westore { i, i + 1 } in L and proceed to sets of cardinality 3, etc.At the very end, U q − opt ( ε ) consists of all subsets u whose values γ u were subtractedfrom A s . The main procedure is outlined in Algorithm 2. In all of the algorithms j max and ℓ max are computational thresholds denoting, respectively, the maximum number ofintervals to be searched through and the maximum allowed cardinality of sets.To search through the sets in a systematic way, we must keep track of the currentset, u , and the index, i , that we are incrementing from. The subroutine increment- u outlined below in Algorithm 1 details how to increment u from index i .11 lgorithm 1 (Subroutine: increment- u ) inputs : u , i output : u u i ← u i + 1 ⊲ updating u i first for r = i + 1 , i + 2 , . . . , | u | do ⊲ incrementing u from index i + 1 u r ← u i + r − i end for return u Algorithm 2 (Constructing the quasi-optimal active set) inputs : ε , p ∗ , s , (¯ γ u ) | u | < ∞ , ( I j ) j max j =1 output : U q − opt ( ε ) L j ← ∅ for all j = 1 , , . . . , j max ⊲ initialising U q − opt ( ε ) ← {∅} T ← A s − ε p ∗ − ¯ γ ∅ ⊲ T tracks difference between A s − ε p ∗ and weights if T ≤ then return U q − opt ( ε ) ⊲ quasi-optimal active set is complete for j = 1 , , . . . , j max do ⊲ looping over intervals ⊲ first handle sets found at previous step ( U q − opt ( ε ) , T, ℓ next , L j +1 ) ← q-opt-search ( U q − opt ( ε ) , T, (¯ γ u ) | u | < ∞ , L j , L j +1 , I j ) if T ≤ then return U q − opt ( ε ) ⊲ quasi-optimal active set is complete for ℓ = ℓ next , ℓ next + 1 , . . . , ℓ max do ⊲ search through unvisited sets u = { , , . . . , ℓ } i ← ℓ ⊲ i keeps track of index to increment u from if ¯ γ u / ∈ I j and ℓ ≥ c then break ⊲ no more u with ¯ γ u ∈ I j while i > do ⊲ when i = 0 there are no more u of cardinality ℓ if ¯ γ u ∈ I j then add u to U q − opt ( ε ) T ← T − ¯ γ u if T ≤ then return U q − opt ( ε ) ⊲ quasi-optimal set is complete i ← ℓ ⊲ continue incrementing from last index else add u to L j +1 i ← i − ⊲ start incrementing from lower index if i = 0 then break ⊲ go to next cardinality end if u ← increment- u ( u , i ) end while end for end for To make the presentation clearer Algorithm 2 is broken into two parts: First, we searchstarting from the sets found in the previous interval, which is handled by the subroutine q-opt-search in Algorithm 3. Then we continue searching through sets in order ofincreasing cardinality (line 9) starting where q-opt-search finished, at cardinality ℓ next .12he basic search structure is the same, however in q-opt-search each set we visit ischecked to reduce multiple visits to a single set and ensure that the same set is not addedto U q − opt ( ε ) more than once.The notation( U q − opt ( ε ) , T, ℓ next , L j +1 ) ← q-opt-search ( U q − opt ( ε ) , T, (¯ γ u ) | u | < ∞ , L j , L j +1 , I j ) , denotes that we call q-opt-search with inputs U q − opt ( ε ), T , (¯ γ u ) | u | < ∞ , L j , L j +1 , I j andthen use the output to update U q − opt ( ε ), T , ℓ next and L j +1 . Algorithm 3 (Subroutine: q-opt-search ) inputs : U q − opt ( ε ), T , (¯ γ u ) | u | < ∞ , L j , L j +1 , I j outputs : U q − opt ( ε ), ℓ next , T , L j +1 ℓ next = 1 for u ∈ L j do i ← | u | while i > do L j ← L j \ u ⊲ reducing the double-handling of sets if ¯ γ u ∈ I j then if u ∈ U q − opt ( ε ) then break ⊲ already visited u and any future increments add u to U q − opt ( ε ) T ← T − ¯ γ u if T ≤ then return ( U q − opt ( ε ) , ℓ next , T, L j +1 ) i = | u | else add u to L j +1 ⊲ u to be checked first in next interval i ← i − if i = 0 then break ⊲ go to next u ∈ L j end if u ← increment- u ( u , i ) end while ℓ next = | u | + 1 ⊲ main search will start at cardinality | u | + 1 end for return ( U q − opt ( ε ) , ℓ next , T, L j +1 )The construction of optimal active sets is very similar. The main difference is that inthe j th step, we first create the list L unsorted j , order its elements u ∈ L unsorted j according todecreasing values of γ u , and next start subtracting the values γ u from A s . Again the lists L j will hold the sets visited in the previous interval.In fact, if we do not care whether or not all of the sets are ordered but only that U opt ( ε ) consists of the sets with the largest weights, then we only need to sort the setswhich come from the final interval. This is because at the previous intervals all of thesets will need to be added to U opt ( ε ), regardless of sorting. To do this in practice, for eachinterval we store the sum of all the weights corresponding to that interval. In Algorithm 4we denote this by T j . At the end of the j th step, we check whether I j is the final interval,13.e., if A s − P ji =1 T i ≤ ε p ⋆ , if so we sort the sets and add them one-by-one until the activeset is complete. Otherwise we add all of the sets in L unsorted j to U opt ( ε ) and go to the nextinterval. For completeness, the construction of optimal active sets is detailed separatelybelow in Algorithm 4 and the subroutine opt-search in Algorithm 5. Algorithm 4 (Constructing the optimal active set) inputs : ε , p ∗ , s , (¯ γ u ) | u | < ∞ , ( I j ) j max j =1 output : U opt ( ε ) T j ← L unsorted j ← ∅ and L j ← ∅ for all j = 1 , , . . . , j max ⊲ initialising U opt ( ε ) ← {∅} T ← A s − ε p ∗ − ¯ γ ∅ ⊲ T tracks difference between A s − ε p ∗ and weights if T ≤ then return U opt ( ε ) ⊲ optimal active set is complete for j = 1 , , . . . , j max do ⊲ looping over intervals ⊲ first handle sets found at previous step ( ℓ next , T j , L unsorted j , L j +1 ) ← opt-search ( U opt ( ε ) , T j , (¯ γ u ) | u | < ∞ , L unsorted j , L j +1 , I j ) for ℓ = ℓ next , ℓ next + 1 , . . . , ℓ max do ⊲ search through unvisited sets u = { , , . . . , ℓ } i ← ℓ ⊲ i keeps track of index to increment u from if ¯ γ u / ∈ I j and ℓ ≥ c then break ⊲ no more u with ¯ γ u ∈ I j while i > do ⊲ when i = 0 there are no more u of cardinality ℓ if ¯ γ u ∈ I j then add u to L unsorted j T j ← T j + ¯ γ u i ← ℓ ⊲ continue incrementing from last index else add u to L j +1 i ← i − ⊲ start incrementing from lower index if i = 0 then break ⊲ go to next cardinality end if u ← increment- u ( u , i ) end while end for if T j ≥ T then ⊲ sorting step, first check if I j is the last interval sort L unsorted j for u ∈ L sorted j do ⊲ add sorted sets until active set is complete add u to U opt ( ε ) T ← T − ¯ γ u if T ≤ then return U opt ( ε ) ⊲ optimal active set is complete end for else ⊲ add all sets for the current interval and continue search add all u to U opt ( ε ) T ← T − T j end if end for lgorithm 5 (Subroutine: opt-search ) inputs : U opt ( ε ), T j , (¯ γ u ) | u | < ∞ , L unsorted j , L j +1 , I j outputs : ℓ next , T j , L unsorted j , L j +1 ℓ next = 1 for u ∈ L j do i ← | u | while i > do L j ← L j \ u ⊲ reducing the double-handling of sets if ¯ γ u ∈ I j then if u ∈ L unsorted j then break ⊲ already visited u and any future increments add u to L unsorted j T j ← T j + ¯ γ u i = | u | else add u to L j +1 ⊲ u to be checked first in next interval i ← i − if i = 0 then break ⊲ go to next u ∈ L j end if u ← increment- u ( u , i ) end while ℓ next = | u | + 1 ⊲ main search will start at cardinality | u | + 1 end for return ( ℓ next , T j , L unsorted j , L j +1 ) In this paper we have introduced the notion of superposition dimension and optimal activesets to be used in the MDM for multivariate integration and presented an algorithmdetailing their construction. We also introduced a second simplified, computationallyless intensive version of the algorithm, which constructs quasi-optimal active sets. Ournumerical results show that the quasi-optimal active sets are of a similar size to the optimalactive sets. Often the two sets are exactly the same. In all of our numerical results theoptimal and quasi-optimal active sets are smaller than, and have superposition dimensionless than or equal to, the active sets using the construction in [7].To observe how different choices of parameters a and c affect our construction, statis-tics on the resulting optimal active sets are given in Tables 1-4. Tables 1 and 2 give,respectively, the size and the superposition dimension of the optimal active set for p = 2and an error request of 10 − . For p = ∞ , ε = 10 − the size and superposition dimensionof the optimal active sets are given in Tables 3 and 4. The results for the quasi-optimalactive set are again very similar and so have not been included here. As expected theseresults demonstrate that as the decay of the weights is slower or the weights become larger( a smaller and c larger) the problem becomes more difficult and the active sets are bynecessity larger. However the superposition dimension remains relatively small, at most 6.15 c | U opt (10 − ) | for p = 2 anddifferent a , c . ac d ( U opt (10 − )) for p = 2 anddifferent a , c . ac | U opt (10 − ) | for p = ∞ anddifferent a , c . ac d ( U opt (10 − )) for p = ∞ and different a , c .Finally, we have constructed the sets U q − opt ( ε ), U opt ( ε ), and U PW ( ε ) for the productweights with c = 1. They are listed explicitly in the Appendix. Acknowledgements
The authors would like to thank the two anonymous referees for their constructive com-ments which helped improve the paper.
Appendix
We list here a selection of the constructed active sets from the previous sections (the verylargest sets have been omitted). To save the space sometimes we write [ ... { x , . . . , x k , x k +1 } ]to denote the sequence of sets,[ ... { x , . . . , x k , x k +1 } ] = { x , . . . , x k , x k + 1 } , . . . , { x , . . . , x k , x k +1 } . For instance [ ... { } ] denotes { } , { } , { } and [ ... { , } ] denotes { , } , { , } , { , } , { , } . Case p = 1 and a = 4 U PW (10 − ) = {∅ , { }} , U PW (10 − ) = {∅ , [ ... { } ] , { , } , { , }} , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] } . Case of p = 1 and a = 3 U PW (10 − ) = {∅ , { } , { } , { , }} , U PW (10 − ) = {∅ , [ ... { } ] , [ .. { , } ] } , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , } , { , } , { , , } , { , , }} . Case p = 1 and a = 2 U PW (10 − ) = {∅ , [ ... { } ] , { , } , { , }} , PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , } , { , } , { , , } , { , , } , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , { , } , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , { , , } , { , , } , { , , } , { , , , } , { , , , }} . Case p = 2 and a = 4 U q − opt (10 − ) = U opt (10 − ) = {∅ , { }} , and U PW (10 − ) = {∅ , { } , { }} , U q − opt (10 − ) = U opt (10 − ) = {∅ , { } , { } , { , }} and U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] } , U q − opt (10 − ) = U opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] } U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , } , { , } , { , , }} . Case of p = 2 and a = 3 U q − opt (10 − ) = U opt (10 − ) = {∅ , { }} and U PW (10 − ) = {∅ , [ ... { } ] , { , }} , U q − opt (10 − ) = U opt (10 − ) = {∅ , [ ... { } ] , { , } , { , }} , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , } , { , , }} , U q − opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , { , , }} , U opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , } , { , } , { , , } , { , , }} , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , { , } , [ ... { , , } ] , [ ... { , , } ] } . Case of p = 2 and a = 2 U opt (10 − ) = {∅ , { } , { } , { , }} and U q − opt (10 − ) = {∅ , [ ... { } ] , { , }} , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , }} U q − opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , { , , }} , U opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , , } ] } , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , { , , } , { , , } , { , , }} . Case p = ∞ and a = 4 U q − opt (10 − ) = U opt (10 − ) = {∅ , { }} and U PW (10 − ) = {∅ , [ ... { } ] , { , } , { , }} . U q − opt (10 − ) = U opt (10 − ) = {∅ , [ ... { } ] , { , }} , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , } , { , } , { , , }} , U q − opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] } , U opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , }} U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , { , } , [ ... { , , } ] , [ ... { , , } ] } Case p = ∞ and a = 3 U q − opt (10 − ) = U opt (10 − ) = {∅ , { } , { }} , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , } , { , } , { , , }} , U q − opt (10 − ) = U opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , }} , PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , { , } , { , } , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , { , , } , { , , } , { , , } , { , , , }} , U q − opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , , } ] } , U opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , , } ] , [ ... { , , } ] } , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , { , } , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , { , , } , { , , } , [ ... { , , } ] , [ ... { , , } ] , { , , } , { , , } , [ ... { , , , } ] , [ ... { , , , } ] , [ ... { , , , } ] , { , , , } , { , , , } , [ ... { , , , } ] , { , , , }} Case p = ∞ and a = 2 U q − opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , { , }} U opt (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , , } ] } , U PW (10 − ) = {∅ , [ ... { } ] , [ ... { , } ] , [ ..., { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , } ] , [ ... { , , , } ] , [ ... { , , , } ] , [ ... { , , , } ] , [ ... { , , , } ] , [ ... { , , , } ] , [ ... { , , , } ] , { , , , } , [ ... { , , , } ] , [ ... { , , , } ] , [ ... { , , , } ] , { , , , } , [ ... { , , , } ] , { , , , } , [ ... { , , , } ] , { , , , } , { , , , , }} . References [1] R. E. Caflisch, W. Morokoff, and A. B. Owen: Valuation of mortgage backed securitiesusing Brownian bridges to reduce the effective dimension. Journal of ComputationalFinance 1: 27–46, 1997.[2] M. Hefter, K. Ritter, and G. W. Wasilkowski, On equivalence of weighted anchoredand ANOVA spaces of functions with mixed smoothness of order one in L or L ∞ norm, J. Complexity (2016), 1-19.[3] P. Kritzer, F. Pillichshammer, and G. W. Wasilkowski, Very low truncation dimen-sion for high dimensional integration under modest error demand, J. Complexity (2016), 63-85. 184] F. Y. Kuo, I. H. Sloan, G. W. Wasilkowski, and H. Wo´zniakowski, Liberating thedimension, J. Complexity (2010), 422-454.[5] R. Liu and A. Owen: Estimating mean dimensionality of analysis of variance decom-positions. J. Amer. Statist. Assoc. 101: 712–721, 2006.[6] A. B. Owen: Effective dimension for weighted function spaces. Technical Re-port, Stanford University, 2012 and in revised version 2014. available under http://statweb.stanford.edu/ ∼ owen/reports/effdim-periodic.pdf [7] L. Plaskota and G. W. Wasilkowski, Tractability of infinite-dimensional integrationin the worst case and randomized settings, J. Complexity (2011), 505–518.[8] I. H. Sloan and H. Wo´zniakowski, When are Quasi-Monte Carlo algorithms efficientfor high dimensional integrals?, J. Complexity (1998), 1-33.[9] X. Wang and K.-T. Fang: Effective dimension and quasi-Monte Carlo integration. J.Complexity 19: 101–124, 2003. Authors’ Addresses:
A. D. GilbertSchool of Mathematics and StatisticsThe University of New South WalesSydney, NSW 2052, AustraliaE-mail: [email protected]
G. W. WasilkowskiDepartment of Computer ScienceUniversity of Kentucky301 David Marksbury BuildingLexington, KY 40506, USAE-mail: [email protected]
Funding: