ConceptioArchivearXiv CS
arXiv CSopen access

Graph Structured Combinatorial Semi-Bandit with Nonlinear Reward Associations through Separable Signals

Unknown · 2026 · arxiv_cs
arXiv CS · Papers · License: Open Access · 2026
Open Source ↗Direct PDF ↓
machine learning, deep learning, neural networks

1

Graph Structured Combinatorial Semi-Bandit with Nonlinear Reward Associations through Separable Signals

arXiv:2606.14650v1 [cs.LG] 12 Jun 2026

Christoph Bauschmann and Setareh Maghsudi Department of Electrical Engineering and Information Technology Ruhr University Bochum Bochum, Germany

Abstract—The identification of optimal structures within vast arrays of interconnected data necessitates significant samplingand computational effort. Learning and leveraging underlying signal dependencies can improve efficiency and predictive capabilities considerably, but the ubiquity of nonlinear statistical relations amplifies the complexity of such undertakings. In this paper, we develop novel generic and adaptive strategies equipped with routines for graph-based causal reward modeling, analytic reproducing kernel methods, and Taylor approximation of functional processes. We establish theoretical performance guarantees sublinear in time and linear in data volume over time. Our analyses cover robustness to a multitude of uncertainties arising from noise interference, gradual model convergence, and solution space mismatch. The framework’s general appeal is substantiated by a minimalistic set of conditions or reliance on prior estimates, while various outlined modifications address specific or extended settings. To demonstrate practical effectiveness, we conduct numerical experiments using both benchmarked synthetic and real-world transportation datasets.

I. I NTRODUCTION EQUENTIAL decision-making processes without influence on the state transition constitute the class of MultiArmed Bandit (MAB) problems. Upon selection, each arm provides a reward, which may follow a stochastic distribution. The agent’s goal is to maximize the reward accumulated over the game horizon by adapting her actions based on the interaction history. The performance of the controlling policy hinges on its ability to balance queries to underexplored distributions in hopes of future payoff against exploitation of sampled answers for immediate gain. Formulated as a minimization objective, the cumulative difference in scoring to an optimal policy is understood to be the agent’s regret. When the agent is allowed to pull multiple (base) arms simultaneously, the viable subsets become the super arms of a combinatorial MAB. The posterior reward distributions of intersecting super arms are associated through their joint dependence on the underlying base arm distributions. If the correlation is linear, the reward is simply a weighted sum of the base signals. In general, causal relationships can take any nonlinear form without ad hoc decomposition. However, an achieving agent will endeavor to make predictive links between similar sets of arms. The potential advantages over conventional MAB techniques that keep individual records are

S

This work has been submitted to the IEEE for possible publication. Copyright may be transferred without notice, after which this version may no longer be accessible.

especially pronounced for large numbers of base arms as their power set grows exponentially. In a semi-bandit setting, the agent’s feedback can entail auxiliary data besides the actual reward. In the case of a combinatorial MAB, natural candidates are the samples from the chosen base arms’ distributions that determine the ultimate reward. If the reward calculation has a certain structure, intermediate values of that computation might become available. In fact, without restrictions to the set of possible functions, a frequentist sub-binomial regret expectation cannot be guaranteed, as each valid base arm combination brings unique coefficients into effect. In a sense, off-policy evaluation shall be supported in the solution space. One such framework to formalize causal interdependencies are Structural Equation Models (SEMs). They work based on a graph representation that aggregates initial node values over weighted edges. In the bandit context, the agent filters the base signal values as initial impulses to the SEM, and the resulting output is accumulated in the overall reward. First introduced in [1] exclusively for linear edge weights, we expand this notion to general activation functions between vertices and by a kernelized optimization that does not rely on perfect feedback for exact weight identification. Literature Overview. Algorithms in the relevant literature need to be distinguished by the assumptions they operate under. There are many ways in which the uncertainty of the reward generation process may be limited. Often, they take the form of prior expertise or access to an oracle to approximate the results, as done for the general analysis carried out in [2]. Moreover, the randomness may be abstracted by concentrating on special signal distributions or feedback that deviates little from its expected value. For instance, [3] have the differences between mean and observation follow a subgaussian distribution homogeneous across super arms and deploy Gaussian processes and likelihood models. The neural-network-based bandits in [4] rely on contextual information that captures all correlations between arms. Additionally, many works do not explicitly account for noise interference. While unstructured approaches may subsume noise independent between arms in the stochastic reward function, structure identification could exploit unrealistic perfect observations. The conditions set for the causality of the process also warrant attention. The control the agent can exert on the system can range from masking signals to soft/hard intervention on their values [5], [6]. Causal graphs as in [7],

2

[8] additionally introduce independence prerequisites such as faithfulness or effect identifiability. Their underlying logic is commonly framed with linear SEMs [9], [10]. Besides linearity, special properties like monotonicity [2] or submodularity [11] may be imposed upon the functions describing the reward associations. In [12], intermediate rewards are sampled along a hierarchy of parametrized distributions rooted in latent priors. The authors develop a Thompson Sampling (TS) approach and investigate specific examples of nonlinearities only for the final accumulating step. The covariance estimations in [13] capture arbitrary associations that precede the selection, but they do not encompass nonlinearities after the action filter. A different graphical framework is that of probabilistically triggered neighboring arms [14], where high connectivity is beneficial rather than a hindrance to the complexity of the learning process. Extensive foundations and an overview of problem classes are supplied by reference books such as [15]. Contribution. For many of these individual scopes, regret bounds that are sublinear in time steps and at most linear in the number of available base arms (for time-independent terms) are established in the respective references (for a tabular comparison, see Appendix A). However, we contribute the first policies for the stochastic structured combinatorial semi-bandit setting to achieve such bounds when the pairwise associations between signals follow arbitrary analytic functions. To simultaneously discover and infer from this flexible description, the routines we develop negotiate the duplicate exploration of inherent distributions and their causal dynamics with the objective. We refrain from prohibitive restrictions to uncertainty and investigate the effect noise, structural misalignments and functional complexity have on convergence. For greater applicability, we delineate alterations that address related or special scenarios. Our theoretical findings and algorithmic analysis are corroborated by numerical tests on benchmarked synthetic and real-world data. Combinatorial MABs accommodate a wide variety of tasks, such as drug discovery [16], online recommender systems [12], neural network design [17] and analysis (e.g., identifying dominant subnetworks or sets of biases), social influence maximization [11] and live routing (e.g., preserving resources for vehicles or wireless communication). The SEM can be most readily applied when inert and aggregated forms of data are apparent or extractable (e.g., from phases of separation or abnormal behavior within the network). It’s capable of describing spatial graphs and signal processing thereon (e.g., of logistical/biological nature like transportation delays respective infection spreads). A motivating example lies in the distribution of variable grid load. Depending on the current local consumption and their proximity, producers can complement or interfere with each other’s efficiency in intricate ways. An agent would have to decide which contractors to enroll so that a sustainable coverage is ensured at the lowest possible cost. Outline. Over the course of this paper, we first provide the formal specifications in Section II, before developing the suggested decision-making policies in Section III. In Section IV, we then prove regret bounds and account for sources of inaccuracy. Subsequently, we analyze the experimental performance in Section V. Finally, we draw conclusions and

discuss future research directions in Section VI. Notation. We use bold lower/upper case letters for vectors and matrices, respectively. the common multiQ We apply Q index notation ab = a[i]b[i] , a! = i∈[N ] i∈[N ] a[i]!, P ⊤ |a| = |a[i]|, and a mod k = (a[i] mod k) i∈[N ] i∈[N ] for a, b ∈ NN , k, N ∈ N. We write a ⊙ b = diag(a)b for the Hadamard product of vectors. II. P ROBLEM S TATEMENT We represent the set of base arms with [N ] = {1, .., N }. Let B1 , . . . , BN be (a priori unknown) probability distributions on closed intervals of R. At each time step t, the vector of independent instantaneous rewards bt is drawn from the composite distribution B on [0, 1]N . To stay concise, we restrict our description to absolutely continuous distributions on [0, 1] that have nonzero lower bounds to their probability density function (p.d.f.). Shifts to other bounded ranges can w.l.o.g. be incorporated into the objective function. Unrestricted distributions and p.d.f.s are discussed at the end of Section IV respective in Appendix B1 (Remark A.2). For the combinatorial bandit problem, an agent is tasked with sequentially selecting subsets of base arms over time. Each such super arm is characterized by the agent’s decision vector xt ∈ {0, 1}N over the inclusion of each base arm, with the base arm i being chosen iff xt [i] = 1. The total cardinality of valid super arms is limited by s ∈ N, thus restricting the set of allowed super arms to X = {x ∈ {0, 1}N | ∥x∥0 ≤ s}, with ∥ · ∥0 counting the nonzero elements of a vector. We use a (sparse) directed non-self-cyclic graph G = ([N ], E, F) on the base arms with edge set E to model the process by which the reward associated with each super arm is calculated. We collectively represent the functional weights for each edge in a matrix-styled operatorPF of dimension (N, N ) N that maps b ∈ RN to F(b) = ( j=1 F[i, j](b[j]))⊤ i∈[N ] . All component functions F[i, j] : R → R are analytic. In this sense, we also understand I = diag(id) as the identity mapping, and conventionalize concatenation as the default connection between operators. More specifically, we assume the causal relationships follow an additively separable SEM [18]. The exogenous input vector to this model is zt = diag(bt )xt = xt ⊙ bt .

(1)

Its endogenous output vector is subsequently acquired as X yt [i] = F[i, j](yt [j]) + G[i, i](zt [i]), ∀i ∈ [N ]. (2) j̸=i

Here, G is a diagonal matrix-styled operator of exogenous activation functions, while F acts endogenously. The solution yt to the above set of equations constitutes the vector of overall rewards in our bandit setup. Hence, each of its entries is conditioned on the associated instantaneous rewards and the overall reward of its immediate neighbors in G. Since diag(F) = 0, the vectorized form of (2) reads as yt = F(yt ) + G(zt ) = F(yt ) + G(xt ⊙ bt ).

3

A. Online Graph Learning

Fig. 1. Causal Reward Relations Following a Separable SEM.

Our algorithms maintain an explicit world model through the approximation of F in a pre-selected space of solutions. For i ∈ [N ], let H(i) be a reproducing kernel Hilbert space (RKHS) with continuous positive-definite kernel, consisting of analytic, zero-true functions from a closed interval in RN to R that are additively separable over their N arguments and zero-imaged on RI[i]⊤ . We use H to denote their composition to a space of operators on RN . For θ > 0, we set (i) Hθ = {u ∈ H(i) | ∥u∥2H(i) ≤ θ}. By Hθ , we summarize the space of operators restricted in each component of H. We collect the feedback up to time t using the matrices Zt = [z1 , . . . , zt ] and Yt = [y1 , . . . , yt ] = F(Yt ) + Zt . Through kernel optimization, we obtain 1 ∥H(Yt ) − Yt + Zt ∥2F H[i]∈H(i) t X ∥H[i]∥2H(i) , +λ

F̂t = arg min When the operator I − F is invertible, this equation can be formulated explicitly as

i∈[N ]

yt = (I − F)−1 (G(xt ⊙ bt )). Unlike the instantaneous rewards, the overall rewards may take on values in vastly different ranges per dimension. After selecting a super arm, the agent receives the strong semi-bandit feedback of the SEM’s exogenous and endogenous vectors for that time step. Following the bandit framework, the extracted payoff ⊤

r(xt , bt ) = rxt (bt ) = h yt = h (I − F)

−1

(G(zt )),

for an operator h of N functions, serves as the agent’s objective function to be maximized. Equivalently, the agent aims to minimize the expected regret R(T ) = T µ(x∗ ) −

T X

(3)

µ(xt )

t=1

when comparing the expected payoff

where ∥ · ∥F and ∥ · ∥H(i) , i ∈ [N ], are the Frobenius and Hilbert space norms, respectively. The factor λ > 0 regularizes the sparsity of the resulting solution. By Lagrangian duality, it can be represented through a convex restriction ∥H[i]∥2H(i) ≤ θ ∀i ∈ [N ].

(4)

The optimization is independent between dimensions i ∈ [N ] and can also be performed in the respective subspaces. In practice, sparsity or bijectivity of I − F̂t can be maintained through edge recognition thresholds. Rather than performing a joint optimization over the entire (o) sample set, we also consider computing a separate Ft for each o ∈ [N ] based solely on the time steps τ in which xτ [o] = 1. Algorithms with this variant then assemble the entries of their overall model recursively for k ∈ [o − 1] as (o)

F̂t [o − k, o] : R → R, b 7→ (I − F̂t )−1 (bI[o]⊤ )[o − k]

µ(x) = Eb∼B [rx (b)] at each time step up to unknown horizon T ∈ N with that of a sought after optimal decision vector x∗ ∈ arg max µ(x). x∈X

The reward structure depicted in Figure 1 displays all potential functional components. For ease of exposition, the main paper’s focus lies on G = I for the preliminary onnode activations, a linear final accumulation h⊤ = (id, . . . , id), which we shall denote simply by rxt (bt ) = 1⊤ yt , and zerotrue F, i.e. F(0) = 0. We generalize these settings and describe further extensions to the framework, including partial endogenous feedback, in Appendix B1.

o−1 X

(5)

j=o−k+1

Should the indexing in (5) not be inherent to the application, edge recognition thresholds should be increased until a DAG is formed. We will elaborate on the chosen method and rationale behind it when presenting the individual algorithms in III-D. B. Internal UCB To handle the exploration-exploitation trade-off for signal passes through the learned model, the algorithms keep track of empirical moments

III. D ECISION -M AKING S TRATEGY In this section, we introduce our policies for the problem described before. Since we do not assume knowledge of the causal dynamics, the agent has to trade off two sources of uncertainty against exploiting its insights. The first relates to the graph functions in F, while the second regards the signal distribution B. Our proposed strategy consists of three distinct parts, which we describe below.

(o)

F̂t [o − k, j]((I − F̂t )−1 (bI[o]⊤ )[j]).

(p)

ϕ̂t [i] =

t X

(zτ [i] − at+1 [i])p 1 {xτ [i] = 1} /mt [i]

(6)

τ =1

about some at+1 ∈ [0, 1]N up to a preselected degree [q] ∋ p, of the observed instantaneous rewards. t P They keep counters mt = xτ for the number of rounds τ =1

each base arm was part of the selected super arms up to time

4

t. After including each base arm at least once, the internal confidence interval has a radius of p Ct [i] = (s + 1) ln(t)/mt [i], ∀i ∈ [N ]. (7) To leverage this confidence, the expected payoff needs to be expressed in terms of the distribution invariants (6). Our approach is based on a representation of the reward function with its Taylor series expansion around any center point a ∈ [0, 1]N for b ∈ [0, 1]N , x ∈ X : ∞ X

rx (b) =

|α|=0 x⊙α=α

Dα r1 (a ⊙ x) (b − a)α . α!

initial derivatives Dα rx = vanish whenever x ⊙ α ̸= α for α ∈ N . The independence of the underlying distributions also allows us to establish q X Dα rx (a) Y (α[i]) ϕ [i] + EB [Rq+1 (b)], µ(x) = α! Several

of

the

∂ |α| r ∂b[1]α[1] ...∂b[N ]α[N ] x N

|α|=0 x⊙α=α

i∈[N ]

for the series approximation of degree q and moments ϕ(p) [i] = EBi [(b − a[i])p ], i ∈ [N ], p ∈ [q]. The Rq+1 is bounded in absolute value by P remainder 1 γ α! max max |D r1 (c)| ( [19]).

|α|=q+1

|γ|=|α|c∈[0,1]N

To boost internal exploration of the current reward model r̂(t−1) = (I − F̂t−1 )−1 , the algorithms optimistically project values within the moments’ confidence interval that maximize the estimated payoff for each arm x ∈ X : Itmax (x) =

max (p)

(t−1)

|α|≤q x⊙α=α

(at ) Dα r̂x α!

The graph learning process must be supplied with sufficient explorative data to gain certainty over the structural dependencies. Outside the model’s reward simulation, our policies’ selection criteria have a confidence addend of p Et (x) = max x[i]4Ct−1 [i] ln(t − 1) i∈[N ]

for x ∈ X , which favors base arms with little observation data. Functional uncertainty is thus captured by a single confidence interval around rx (a), and required, as derivatives can be arbitrarily close to 0. Introducing intervals for every derivative value Dα rx (a) would risk instability or high computation load for tight estimates. Under the assumption of an underlying SEM, the presence of a base arm in the played super arm will already suffice to build up estimation certainty over the associated outward edge weights. The impact of this rewardindependent bound can be fine-tuned by multiplication with a hyperparameter. The logarithmic growth ensures the influence of the external exploration on the decision making process over time, independent of scaling factors that would otherwise be critical to surpass the sample order dictated by kernel optimization accuracy guarantees. P For short timephorizons, an alternative external UCB of x[i]4Ct−1 [i] ln(t − 1)/s i∈[N ]

Y 

 (α[j]) (α[j]) ϕ̂t−1 [j] + Wt−1 [j] .

j:x[j]=1

t X

zτ [i]/mt [i], i ∈ [N ].

facilitates broader overall information gain, without displacing any T -dependent terms in the theoretical bounds of Section IV. D. SSEM Algorithms

The specific growth behavior of the functions from the RKHS can be used to narrow the range of potentially maximizing values within the confidence interval. For example, if the functions have only non-negative derivatives and all moments are (p) non-negative, one can substitute Wt = 2p mod 2 Ct , ∀p ∈ [q]. Under these circumstances, only super arms of cardinality s would be candidates for maximizing the expected payoff. To center the Taylor approximation around relevant values, at+1 can be chosen as the empirical means β̂ t [i] =

C. External UCB

(p)

Wt−1 ∈RN : Wt−1 [i] =2p mod 2 Ct−1 [i] ∀p∈[q],i∈[N ]

X

lies in extrapolating the moments for replacement with their approximations, while the partial derivatives allow for pointwise evaluation. Otherwise, the entire function chain would have to be resolved symbolically to identify the occurrences of each power in the payoff, or be subject to another kernel estimation.

(8)

τ =1

In combination, one can also track unbiased estimates of the central moments (as is the case in Algorithm 3 in Appendix B2). A fixed at on the other hand eases the empirical moment update and potentially stabilizes the learning progress. Especially at = 0 simplifies derivative evaluation to a single point across all super arms. Many of our results would also hold for other multivariate polynomial approximation techniques with bounded coefficients. A practical advantage of working with the Taylor series

All variants of our algorithm are initialized with a series of ⌈N/s⌉ super arms covering the entire base arm set to obtain first estimates of the moments of B[i], ∀i ∈ [N ]. For remaining time steps t up to hidden time horizon T ∈ N, they then iterate over the mechanics described in the previous subsections. First, F̂t−1 is learned and used to conjecture a reward function approximation r̂(t−1) . Subsequently, the unbiased moment estimates and confidence bounds are updated with the most recent observations. A super arm x ∈ X maximizing a weighted sum of the internal and external UCBs, Itmax (x) + ∆It Et (x),

(9)

is selected for play, with ties broken at random, and new observations of resulting endo- and exogenous vectors are stored. To strike a balance between the UCBs, we choose the internal span ∆It = ′max Itmax (x′ ) − Itmax (x′′ ) as the ′′ x ,x ∈X

dynamic weighting factor. Alternatively, the internal UCB can be limited to an estimate of the maximal true payoff. The variants differ as follows: SSEM-UCB (Alg. 1) learns graphs from the samples associated with each base arm node individually, and combines them suitable for a Directed

5

Algorithm 1 SSEM-UCB: Separable Structural Equation Model - Upper Confidence Bound Input: Taylor degree q, centers (at )t∈[T ] ; separable RKHSs (H(i) )i∈[N ] , norm bound θ. 1: for t = 1, . . . , ⌈N/s⌉ do 2: Select xt with xt [i] = 1 ⇐⇒ (t − 1)s < i ≤ ts. 3: Observe zt and yt . 4: end for 5: for t = ⌈N/s⌉ + 1, . . . , T do 6: for o = 1, . . . , N do (o) 7: Obtain F̂t−1 as stated in (3) and (4) but exclusively from feedback during {τ ∈ [t − 1] | xτ [o] = 1}. 8: Define F̂t−1 [o − k, o] for k ∈ [o − 1] following (5) as entries of upper triangular operator F̂t−1 . 9: end for 10: Set r̂(t−1) = 1⊤ (I − F̂t−1 )−1 . (p) 11: Calculate Ct−1 and (ϕ̂t−1 )p∈[q] according to (7), (6). 12: Select decision vector xt that maximizes (9). 13: Observe zt and yt . 14: end for Algorithm 2 SSEM-UCB-JO: Separable Structural Equation Model - Upper Confidence Bound - Joint Optimization Exactly as SSEM-UCB, only collapsing the loop over o ∈ [N ] to a single optimization, and exchanging a timed factor in the external UCB. 6: Obtain F̂t−1 as stated in (3) and (4). 12: Select decision vector xt that instead of (9) solves q√ max I arg max It (x) + ∆t Et (x) t − 1/ ln(t − 1). x∈X

Acyclic Graph (DAG). This construction will support our theoretical regret analysis for DAGs, which will argue over this very prediction quality. Less computation intensive and unpresumptuous in structure, SSEM-UCB-JO (Alg. 2) optimizes for a single joint graph per time step, at the cost of necessitating an ordinal increase in T to its external UCB. Due to its greater generality and computational efficiency, we will employ this version for empirical deployment. In Appendix B2, we outline an alteration tailored to normal instantaneous reward distributions in particular, alongside a heuristic for choosing the cut-off degree q. Further, we analyze algorithmic complexity and potential failure modes in Appendix C.

by backward substitution, and could more generally be relaxed through pseudo-inversion. Assumption IV.2. To avoid asymptotic approximations, for Alg. 1 we presume the RKHS of H to be chosen in such a way that the possible entries of F̂t−1 as in (5) are from a finite dimensional function space (e.g., polynomials up to a degree). We use the following definitions in our regret analysis. For any decision vector x ∈ X , let ∆(x) = µ(x∗ ) − µ(x) be the suboptimality gap. We define ∆max = max ∆(x), ∆min = x∈X

min

x∈X :µ(x)<µ(x∗ )

∆(x) and ∆I = max ∆Iτ . Moreover, let τ ∈[T ]

−1) Dα r̂(τ (aτ ) , xτ X −1) Γ(q) Dγ r̂(τ (c) /α!, max = max max max max x

wmax = max

max

τ ∈[T ] |α|≤q+1

τ ∈[T ]|γ|=q+1 x∈X c∈[0,1]N



(q,s)

ψmin = min

|α|=q+1

− min{q,s}

min{q, s}s wmax (q + 1)s

,

1 ∆I



(q)

∆min − 4Γmax . 8

(q)

Theorem IV.3. Let G be a DAG and Γmax < ∆min /4 for Taylor approximation degree q. Under Asm. IV.1 and, for Alg. 1, Asm. IV.2 , the expected regret of the algorithms is upper bounded by "

" TB,Hθ ,N +

4(s + 1)KT ln T (q,s) 2

ψmin

# # π2 (2qs + 1) N ∆max , +1+ 3

with KT = ln(T ) for Alg. 1 (SSEM-UCB) and KT = Alg. 2 (SSEM-UCB-JO). [E]

T for

Here, TB,Hθ ,N ∈ N corresponds to the length of respective initialization phases dictated by the deployed kernel optimization. As our method of choice involves various invariants of B and H, its value is left for specification in the proof. Dependence on N may also be reallocated to the maximal graph degree and path length (see Supplementary E4). Section D furthermore offers different intuitions through alternative expressions of involved problem invariants. A version of the results from Thm. IV.3 for a confidence interval length and regret bound in O(ln(T )N ) is possible, but would build on a priori unknown B-dependent factors for the selection criterion, which could in principle be approximated over time. That regret bounds of this ordinality are close to optimal is apparent from the problem-dependent lower bounds Ω(ln(T )N/∆min ) in [20], which apply even to a priori known reward functions and thereby obsolete endogenous feedback. Conversely, Alg. 2 carries less computational expense for optimization and evaluation, and has more general appeal. (q)

IV. T HEORETICAL A NALYSIS In this section, we prove upper bounds on the expected regret of the SSEM algorithms. Each theorem closes by listing the supplementary material section of its proof. Assumption IV.1. The true operator is realizible, i.e., F ∈ Hθ , and unique, i.e., (I − F)−1 exists. Realizibility can subsequently be foregone for certain choices of dense RKHSs, such as linearized polynomials (see Rem. IV.5). Uniqueness is for example guaranteed for DAGs

Proposition IV.4. For any G, provided Γmax < ∆min /4 and Asm. √IV.1, the time-dependent expected regret of Alg. 2 is O(N T ln(T )). [F2] (q)

Instead of preconditioning Γmax < ∆min /4, in Cor. A.3 (Appendix D) we restate the regret bound of Thm. IV.3 (q) for a 4Γmax -approximate optimal solution. There, we also derive specialized bounds by assuming normally distributed instantaneous rewards (Cor. A.4) or relying on the absence of uncertainty (Cor. A.5). The rationale of Thm. IV.3 and SSEMUCB however is robust in the following ways.

6

Remark IV.5. To estimate the performance on reward functions outside the solution space described by H, let q φHθ = max min max E[(F − F̂)(yx (b))[i]2 ] i∈[N ]F̂[i]∈H(i) x∈X θ

with yx (b) = y(x ⊙ b) denote the optimal obtainable objective gap. The proof of Thm. IV.3 actually tracks this error in TB,Hθ ,N and shows that its growth in q is at most exponential. For a multivariate linearized polynomial     (i)   N dim H P P RKHS in RN → R, z 7→ rj,k z[j]k rj,k ∈ R ,   j=1 k=1   j̸=i

the selection error can be expressed based on the qth Taylor remainder of the payoff as well. Since the remainder bound falls in a factorial of q, matching q-exponential bounds on the derivatives of r would suffice to validate Thm. IV.3 even without assuming realizibility. Further, Thm. IV.3 can easily be adapted to cover multiple types of noise interference, as illustrated next. Theorem IV.6. Suppose that the exogenous vector is subject to bounded white observation noise nt ∼ N with B + N supported on [0, 1]N . Additionally, the SEM may feature unobserved white model noise et ∼ M on RN so that yt = F(yt )+zt +et , ∀t ∈ [T ]. Let nt , mt , bt be independent across dimensions and among each other. Then, the regret bound of Thm. IV.3 holds under the same conditions. [F6] The noise directly affects the range of the optimization target and thus both wmax and TB,Hθ ,N . The implicit rescaling of the reward function to keep the noisy exogenous observation in [0, 1]N adds to this effect. Variants of Thms. IV.3 and IV.6 for subgaussian signals respective noise can be derived from martingale forms of Azuma’s inequality [21] in E2. For endogenous observation noise, see Remark A.2. V. N UMERICAL A NALYSIS This section is dedicated to the experimental demonstration and performance elaboration of the proposed policy in suitable scenarios. We discuss implementation details and comparisons to a conceptually diverse assortment of benchmark algorithms. Our analysis extrapolates problem quantities indicative for the success and robustness of our approach on noisy synthetic as well as real-world datasets. Benchmarks. We test SSEM-UCB-JO with multiple choices for the degree of approximation q. For its lowest value 1 and linear function estimates, the method corresponds to a slightly modified version of the SEM-UCB algorithm in [1]. Relying on the similar assumption of a separable action space, SGB [11] accumulates its super arm through gradual inclusion of base arms with the best sample average over a preset number of rounds, before committing and exploiting the final super arm for the remainder of the known time horizon. The structural bandit algorithm METS [12] expects a joint effect prior from which the parameters of action prior distributions and subsequently rewards are sampled. For fairness, we provide it with the distribution type and parameter ranges of the exogenous distributions and have its computation

Fig. 2. Policies’ Regret Confidence Across Runs.

tree form an SEM. For a kernelized bandit comparison, we include IGP-UCB [3] which calculates reward estimates from kernel matrices between previous and suggested super arms following Gaussian likelihood models with hyperparameters as specified in Section G1. A. Synthetic Tests For the synthetic setup we allow super arms of cardinality up to s = 6 and generate DAG adjacency matrices with an edge density of 0.2. We conduct trials for 10 base arms and a mixture of linear and periodic endogenous functions with coefficients in [0.4, 0.7]. The exogenous inputs are sampled from normal distributions with means in range [0.3, 0.7] and standard deviations in [0.1, 0.4], rescaled to be supported almost entirely on [0, 1]. Two further sets of trials with different distributions, component functions, parameter ranges, or more base arms are compared in G1. The observations are obfuscated by white noise from standard normal distributions, exogenously scaled by 0.02, and endogenously by 2 % of the average reward. At each time step, we provide each algorithm with the same exogenous sample masked by its individual super arm choice. The regularization weight of the kernel optimization is tuned adaptively during runtime so that the empirical graphs are appropriately sparse. Without access to density estimates, a computationally more expensive validation split could determine the best value per step. Recommendations towards an efficient implementation are given in Section G1. We determine the optimal reward empirically and measure the time-averaged regret each algorithm incurs in total up to the current step. In this metric, the algorithms eventually converge to the regret of the super arm they predict to be optimal. By incorporating the problem structure into their arm search, the SSEM instances and METS incur much less initial regret than the SGB and IGP-UCB benchmarks and settle for closer to optimal arms sooner. Although the highest derivative approximation of the analytic reward function performs the best among them, the increase is not necessarily steady between degrees, as higher Taylor remainders are not always smaller in absolute expectation.

7

In Figure 2, only SSEM-UCB-JO of degree 4 and METS achieved near-zero regret average in time, the latter of which had additional access to bounds on all synthetic distribution parameters. Lower degrees of SSEM-UCB-JO and the remaining benchmarks do not reliably converge to the true optimal arm due to violations of their assumptions. In the case of SSEM-UCB-JO, the more simplistic Taylor approximations are not sufficiently close to the reward series, whereas IGP-UCB and SGB require the expected reward to be in subgaussian distance respective separable. The unstructured benchmarks also experience only the overall effect of the exogenous noise, which affects each super arm differently and in ways not covered by their theoretical guarantees. In addition, the 95% confidence intervals across 10 runs exhibit a lower diameter and hence greater certainty for SSEM-UCBJO of higher degrees. B. Train Delay Dataset To demonstrate the application of SSEM-UCB-JO with degree q = 3 in real-world networks, we work with a dataset on train delays in Germany [22]. The set contains data points for nearly 2, 000 stations taken at an hourly rate on 16 days between 23.06.2024 and 14.07.2024. We extract a local subset of 20 stops that are at least categorized as urban hubs (category 3 or less) in geographic locations closest to 10.5◦ N, 51.2◦ E, covering a circular area ca. 95 km in radius. Each data entry records the deviation of a train from its current schedule with its complete track record. Judging by the direct connections featured on the routes, the physical railway system in the chosen region exhibits a graph density close to 17 %, with various factors potentially inducing additional associations between stations. We consider the time difference between planned and actual stay, averaged over all trains at a stop per hour, to be the independent exogenous signal inherent to that node in the transportation network. The train-averaged total delay per station accumulated over the course of each hour is understood as the associated endogenous signal under the complete super arm encompassing the entire set of vertices. Endogenous feedback for other actions is generated according to the SEM computation structure after fitting a graph to the aforementioned data. Projected graph adjacencies at the very least have to include direct connections listed regularly in the track histories and can form cycles. The model is allowed to feature component functions of linear, cubic or (scaled and shifted) hyperbolic tangent variety. Figure 3 displays the learned model’s progress when reconstructing the overall delay from the immediate delay based on earlier observations. For full names as well as details to the implementation and further visualization, we refer to Section G2. For the agent, we permit super arms of up to 6 stations with the goal of identifying systemic problems in clusters of the railroad infrastructure. Our simplified test setting does not account for periodic differences between days, nor does it involve the accompanying information on causes of obstruction. After the initialization phase, we find the rate of change in super arm selection to be generally decreasing, with shifts

Fig. 3. Computational Graph Learning Progress.

occurring whenever a better performing combination was explored. The resulting cluster of stations is predicted as chief in overall delay when only selected stations contribute their inherent delay. This subnetwork analysis differs from the apparent set of stations with top overall delay in that it seeks to account for the impact of each station’s delay reverberated throughout the system. The final super arm does contain 3 stations with high inherent delay, and 3 city stops in close proximity of each other. VI. C ONCLUSION By only relying on very elemental properties, the methodologies we introduced are highly customizable through kernel selection and rebalancing of trade-offs befitting specific fields. Expert knowledge can optionally impose additional constraints and guide initial guesses, while the model’s causal coherence renders it resilient to a fair amount of inaccuracies and interference. As evidenced by the regret bounds we obtained, in time the algorithms greatly condense the combinatorial intricacy to a linear growth in base arm set size whilst staying sublinear in step count. We exemplified the transfer to data analysis in realistic settings, where they successfully competed against several benchmarks. Subsequent research could further improve computational scalability through a TS variant that circumvents the determination of optimistic scores by querying instantaneous reward distribution replicas. The parametrization of our model already allows for dimensionality reduction by associating exogenous vectors to multiple nodes and thus reducing the ratio to endogenous outputs. Infinite action spaces can also be covered by contextual or continuously budgeted versions. One way to capture the full reinforcement learning problem beyond bandits would be to extend the fundamental SEM recurrently with time-lagged adjacencies [18].

8

R EFERENCES [1] B. Nourani-Koliji, S. Ghoorchian, and S. Maghsudi, “Linear combinatorial semi-bandit with causally related rewards,” in Proceedings of the Thirty-First International Joint Conference on Artificial Intelligence, IJCAI-22, L. D. Raedt, Ed. International Joint Conferences on Artificial Intelligence Organization, 7 2022, pp. 4878–4884, main Track. [2] W. Chen, W. Hu, F. Li, J. Li, Y. Liu, and P. Lu, “Combinatorial multiarmed bandit with general reward functions,” in Advances in Neural Information Processing Systems, D. Lee, M. Sugiyama, U. Luxburg, I. Guyon, and R. Garnett, Eds., vol. 29. Curran Associates, Inc., 2016. [3] S. R. Chowdhury and A. Gopalan, “On kernelized multi-armed bandits,” in Proceedings of the 34th International Conference on Machine Learning, ser. Proceedings of Machine Learning Research, D. Precup and Y. W. Teh, Eds., vol. 70. PMLR, 06–11 Aug 2017, pp. 844–853. [4] T. Hwang, K. Chai, and M.-H. Oh, “Combinatorial neural bandits,” in Proceedings of the 40th International Conference on Machine Learning, ser. Proceedings of Machine Learning Research, A. Krause, E. Brunskill, K. Cho, B. Engelhardt, S. Sabato, and J. Scarlett, Eds., vol. 202. PMLR, 23–29 Jul 2023, pp. 14 203–14 236. [5] Z. Yan, D. Wei, D. A Katz, P. Sattigeri, and A. Tajer, “Causal bandits with general causal models and interventions,” in Proceedings of The 27th International Conference on Artificial Intelligence and Statistics, ser. Proceedings of Machine Learning Research, S. Dasgupta, S. Mandt, and Y. Li, Eds., vol. 238. PMLR, 02–04 May 2024, pp. 4609–4617. [6] S. Feng and W. Chen, “Combinatorial causal bandits,” Proceedings of the AAAI Conference on Artificial Intelligence, vol. 37, no. 6, p. 7550–7558, Jun. 2023. [7] Y. Lu, A. Meisami, and A. Tewari, “Causal bandits with unknown graph structure,” in Advances in Neural Information Processing Systems, M. Ranzato, A. Beygelzimer, Y. Dauphin, P. Liang, and J. W. Vaughan, Eds., vol. 34. Curran Associates, Inc., 2021, pp. 24 817–24 828. [8] M. Konobeev, J. Etesami, and N. Kiyavash, “Causal bandits without graph learning,” in Proceedings of the Fourth Conference on Causal Learning and Reasoning, ser. Proceedings of Machine Learning Research, B. Huang and M. Drton, Eds., vol. 275. PMLR, 07–09 May 2025, pp. 31–63. [9] B. Varici, K. Shanmugam, P. Sattigeri, and A. Tajer, “Causal bandits for linear structural equation models,” Journal of Machine Learning Research, vol. 24, no. 297, pp. 1–59, 2023. [10] Z. Yan and A. Tajer, “Linear causal bandits: Unknown graph and soft interventions,” in Advances in Neural Information Processing Systems, A. Globerson, L. Mackey, D. Belgrave, A. Fan, U. Paquet, J. Tomczak, and C. Zhang, Eds., vol. 37. Curran Associates, Inc., 2024, pp. 23 939– 23 987. [11] F. Fourati, C. J. Quinn, M.-S. Alouini, and V. Aggarwal, “Combinatorial stochastic-greedy bandit,” Proceedings of the AAAI Conference on Artificial Intelligence, vol. 38, no. 11, 2024. [12] I. Aouali, B. Kveton, and S. Katariya, “Mixed-effect thompson sampling,” in Proceedings of The 26th International Conference on Artificial Intelligence and Statistics, ser. Proceedings of Machine Learning Research, F. Ruiz, J. Dy, and J.-W. van de Meent, Eds., vol. 206. PMLR, 25–27 Apr 2023, pp. 2087–2115. [13] J. Zhou, P. Gaillard, T. Rahier, H. Zenati, and J. Arbel, “Towards efficient and optimal covariance-adaptive algorithms for combinatorial semi-bandits,” in Advances in Neural Information Processing Systems, A. Globerson, L. Mackey, D. Belgrave, A. Fan, U. Paquet, J. Tomczak, and C. Zhang, Eds., vol. 37. Curran Associates, Inc., 2024, pp. 30 474– 30 521. [14] T. Kocák and A. Carpentier, “Online learning with feedback graphs: The true shape of regret,” in Proceedings of the 40th International Conference on Machine Learning, ser. Proceedings of Machine Learning Research, A. Krause, E. Brunskill, K. Cho, B. Engelhardt, S. Sabato, and J. Scarlett, Eds., vol. 202. PMLR, 23–29 Jul 2023, pp. 17 260–17 282. [15] T. Lattimore and C. Szepesvári, Bandit Algorithms. Cambridge University Press, 08 2020. [16] A. Durand, C. Achilleos, D. Iacovides, K. Strati, G. D. Mitsis, and J. Pineau, “Contextual bandits for adapting treatment in a mouse model of de novo carcinogenesis,” in Proceedings of the 3rd Machine Learning for Healthcare Conference, ser. Proceedings of Machine Learning Research, F. Doshi-Velez, J. Fackler, K. Jung, D. Kale, R. Ranganath, B. Wallace, and J. Wiens, Eds., vol. 85. PMLR, 17–18 Aug 2018, pp. 67–82. [17] H. Huang, X. Ma, S. Erfani, and J. Bailey, “Neural architecture search via combinatorial multi-armed bandit,” in 2021 International Joint Conference on Neural Networks, 07 2021, pp. 1–8.

[18] G. B. Giannakis, Y. Shen, and G. V. Karanikolas, “Topology identification and learning over graphs: Accounting for nonlinearities and dynamics,” Proceedings of the IEEE, vol. 106, no. 5, pp. 787–807, 2018. [19] K. Königsberger, Analysis 2. Springer Berlin Heidelberg, 2004, p. 64 ff. [20] N. Merlis and S. Mannor, “Tight lower bounds for combinatorial multiarmed bandits,” in Proceedings of Thirty Third Conference on Learning Theory, ser. Proceedings of Machine Learning Research, J. Abernethy and S. Agarwal, Eds., vol. 125. PMLR, 09–12 Jul 2020, pp. 2830–2857. [21] K. Azuma, “Weighted sums of certain dependent random variables,” Tohoku Mathematical Journal, vol. 19, no. 3, pp. 357 – 367, 1967. [22] DB, “Api marketplace data,” https://github.com/Nokkyuu/Fantastic trains and when to find them, 2024, accessed: 2025-06-30. [23] B. Nourani-Koliji, S. Bilaj, A. Balef, and S. Maghsudi, PiecewiseStationary Combinatorial Semi-Bandit with Causally Related Rewards. IOS Press, 09 2023, pp. 1787–1794. [24] S. Ghoorchian, S. Bilaj, and S. Maghsudi, “Non-stationary delayed combinatorial semi-bandit with causally related rewards,” IEEE Open Journal of Signal Processing, vol. 6, pp. 369–384, 2025. [25] D. Russo and B. Van Roy, “Learning to optimize via informationdirected sampling,” in Advances in Neural Information Processing Systems, Z. Ghahramani, M. Welling, C. Cortes, N. Lawrence, and K. Weinberger, Eds., vol. 27. Curran Associates, Inc., 2014, pp. 1583– 1591. [26] Y. Shen, G. B. Giannakis, and B. Baingana, “Nonlinear structural vector autoregressive models with application to directed brain networks,” IEEE Transactions on Signal Processing, vol. 67, no. 20, pp. 5325–5339, 2019. [27] Y. Gai, B. Krishnamachari, and R. Jain, “Combinatorial network optimization with unknown variables: Multi-armed bandits with linear rewards and individual observations,” IEEE/ACM Transactions on Networking, vol. 20, no. 5, pp. 1466–1478, 2012. [28] D. Bertsimas and N. Koduri, “Data-Driven Optimization: A Reproducing Kernel Hilbert Space Approach,” Operations Research, vol. 70, no. 1, pp. 454–471, 2021.

9

TABLE I O RDINAL C OMPARISON OF R EGRET B OUNDS Work

Regret

Key Assumptions

[1]

O(N ln T ) √ O( N T ln T ) p O( T γT (γT + ln(1/δ))), where γT = O((ln T )N +1 ) for squared exponential kernels  q ˜ s}T , d˜max{d, Õ

Linear SEM, noiseless observations, DAG

where d˜ is the effective neural tangent dimension p O(dL−1 T ln(N T )) √ O(N dT log T ) √ Õ(βTL+1 (d + 1)L/2 p N T ), where βT = Õ(d ln(N T )) √ Õ(dL−0.5 T ) √ O(( N T ln T )2/3 ) √ Õ( N T ) √ Õ( N T )

Contexts encode arm dependencies, monotone score accumulation

[2] [3]

[4]

[5] [6] [9] [10] [11] [12] [13]

Oracle access, monotone reward function Expected scores observed with homogeneous noise across arms

Free selection of ingoing signal values, no action cardinality restriction, known DAG structure Monotonic generalized linear functions, endogenous interventional control, known DAG structure Linear SEM, DAG Linear SEM, known graph degree, DAG Submodular and monotonic expected reward Exemplary nonlinearities only for final accumulation, DAG No nonlinear transformation after action selection

2

[ours]

O(N √ ln (T )) for Alg. 1, O(N T ln(T )) for Alg. 2

Separable SEM, DAG for Alg. 1

A PPENDIX A. COMPARISON In Table I, we detail the time-dependent regret bounds from related papers. Alongside each, we list a few key assumptions in which they differ from the work conducted here. Some framworks defer the regret-impact of base arm number N to graph properties such as the maximal degree d and maximal path length L, which without further assumptions can be of order N . While frameworks and semi-bandit observations vary, intervened node values or prior contexts may be compared to filtered feedback by the volume of data available for learning. Compared to general SEMs, separable SEMs are less computational expansive for high-dimensional data, carry less risk of overfitting, and interactions associated to each edge are more easily interpretable. For nonlinear component functions, the inverse (I−F)−1 still captures more complex interactions, as it is in general not separable, nor are equilibrium solutions of cyclic SEMs. B. ADAPTATIONS 1) Extensions: For simplicity, we fixed G = I for the preliminary on-node activations. Generic exogenous transformations could be identified by optimization together with the entries of F, effectively leading to regret bounds as if another base arm were involved. Our proofs will still specify the step in which this identification would take place. Besides, we focused on zero-true F. Otherwise, non-zero N P F0 : [0, 1]N → RN , b 7→ ( F[i, j](0))⊤ i∈[N ] can be linearly j=1

transformed as I − F + F0 = (I + F0 )(I − F) towards a

zero-true operator F−F0 covered in our analysis. Uncertainty over the entries of F0 would have the same effect on the regret bounds as another exogenous operator since (I−F)−1 = (I − (F − F0 ))−1 (I + F0 ), assuming the inverses exist. Although we accumulated the final payoff from the overall rewards in a strictly linear manner, one could frame it through a function operator h of dimensions (N, 1), i.e., P rxt (bt ) = h⊤ yt = h[i](y[i]). If h is known, it is i∈[N ] straightforward to combine our theoretical results with known bounds for the expectation of the accumulating functions; Otherwise, their learning process would require the payoff as additional feedback and increase TB,Hθ ,N according to the complexity of the candidate function space in the regret bounds of Section IV. Remark A.1. i) Setting rx (b) = r(diag(x)b) implies that only the instantaneous rewards of the selected arms impact the reward. If the selection only altered the signal strength, e.g., rx (b) = r((diag(x) + kI)b), k ∈ R, we would need to treat s as N for most of our theoretical analysis but drastically reduce TB,Hθ ,N , as every base distribution would be involved no matter the agent’s decision. Separate sets of arms completely outside the agent’s control but not perception could be treated similarly. ii) In our setup, the agent selects exogenously, as its decision vector masks the independent instantaneous rewards and their observability directly. Were we to restate the problem by shifting the mask to the associated endogenous vectors, i.e., rxt (bt ) = 1⊤ (xt ⊙ yt ) = 1⊤ (xt ⊙ (I − F)−1 (G(bt ))), we could replicate the analyses of this paper with appropriate adjustments to the selection strategies and only minor changes

10

in theoretical variables. In either selection scenario, if both vectors’ observability is restricted to the chosen nodes, the external UCB should track instances of pairwise base arm occurrence, leading to a replacement of N with approximately N 2 in regret bounds. Partially masked feedback vectors could also be covered by the objective gap. iii) Alternative uniquely determined SEMs of the form y = F(y + z) require only adjustments of signal ranges as they are equivalent to y = (I − F)−1 (z) − z. Structural input to the agent at each time step, i.e. contexts, for example as exogenous addends, likewise would not affect ordinality. Remark A.2. i) Were we to forego the lower bounds to the distributions Bi , the continuity of their p.d.f.s would allow for an asymptotic argument in expected difference to a truncation from below. This argument would validate the regret bound of Theorem IV.3 with another invariant of B involved in the calculation of TB,Hθ ,N . Moreover, discrete distributions would require the RKHS to be of sufficiently low dimension to determine component functions from occupied values. ii) Ordinally, Theorem IV.6 applies to i.i.d. samples from white endogenous N2 instead of exogenous observation noise as well. In conjunction, bounds on the conditional denoising error E b∼B [|F(y(b)) − E b̃∼B [F(y(b̃)) | y(b̃) + ϵ̃ = ϵ∼N2

ϵ̃∼N2

y(b) + ϵ]|] can be utilized in the general regret analysis without imposing additional properties on the distributions or the functions. One could seamlessly employ the subroutines of [23] or [24] to better cope with structurally nonstationary environments, shifting distributions or feedback delays. Besides its expressiveness, our choice of causal model class mostly pertains to its sample efficient identification, and only the related parts would require major overhauls for different representations. In particular, general bounds can be derived from our results for any analytic objective function if known to the agent. Other signal transformations, such as graph neural networks, could be combined with matching confidence schemes (e.g., neural bandit [4]) to encompass an even wider range of problem setups categorically. On a DAG, the SEM already resembles the forward propagation of a neural network that outputs all its intermediate node values, with the instantaneous rewards acting as biases and function activations associated to edges rather than nodes. Though not studied here, ideally the external UCB would take the postulated graph adjacencies into account when evaluating the impact of optimization certainty improvements gained from inclusion of each base arm, reminiscent of probabilistic trigger settings [14]. From a Bayesian perspective, replacement with a more targeted parameter informationdirected sampling design [25] could aid the graph learning. 2) Specializations: Tailored to normal instantaneous reward distributions in particular, SSEM-UCB-Norm (Algorithm 3) decreases the number of moments that go into the internal UCB calculation to two regardless of q. In turn, it also needs to apply its confidence to the derivatives’ argument. The maximization over the confidence interval is well-defined due to the continuity of the underlying functions and their derivatives.

Algorithm 3 SSEM-UCB-Norm: Separable Structural Equation Model - Upper Confidence Bound - Normal Exactly as SSEM-UCB, but running the initialization loop twice, and replacing empirical moment calculation and selection criterion as follows. (2) 1: Calculate Ct−1 , β̂ t−1 = at and ϕ̂t−1 according to (7), t−1 [j] (8) and (6), respectively, and set m̃t−1 [j] = mm for t−1 [j]−1 j ∈ [N ]. 2: Select decision vector xt that solves arg max Et (x) + x∈X

X |α|≤q α mod 2=0

Y 

max Q Ct−1 ∈ i∈[N ] [−Ct−1 [i],Ct−1 [i]] Q (2) Wt−1 ∈ i∈[N ] {−Ct−1 [i],Ct−1 [i]}

Dα r̂(t−1) (x ⊙ (β̂ t−1 + Ct−1 )) α!! (2)

m̃t−1 [j]ϕ̂t−1 [j] +

p

(2)

m̃t−1 [j]Wt−1 [j]

α[j]/2

.

j∈I(x)

Algorithm 4 q-search Input: Time step t, input of SSEM-UCB except for q. ˆ min = 0. 1: Initialize q = −1, Γ = 0, ∆ ˆ min /4 and q < t do 2: while Γ ≥ ∆ 3: Increase q ← q + 1. ˆ min ← 4: Calculate ∆ min

Smax − S(x) for

x∈X Smax −S(x)>0

Smax ← max S(x′ ) based on the selection strategy ′ x ∈X

S and reward estimate r̂(t−1) of SSEM-UCB with parameter q. 5: Set Γ to a q-remainder bound of r̂(t−1) (either as provided in III-B or as an empirical estimate). 6: end while 7: return q

Since the distance of the deployed Taylor series approximation to the actual reward function is problem-dependent, the chosen cutoff degree q can be gradually increased over the course of the algorithm’s runtime until the diminishing analytical remainder falls within bounds proportional to the minimal empirical payoff discrepancy of suboptimal super arms, as outlined in Algorithm 4. C. PERFORMANCE Complexity. The computational expense of each time step t is split between kernel optimization, differentiation and arm comparison. The former two strongly depend on the method employed. For example, the Alternating Direction Method of Multipliers (ADMM) identifies the P best sample fit incrementally in O(t2 N ) for Alg. 2 and O( i∈[N ] xt [i]mt [i]2 N ) for Alg. 1 [26]. For at = 0, very basic numerical differentiation can be performed in O(w eval(F̂t−1 )), where w = (qN )min{s,q/2} and eval(F̂t−1 ) is the complexity of evaluating F̂t−1 . With all partial derivative values in memory, the

11

 Ps optimal combinatorial action amongst |X | = k=0 N k can always be picked in O(w + N s s2 ) by forming a computation tree of viable combinations. These most general bounds can be drastically decreased by graph sparsity (enforced through regularization or edge recognition thresholds), optimized block matrix inversion, refitting the graph on demand or in parallel, exploiting further structure or restrictions on the action space, randomizing challengers to the leading arm, and other adjustments described in Sections V and G. Failure modes. Since the external UCB ensures functional certainty through minimal base arm count, if all super arms containing a specific base arm perform exceptionally subpar, its cost might become unduly dominant, as is evidenced by the regret bounds linearity in ∆max . A successive arm elimination subroutine, where discrepancies in experienced performance are compared to a threshold for further inclusion, helps the algorithms to cope with this scenario. When dealing with very complex or extensive graphs, unless there are substantial risk factors, it is also advisable to spend additional iterations on the initial pure exploration phase to build up sufficient structural knowledge sooner. Further, for problems in which the operator I − F is singular, the solution to the SEM is not unique and structural learning instable. Specifying endogenous ranges, or secondary constraints such as endogenous norm minimization, can facilitate stable pseudo-inversion of the SEM while reducing the computational complexity of exact inversion. Greater detail on such generalizations and their integration into the algorithms is left to future work. D. COROLLARIES An alternative expression for the regret bound of Theorem IV.3 presents itself through Cauchy’s estimate and Hartogs’s theorem on the analytic continuations rτ , τ ∈ [T ], of the empiric reward functions to the complex numbers. For any l ∈ R>0 within the radius of converP (τ −1) |Dα rxτ (aτ )|/α! ≤ gence, we thereby have max τ ∈[T ] |α|≤q P sup |r(τ −1) |/l|α| , where Pl (aτ ) is the polydisc max τ ∈[T ] |α|≤q Pl (aτ )

of radius l in all directions around aτ in C. For example, for l = s + q − 1 we get e max sup |r(τ −1) | as a limit to τ ∈[T ] Ps+q−1 (aτ )

(q) the absolute sum. Multiplied by (∆min − 4Γmax )/8, this term (q,s) can be substituted for ψmin in Theorem IV.3.

For an approximate formulation of Theorem IV.3, we (q) (q) set ∆max = max ∆(x) and ∆min = (q)

x:µ(x)+4Γmax <µ(x∗ ) (q) min ∆(x) ≥ 4Γmax . (q) ∗ x:µ(x)+4Γmax <µ(x )

(i) Provided that F[i] ∈ Hθ ∀ i ∈ [N ], then,

Corollary A.3. for any degree q ∈ N, the expected regret of the SSEM-UCB algorithm is upper bounded as   2 2 4(s + 1) ln T π  +1+ (2qs + 1) N ∆(q) max (q,s) 2 3 ψ min

(q) + TB,Hθ ,N ∆(q) max + 4Γmax T.

Proof. See Section F3.

Next, we state the regret bound associated with the specialized Algorithm 3 in Section B2. Corollary A.4. In the setting of Theorem IV.3, if the instantaneous reward distributions Bi are all normal and all T samples bt lie in [0, 0.5]N , the expected regret of the SSEMUCB-Norm algorithm is upper bounded as 

2

2

 9(q − 1)!!q(s + 1) ln T + 2 + π (4s + 1) N ∆max (q,s) 2 3 ψmin + TB,Hθ ,N ∆max . □

Proof. See Section F4.

Introducing noise would lead to results analogous to Theorem IV.6. The benefits of less averages to track lie in memory usage as well as computation effort required to identify each super arm’s most optimistic conjecture. Their impact on performance would be more noticeable in Bayesian regret analysis than the absolute guarantees studied here, which necessarily involve many worst-case problem instance evaluations. A Bayesian view would likely improve all acquired bounds in their ordinality, especially for sparse graphs, which are the recommended application of SSEM-UCB due to the computation load associated with long paths and, to a lesser extent, high degrees. Corollary A.5. Let M ∈ [0, 1]N ×N be an invertible upper triangular matrix. In the idealized setting of Theorem IV.3, an initialization phase to SSEM-UCB that plays the ith column vector of M a number dim(H(i) ) times leads to a regret bound of 

" X

dim(H

(i)

)+

4(s + 1) ln T (q,s) 2

ψmin

i∈[N ]

(q)

#  2π 2 qs N  ∆max +1+ 3

(q,s)

without the term −4Γmax in ψmin . Proof. See Section F5.

In contrast to Theorem IV.6, the preceding corollary is valid only in the absence of noise or when the graph is known.

E. PROOF OF THEOREM IV.3 1) Notations: Before proceeding to the proof itself, we introduce some important notations together with their definitions. We define the index set of a decision vector x ∈ X as I(x) = {i | x[i] ̸= 0, ∀i ∈ [N ]}. The confidence q bound of (s+1) ln t base arm i at time t was declared as Ct [i] = mt [i] . At each time t, we collect the unbiased empirical estimates of (k) moments ϕ̂t [i] about the center point at+1 ∈ [0, 1]N and the calculated confidence bounds Ct [i] of all base arms i ∈ [N ] (k) (k) in vectors (ϕ̂t )k∈[q] and Ct , respectively. We write ϕt for the actual moments about at+1 of the instantaneous rewards.

12

For time step t + 1, with kp = 22 mod p , we set the internal selection index for a decision vector x ∈ X to max It+1 (x) =

max (p)

Wt

Q

(m)

(m)

},min{kp Ct [i],1−ϕ̂t

}}

(α[j]) (α[j]) ϕ̂t [j] + Wt [j]



{max{−kp Ct [i],1−kp −ϕ̂t

Lemma A.6. [21] Let z1 , z2 , . . . , zm be random variables, c, d ∈ R and zi ∈ [c, d], ∀i. Moreover, E[zt |z1 , . . . , zt−1 ] = γ, for all t = 1, . . . , m. Then, for all D ≥ 0, # " m 2 X − 2D zi − mγ ≥ D ≤ 2e m(d−c)2 . P (10) i=1

i∈[N ]

X Dα r̂(t) Y  x (at+1 ) |α|≤q

α!

j∈I(x) (k)

min and It+1 (x) analogously. Here, we w.l.o.g. restricted ϕ̂t [i]+ (k) Wt [i] to the range [−1, 1] of the actual moments. Let (p) (p) Wt−1max (x) and −Wt−1min (x), p ∈ [q], be confidence scores that realize Itmax (x) and Itmin (x) respectively. We call p Et+1 (x) = ∆It+1 max x[i]4Ct [i] ln(t) i∈[N ]

the (weighted) external selection index. At time t, the algorithm selects a decision vector xt that maximizes Itmax + Et . For ease of presentation, in the sequel, we use the equivalence 1⊤ (I − F̂t−1 )−1 Ĝt−1 (diag(b)xt ) = 1⊤ (I − F̂t−1 )−1 Ĝt−1 diag(xt · id)(b) for all b ∈ [0, 1]N . To simplify the notation, sometimes we drop the time index t in mt [i] and use m[i] to denote the number of times that the base arm i has been observed up to the current time instance. Since the choice of center point for the Taylor series expansion at each time step is free, we also use a and ϕ(k) in lieu of at+1 (k) and ϕt , respectively. For any x ∈ X , we use the counter Tx (t) to represent the total number of times the decision vector x has been selected up to time t. Finally, for each base arm i ∈ [N ], we define a counter Ti (t) which is updated as follows. At each time t > TB,Hθ ,N , after a starting point TB,Hθ ,N ∈ N to be specified later on, that a suboptimal decision vector xt is selected, we have at least one base arm i ∈ [N ] such that i = arg min mt [i]. In this case, if the base arm i is unique, i∈I(xt )

we increment Ti (t) by 1. If there is more than one such base arm, we break the tie and select one of them arbitrarily to increment its corresponding counter. For a distribution D on a Lebesgue measurable space (S, Σ, µ) and analytic functions h, g : S → R, we call p dD (h, g) = Ed∼D [(h(d) − g(d))2 ] sZ = p.d.f.D (d)(h − g)(d)2 dd S

the distance between h and g respective D. We note that on finite dimensional spaces of analytic functions and for µ(supp(D)) > 0, the corresponding weighted L2 -seminorm is a true norm, because the intersection of distinct analytic functions has measure 0. In principle, this would allow for norm equivalence arguments in the optimization part of the proof. Since we need to track dependencies on T or N however, we will state equivalence factors explicitly. 2) Auxiliary Results : We use the following lemmata in the proof of Theorem 1.

Lemma A.7. Let H ⊆ L2 (S, µ) be a vector space of analytic functions and h ∈ H map from the set of a Lebesgue measure space (S, Σ, µ) to R. Let B1 and B2 be distributions on subsets S1 and S2 of S, respectively. Further, let lB1 be a lower bound to the existing p.d.f. of B1 , µ(S1 ) > 0 and (e1 , . . . , el ) be an orthonormal basis of H relative S1 . Then Pl p EB2 [h2 ] ≤ EB1 [h2 ]( i=1 EB2 [e2i ])2 /lB1 holds true. Pl Proof. For the basis representation h = i=1 ki ei we have !2 Z l X 2 (a) EB1 [h(b) ] = ki ei (b) p.d.f.B1 (b)db S1

i=1 l X

Z ≥ lB1 S1 (c)

= lB1

!2 ki ei (b)

db

i=1

l X

ki2 ≥ lB1 max ki2 , i∈[l]

i=1

(11)

where (c) is due to the orthonormality of the basis over S1 and (a) obeys the law of the unconscious statistician. This bound to the maximum norm allows us to conclude p

l (a) X

EB2 [h(b)2 ] ≤

p |ki | EB2 [ei (b)2 ]

i=1 (11)

s

1 lB1

EB1 [h(b)2 ]

l p X

EB2 [ei (b)2 ],

i=1

with (a) derived from the Minkowski inequality. □ Lastly, we note that, since summation and concatenation preserve the analytic property, it is maintained by the vector functions I − F̂t on [0, 1]N . Their Jacobians are upper triangular matrices with ones on their diagonal and thus of full rank. By Lagrange’s inversion theorem, the analytic property propagates to the inverses (I − F̂t )−1 and hence the empiric (t) reward functions r̂x = 1⊤ (I − F̂t )−1 Ĝt diag(x · id) as well. (t) In particular, each partial derivative Dα r̂xt+1 with α ∈ NN is uniformly continuous on the compact N -dimensional interQ val i∈[N ] [−Ct [i], 1 + Ct [i]], a fact we will use implicitly throughout the proof. 3) Proof: The proof follows the same basic outline as in [1], but the analysis of uncertainties and nonlinearities requires substantial extensions and modifications. Here, we will focus on SSEM-UCB, and defer the changes required for SSEMUCB-JO to Section F1. We start by rewriting the expected regret as R(T ) = T µ(x∗ ) −

T X

µ(xt )

t=1

=

X x:µ(x)<µ(x∗ )

∆(x)E[Tx (T )].

(12)

13

Based on the definition of the counters Ti (t) for the base arms i ∈ [N ], at each time t that a suboptimal decision vector is selected, only one of such counters is incremented by 1. Thus, we have [27]   "N # X X E Tx (t) ≤ TB,Hθ ,N + E Ti (t) , x:µ(x)<µ(x∗ )

i=1

which implies that X

N X

E [Tx (t)] ≤ TB,Hθ ,N +

x:µ(x)<µ(x∗ )

⊤ Let vt+1 = 1⊤ (I − F̂t )−1 Ĝt . We order the elements in sets ∗ I(x ) and I(xt+1 ) arbitrarily. In the following, our results are independent of the way we order these sets. Let vk , k = 1, . . . , |I(x∗ )| ≤ s, represent the kth element in I(x∗ ) and uk , k = 1, . . . , |I(xt+1 )| ≤ s, represent the kth element in I(xt+1 ). Explicitly, we have ( T X Ti (T ) ≤ l + 1 min 0<m[v1 ],...,m[v|I(x∗ )| ]≤t

t=TB,Hθ ,N

X

E [Ti (t)] .

i=1

|α|≤q x∗ ⊙α=α

Therefore, we observe that X R(T ) = ∆(x)E[Tx (T )]

⊤ Dα vt+1 diag(x∗ · id) (at+1 )) α!

|I(x∗ )| (α[vj ])

Y

(ϕ̂t

(α[j])max

[vj ] + Wt

(x∗ )[vj ]) + Et+1 (x∗ )

j=1

x:µ(x)<µ(x∗ ) (∗)

≤ ∆max TB,Hθ ,N + ∆max

N X

≤ Et+1 (xt+1 ) + E[Ti (T )],

(13)

X

max l≤m[u1 ],...,m[u|I(xt+1 )| ]≤t

i=1

⊤ Dα vt+1 diag(xt+1 · id) (at+1 ) α! ) |I(xt+1 )| Y (α[uj ]) (α[j]) max (ϕ̂t [uj ] + Wt (xt+1 )[uj ])

where (∗) follows from the definition of ∆max . Let Ii (t) denote the indicator function which is equal to 1 if Ti (t) is increased by 1 at time t, and is 0 otherwise. Consequently, T X

Ti (T ) =

j=1

1 {Ii (t) = 1} .

If Ii (t) = 1, it means that a suboptimal decision vector xt is selected at time t. In this case, mt [i] = min {mt [j]|j ∈ I(xt )}. Let −1) max |Dα r̂(τ (aτ )|, xτ

τ ∈[T ] |α|≤q+1

ψmin = min

 wmax

min{q,s} , 1 (q+1)smin{q,s}+1 ∆I

∆

∞ X

≤l+

t=TB,Hθ ,N +1

wmax = max

|α|≤q xt+1 ⊙α=α

(q) min − 4Γmax

8

t X

t=TB,Hθ ,N mv1 =1 ⊤ Dα vt+1

(

t X

t X

...

mv|I(x∗ )| =1 mu1 =l mu|I(x ∗ |I(x P )|

t+1 )|

=l

! (at+1 [vj ])I[vj ]⊤

j=1

X

1

t X

...

α!

|α|≤q x∗ ⊙α=α |I(x∗ )|

Y

(α[vj ])

(ϕ̂t

(α[j]) max

[vj ] + Wt

(x∗ )[vj ]) + Et+1 (x∗ )

j=1

and

4(s + 1) ln2 T l= 2 ψmin 

 .

Ti (T ) =

T X

1 {Ii (t) = 1}

T X

1 {Ii (t) = 1 & Ti (t − 1) ≥ l}

≤l+

1{Itmax (x∗ ) + Et (x∗ )

t=TB,Hθ ,N +1

≤ Itmax (xt ) + Et (xt ) & Ti (t − 1) ≥ l} =l+

T X

α!

|I(xt+1 )|

(α[uj ])

(ϕ̂t

(α[j]) max

[uj ] + Wt

(xt+1 )[uj ])

j=1

) + Et+1 (xt+1 ) .

t=TB,Hθ ,N +1 T X

(at+1 [uj ])I[uj ]

j=1

|α|≤q xt+1 ⊙α=α

Y

t=TB,Hθ ,N +1

!

|I(x t+1 )| P

X

Then,

≤l+

⊤ Dα vt+1

max ∗ 1{It+1 (x ) + Et+1 (x∗ )

t=TB,Hθ ,N max ≤ It+1 (xt+1 ) + Et+1 (xt+1 ) & Ti (t) ≥ l}.

Based on the definition of Ti (t), we have Ti (t) ≤ mt [i], ∀i ∈ [N ]. Hence, when Ti (t) ≥ l, we know that [27] l ≤ Ti (t) ≤ mt [j],

∀j ∈ I(xt+1 ).

(14)

Given mv1 , . . . , mv|I(x∗ )| , mu1 , . . . , mu|I(xt+1 )| , we define event P as max ∗ max It+1 (x ) + Et+1 (x∗ ) ≤ It+1 (xt+1 ) + Et+1 (xt+1 ).

Event P would necessitate at least one of the following events to apply: max ∗ ∗ It+1 (x ) + ∆min /4 + Γ(q) max < µ(x ),

(15)

min It+1 (xt+1 ) − ∆min /4 − Γ(q) max > µ(xt+1 ), ∗ max min µ(x ) − µ(xt+1 ) ≤ It+1 (xt+1 ) − It+1 (xt+1 ) (q) + Et+1 (xt+1 ) + ∆min /2 + 2Γmax .

(16) (17)

14

P(Z = z|t̂ ∼ Uni({t ∈ [t]|o ∈ I(xt )}, z ∼ B xt̂ ). For i ∈ [N ] we write

First, we consider (15). We set X Dα r̂(t−1) (at+1 ) Y  x

Jt (x) =

α!

|α|≤q

ϕ(α[j]) [j]



(i) h∈Hθ

j∈I(x)

and investigate the event V

max ∗ It+1 (x ) < Jt+1 (x∗ ), | {z }

ϱt ≥ P(EB(o) [(F̂t (y(b))[i] − F(y(b))[i])2 ]

∗ Jt+1 (x ) + ∆min /4 + Γ(q) max < µ(x ) .

{z

(o)

− EB(o) [(F(y(b))[i] − FB(o) (y(b))[i])2 ] > Ht ) t

Therefore, we have

t

with FB(o) [i] = arg min dB(o) (hy, Fy[i]), ϱt = 1 − t

P [V] ≤ P [I] + P [II] .

(19)

(m) max ∗ For event I to be true, by selection of Wt (x ), there ∗

would have to be k ∈ [|I(x )|] and m ∈ [q] with (m)

22 mod m Ct [vk ] < ϕ(m) [vk ] − ϕ̂t

[vk ] .

(20)

Considering each case, we have   (m) P 22 mod m Ct [vk ] < ϕ(m) [vk ] − ϕ̂t [vk ]  = P 22 mod m mt [vk ]Ct [vk ]  (m) < mt [vk ] ϕ(m) [vk ] − ϕ̂t [vk ] (a)

2

≤ 2e−(2/mt [vk ])mt [vk ] Ct [vk ]

= 2t

t

1 − t−2(s+1) , s

ln(2/ϱt ) , 2mt [o] mt [o] √ √ √ where M1 = 4(κ θ + max r(b))κ θ, M2 = 2(κ θ + N r b∈[0,1] max r(b))2 , κ = sup K(b, b) and K is the positive (o)

Ht

= M1 p

1

+ M2

b∈[0,1]

i∈[N ]F̂[i]∈H(i) x∈X θ

−2(s+1) ln t

,

(21)

(22)

Now, we investigate Event II:   P |µ(x∗ ) − Jt+1 (x∗ )| > ∆min /4 + Γ(q) max   (a) (t) ≤ P E[rx∗ (b) − r̂x∗ (b)] > ∆min /4   2 (t) 2 = P E[rx∗ (b) − r̂x∗ (b)] > ∆min /16 ≤P

(i)

h∈Hθ

semi-definite kernel associated with the RKHS H(i) . Since K is continuous, a finite supremum exists over the closed interval domain of functions in H(i) . Even though the objective gap q φHθ = max min max E[(F − F̂)(yx (b))[i]2 ]

2

P [I] ≤ q|I(x∗ )|2t−2(s+1) ≤ 2qst−2(s+1) .



√ N6

b∈[0,1]N

where (a) follows from Lemma A.6, and (b) from the definition of Ct . Since this bound applies to each case’s probability, for Event I, we conclude that

(b)

(24)

t

}

II

−2(s+1)

τ =1

(o)

I

= 2e

(h(y(zτ )) − y(zτ )[i] + zτ [i])2

for the solution to the minimization problem with norm constraint θ (see Algorithm 1). A Rademacher complexity analysis of our kernel optimization (as described in Prop. 2 of [28]) results in an accuracy bound of

(b)

t X

(i) h∈Hθ

If V is true, then at least one of the following must hold.

(h(y(zτ )) − F(y(zτ ))[i])2

τ =1

= arg min

max ∗ It+1 (x ) + Jt+1 (x∗ ) + ∆min /4 + Γ(q) max < Jt+1 (x∗ ) + µ(x∗ ).

|

t X

(o)

F̂t [i] = arg min

(18)

(t) E[(rx∗ (b) − r̂x∗ (b))2 ] > ∆2min /16 (q) Γmax ,



,

is zero as by the assumption of Theorem IV.3, we treat it as a variable to support the generalizations in Theorem IV.6. For any x we can represent B x as the distribution induced by B over [0, 1]x = [0, 1]|I(x)∩{i}| , which is nonzero

×

i∈[N ] ⊤

almost everywhere. In particular, B˜o := B I[o] assigns 0 in every dimension other than the oth. Taking into account that dB(o) (Fy[i], FB(o) y[i]) ≤ dB(o) (Fy[i], FHθ y[i]) ≤ φHθ for a t t t realization FHθ of φHθ , we deduce that (24)

ϱt ≥ P(

min

(o)

EBx [(F̂t (y(b))[i] − F(y(b))[i])2 ]

x∈X :o∈I(x) (o)

> Ht

+ φ2Hθ )

(o)

(23)

where (a) results from the definition of and (b) is implied by Jensen’s inequality. Next, we will derive an ordinal bound to this probability from our optimization guarantees. We interpret the exogenous vectors observed at time steps before t when base action o ∈ [N ] was part of the super arm (o) as sampled from the empirical distribution Bt , defined by

≥ P(EB̃o [(F̂t (y(b))[i] − F(y(b))[i])2 ] uB (o) > s (Ht + φ2Hθ )), lB

(25)

where uB , lB > 0q are an upper respective lower bound on Pl the factors ( i=1 EB̃o [e2i ])2 , ly(Bx ) in Lemma A.7 for all distributions y(B x ), x ∈ X , o ∈ [N ]. Had we not fixed the exogenous operator G = I, we would jointly optimize for its prediction alongside F and directly gain a prediction certainty for G[o, o] regarding B̃o , which we would utilize analogously in subsequent steps.

15

Since F is not cyclic, the SEM describes the recursive formula of backward substitution: N X  −1 (I − F) (b)[i] = b[i] + fij (I − F)−1 (b)[j] , (26) j=i+1

(o)

dB̃o (ŷt [i], y[i])

(27)

(o) (o) = dB̃o (F̂t ŷt [i], Fy[i]) (o) (o) (o) (o) ≤ dB̃o (F̂t y[i], Fy[i]) + dB̃o (F̂t ŷt [i], F̂t y[i]) o−1 X (a) (o) (o) ≤ dB̃o (F̂t y[i], Fy[i]) + η dB̃o (ŷt [j], y[j]) j=i+1 (o)

= dB̃o (F̂t y[i], Fy[i]) + η

=

o−2−i X

(η + 1)k

k=0

FB̃o )−1 ,

(I =

arg min (i)

FB̃o [i]

=

dB̃o (hy, Fy[i]),

h∈Hθ

K(x, y)

=

Pdim(F (j) ) p

max

i=1

j∈[N ]

where i = N, . . . , 1, b ∼ B. As such, we can express the (o) (o) expectation for ŷt = (I − F̂t )−1 with

(b)

yB̃o (B̃o ) (fi,j )j∈[N ]

Here,

EBx [ei (y(b))2 ]/lB , and

(e1 , . . . , edim(F (j) ) ) is a basis of the space F (j) containing the algorithm’s possible values for the entries of F̂t that is orthonormal over [0, 1] in accordance to Lemma A.7. (For SSEM-UCB-JO, F (j) = H(j) .) Further, g1 (η) with g1 : R≥0 → R≥0 exists as a gradient bound on each component function. A numerical value of this bound can be determined by replicating the transformation in (27) (o) (o) for |ŷt (bI[o]⊤ )[j] − ŷt (cI[o]⊤ [j])| ∀j ∈ [N ], yielding N N g1 (η) = (1 + η(1 + η) ) ). (For SSEM-UCB-JO, g1 (η) = η would suffice as each component function lies in one of the RKHSs.) Replicating (27) also reveals

(o)

dB̃o (F̂t y[i + 1 + k], Fy[i + 1 + k]), where (b) transits to the explicit form of the recursive formula (o) and (a) utilizes the continuity of F̂t [i]. Here, η is an upper bound on all the B̃o , o ∈ [N ], associated distances of the gradient for all the bounded functions in Hθ , which exists since derivation constitutes a linear map. This gives rise to (o)

P(EB̃o [(ŷt (b)[i] − y(b)[i])2 ] (25) uB

>

(o)

s (Ht

lB

+ φ2Hθ )(1 + η

o−2−i X

(η + 1)k ))2

k=0

≤ 1 − (1 − ϱt )o−1 .

(28)

The prediction accuracies over the separated B̃o ’s will provide an accuracy bound for each component function of (t) F̂t = (fˆij )i,j∈[N ] , as defined in Algorithm 1, by consulting (t) (26) again. To this end, batches (fˆo−k,o )o∈{k+1,...,N } can be assessed in chronological order k = 1, . . . , N − 1 as (t) (B̃o ) dBo (fˆo−k,o , fo−k,o )   (o) ≤ dB̃o ŷt [o − k], yB̃o [o − k]

+

  (o) dB̃o ŷt [j], yB̃o [j]    (o) ≤ dB̃o ŷt [j], y[j] + dB̃o yB̃o [j], y[j]   (o) ≤ dB̃o ŷt [j], y[j] + g2 (φHθ )

o−1 X

dB̃o



with g2 (φHθ ) = φHθ (1 + η

o−2−j P

recursive form of (29), we then find that

(29)

(t) (o) (B̃o ) fˆo−k,j (ŷt [j]), fo−k,j (yB̃o [j])



 P EBo

(t) (B̃o ) fˆo−k,o (b) − fo−k,o (b)

j=o−k+1

≤ dB̃o +





(o) ŷt [o − k], yB̃o [o − k]

o−1 X

> (1 + K(I[o] , yB̃o )

+

k−2 X

2 

(K(I[o]⊤ , yB̃o ) + 1)j )2

j=0

  (t) (B̃o ) dB̃o fˆo−k,j (yB̃o [j]), fo−k,j (yB̃o [j])

2

(1 + (k − 1)g1 (η)) (

j=o−k+1 o−1 X

(η + 1)k ). Resolving the

k=0

  (t) (o) (t) dB̃o fˆo−k,j (ŷt [j]), fˆo−k,j (yB̃o [j])

q

! (o) Lt + g2 (φHθ ))2

(28)

≤ 1 − (1 − ϱt )(o−1)(k−1)k ,

(30)

j=o−k+1

  (o) ≤ dB̃o ŷt [o − k], yB̃o [o − k] +

o−1 X

  (t) (B̃o ) K(I[o]⊤ , yB̃o )dBj fˆo−k,j , fo−k,j

j=o−k+1

+

o−1 X j=o−k+1

(o)

abbreviating Lt 

(o)



g1 (η)dB̃o ŷt [j], yB̃o [j] .

(o,k)

(o)

= ulsB (Ht B

+ φ2Hθ )(1 + η

o−2−i P

(η + 1)k ))2

k=0

and Mt for the lower side of the inequality in the above probability.

16

(o)

[N ]. The negative derivative of Ht indicates that Gt is monotonically falling as well. Thus when φHθ is sufficiently small, there is TB,Hθ ,N ∈ N such that for t > TB,Hθ ,N ,

We can now show that dBx∗ (ŷt [i], y[i]) ≤

N X

(t)

dBx∗ (fˆij (ŷt [j]), fij (y[j]))

j=i+1

N  X

(BI[j] )

(t)

dBx∗ (fˆij (y[j]), fij

4Gt < 4GTB,Hθ ,N < ∆min .

(y[j]))

j=i+1 ⊤

(BI[j] )

+ dBx∗ (fij (y[j]), fij +

N X

(y[j]))

 To subsume the technical initialization phase of the algorithm, in which a first sample of each instantaneous reward was collected, we also demand TB,Hθ ,N ≥ ⌈N/s⌉.

(t) (t) dBx∗ (fˆij (ŷt [j]), fˆij (y[j]))

j=i+1

N X

  (t) K(x∗ , y) dBj (fˆij , fij ) + g3 (φHθ )

Overall, together with (23) and (22), we find that

j=i+1

+

N X

g1 (η)dBx∗ (ŷt [j], y[j]).

P[V] ≤ (2qs + 1)t−2(s+1) .

j=i+1

The numerical value of g3 (φHθ ) with zero-true g3 : R≥0 → (o,k) R≥0 can be stated along the lines of (29) and similar to Mt , which we forego here. Finally our analysis yields   P EBx∗ [(r̂(t) (b) − r(b))2 ] > G2t (30)

N6

≤ 1 − (1 − ϱt )

−2(s+1)

=t

(32)

The same upper bound can be derived for the probability of (16) in a similar fashion. Finally, we consider (17). For every m ∈ [q], j ∈ I(xt+1 ), we have

for Gt =(1 + g1 (η)

N −1 X

√ 3 (g1 (η) + 1)k ) N

k=0 r  (o,k) K(x , y) max Mt + g3 (φHθ ) . ∗

o,k∈[N ]



(31)

From our selection criterion we know that v u u 16(s + 1) ln2 (t) max 0 ≤ It+1 (xt+1 ) + ∆It+1 t min mt [i] i∈I(xt+1 ) v u u 16(s + 1) ln2 (t) max − It+1 (I[arg min mt [i]]⊤ ) − ∆It+1 t min mt [i] i∈[N ] i∈[N ]   

(m)

 (m) [j] + Wt max (xt+1 )[j]   (m) (m) − ϕ̂t [j] − Wt min (xt+1 )[j] s r (s + 1) ln t (s + 1) ln t (a) ≤4 ≤ 4Ct [j] = 4 mt [j] l

ϕ̂t

(q) min − 4Γmax ) , (q + 1)smin{q,s}+1

(b) min{q, s}(∆

4wmax

(33)

where in (a) and (b), we substituted the value for Ct [j] and l, respectively. On the same note, we find that

  1  √ 1   , ≤ ∆It+1  1 + 4 s + 1 ln(t)  √l − r  min mt [i] i∈[N ]

v

and, as ties are broken randomly, w.l.o.g. ∆It+1 > 0, hence −2 1 min mt [i] ≥  q +√  i∈[N ] 2 l 16(s + 1) ln (t) 

u u 16(s + 1) ln2 (t) Et+1 (xt+1 ) = ∆It+1 t min mt [i] i∈I(xt+1 )

1

≤ ∆It+1

4(s + 1) ln2 (t) . (1/2 + ψmin )2

16(s + 1) ln2 (T ) l

≤ 2∆It+1 ψmin ≤ ∆min /4 − Γ(q) max .

Through appropriate change of variable and application of L’Hôpital’s rule, it is straightforward to verify that (o) lim ln(ϱτ )/ ln2 (τ ) = 0 and therefore lim Hτ = 0 ∀o ∈

τ →∞

s

τ →∞

(α)

Setting Zmax [j] =



(α[j])

ϕ̂t

(α[j]) max

[j] + Wt

(xt+1 )[j]

(34)



and

17

(α)

Zmin [j] accordingly, we can conclude max min It+1 (xt+1 ) − It+1 (xt+1 ) + Et+1 (xt+1 ) (t) X Dα r̂xt+1 (at+1 )

=

α!

|α|≤q



Y j∈I(xt+1 )

(c)

j∈I(xt+1 )

X |Dα r̂(t) xt+1 (at+1 )| α!

|α|≤q

 (α) Zmin [j] + Et+1 (xt+1 )

Y

(α) Zmax [j] −

(α)

(α) Zmax − Zmin

1

+ Et+1 (xt+1 ) (q) X (33) wmax min{q, s}(∆min − 4Γmax ) < s α! 4wmax (q + 1)smin{q,s}+1 |α|≤q xt+1 ⊙α=α

+ Et+1 (xt+1 ) q (q) X sk min − 4Γmax ) ≤ + Et+1 (xt+1 ) k! 4(q + 1)smin{q,s} k=0 ≤ ∆min /4 − Γ(q) max + Et+1 (xt+1 ) (34) ≤ ∆t+1 − ∆min /2 − 2Γ(q) max = µ(x∗ ) − µ(xt+1 ) − ∆min /2 − 2Γ(q) max ,

(d) min{q, s}(∆

(35)

where (d) is a combinatorial result from induction over s and application of the binomial theorem. Moreover, (c) (α) (α) results from Zmax [j], Zmin [j] ∈ [−1, 1] by decomposing the difference of products into a sum of products of differences. Hence, we conclude that (17) never happens. By using (14), (32) and (35), we achieve   4(s + 1) ln2 T E[Ti (T )] ≤ 2 ψmin t t t ∞ t X X X X X ··· ··· + t=TB,Hθ ,N mv1 =1

2(2qs + 1)t

mvs =1 mu1 =l

4) Dependencies: i) The maximum norm in regards to an orthonormal basis of H(i) (w.r.t. the inner product on (i) H(i) ), i ∈ [N ], is bounded by θ for functions from Hθ similar to Lemma A.7. Maximum norms in other bases are therefore bounded on a combination of θ and dim(H(i) ). The number of nonzero coefficients in the basis representation is proportional to the highest degree among the empirical graphs corresponding to (F̂t )t∈[T ] . For a polynomial Hilbert space, the coefficients of any r̂(t) are then limited by a power of this quantity determined by q, as are the linearly mapped partial derivatives of r̂(t) subsumed by wmax . Instead of q, for any RKHS, an upper bound can be stated in terms of the maximum path length among the empirical graphs, which equals the longest function chain in the inverse (I − F̂t )−1 . The latter should be combined with a norm-threshold for edge recognition in the algorithm. Both upper limits on wmax are independent of T and grow with the density of the graphs. ii) Similar to i), the gradient bound η introduced in (27) depends on θ as well as the maximum path length p and maximum degree d of the empirical graphs belonging to (F̂t )t∈[T ] . The number of relevant indices in all summations in (27)-(31) is at most d. We can then find TB,Hθ ,p,d ∈ N independent of N for which (32) is replaced by P[V] ≤ (2qs + N )t−2(s+1) to accommodate the summation of the N overall rewards. F. REMAINING PROOFS 1) Proof of Theorem IV.3 for SSEM-UCB-JO: We only need to update the kernel optimization guarantees utilized in proof E. To this end, we √ use the same variable terms with l changed by a factor of T / ln(T ). We write ŷt+1 = (I − F̂t+1 )−1 . We understand the exogenous vectors observed up to time t as sampled from the empirical distribution Bt . Optimizing over this entire data set results in a complexity bound of ϱt ≥ P(EBt [(F̂t (y(b))[i] − F(y(b))[i])2 ] − EBt [(F(y(b))[i] − FBt (y(b))[i])2 ] > Ht ) for every i ∈ [N ], with FBt [i] = arg min dBt (hy, Fy[i]) and

mus =l

(i)

h∈Hθ

−2(s+1)

r

∞ X 4(s + 1) ln2 T + 1 + (2qs + 1) 2t−2 2 ψmin t=1

4(s + 1) ln2 T π2 +1+ (2qs + 1). 2 ψmin 3

(36)

1 ln(2/ϱt ) Ht = M1 √ + M2 2t t as opposed to (24). In combination with the definition of FBt , for every o ∈ [N ] we know that ϱt ≥ P(EBt [(F̂t (y(b))[i] − F(y(b))[i])2 ]

Therefore, the expected regret is upper bounded as

− EBt [(F(y(b))[i] − FB(¬o) (y(b))[i])2 ] > Ht ) t

R(T ) ≤ TB,Hθ ,N ∆max + ∆max

N X

≥ P(EB(o) [(F̂t (y(b))[i] − F(y(b))[i])2 ]

E[Ti (T )]

t

− EB(o) [(F(y(b))[i] − FB(¬o) (y(b))[i])2 ] >

i=1

≤ TB,Hθ ,N ∆max  N  X 4(s + 1) ln2 T π2 + + 1 + (2qs + 1) ∆max 2 ψmin 3 i=1

t

t

t Ht ) mt [o]

≥ P(EB(o) [(F̂t (y(b))[i] − F(y(b))[i])2 ] t

>

≤ TB,Hθ ,N ∆max   4(s + 1) ln2 T π2 +1+ (2qs + 1) N ∆max . + 2 ψmin 3

t Ht + φ2Hθ ), mt [o]

(¬o)

with Bt formed from samples in time steps where o was not selected, and altered objective gap q φHθ = max max E[(F − FBx′ )(yx (b))[i]2 ]. ′ i∈[N ]x,x ∈X

18

Following the proof of Theorem IV.3, we only need to establish that lim mττ [o] Hτ = 0 to be able to infer a sufficient τ →∞ TB,Hθ ,N ∈ N accordingly. From the algorithm’s selection √ 4(s+1) t ln t criterion, we see that min mt [i] ≥ (1/2+ψmin )2 . Combined i∈[N ] p √ with lim τ ln(ϱτ )/( τ ln(τ )) = 0, the result follows. ■ τ →∞ 2) Proof of Proposition IV.4: For a cyclic SEM, convergence of F̂t towards F on the sample measure EBt implies convergence in an associated weighted L2 -norm since the operators are analytic. As all H(i) , i ∈ [N ], are finitedimensional, all norms on them are equivalent and the operators in particular converge with respect to their Hilbert space norm. For any RKHS, the convergence must be uniform, as |(F̂t − F)(y)[i]| = |⟨F̂t [i] − F[i], Ky(i) ⟩H(i) | ≤ ∥Ky(i) ∥H(i) ∥F̂t [i] − F[i]∥H(i) (i)

by Cauchy-Schwarz, where Ky = K (i) (y, ·) is the kernel function of H(i) and y is a solution to the SEM of F for any b ∈ [0, 1]N . The inverse function theorem then implies convergence of the inverses (I − F̂t )−1 at the same timerate, although with different constant factors affecting only THθ ,B,N . Past the operator accuracy, the remainder of proof E does not exploit acyclicity and can be replicated here. ■ 3) Proof of Corollary A.3: We have R(T ) =

T X

∆(xt )

t=1

X

=

∆(xt ) +

t∈[T ] 4Γ(q) max <∆(xt )

 (Thm. IV.3)

2  4(s + 1) ln T (q,s) 2

ψmin

X

∆(xt )

t∈[T ] 4Γ(q) max ≥∆(xt )

 π2 +1+ (2qs + 1) N ∆(q) max 3

(q) + TB,Hθ ,N ∆(q) max + 4Γmax T.

■ 4) Proof of Corollary A.4: All central moments ϕ(p) [j] = EBj [(b − EBj [c])p ] of normal distributions can be calculated from their respective variance: ∀p ∈ N, j ∈ [N ] : ϕ(2p) [j] = (2p − 1)!!ϕ(2) [j]p and ϕ(2p−1) [j] = 0. Replicating the proof of Theorem IV.3, only the unbiased estimates of mean and variance need to be taken into account in (20). Further, (33) needs to be extended by the argument that for |a|, |b|, |c| ≤ 0.5, p ∈ N, by induction we have |(b − a)p − (c − a)p | ≤ p|(b − a) − (c − a)| = p|b − c| (t)

and that all Dα r̂xt−1 for |α| ≤ q, t ∈ [T ], are swmax -Lipschitz continuous. ■ 5) Proof of Corollary A.5: The dim(H(i) ) samples obtained while playing column vector M[i] are with probability 1 in general position and allow for exact identification of F[i] for any i ∈ [N ]. Thus, the terms inherited from the optimization uncertainty after the initialization phase can be omitted from the bound established in Theorem IV.3. ■

(m)

6) Proof of Theorem IV.6: We set ϕ̃ = Eb∼B,ϵ∼N [(b + ϵ − at+1 )m ] for the noisy moments about at+1 of the exogenous signals. Only definition (18) (and its subsequent appearances) in the proof of Theorem IV.3 has to be altered to feature the noisy moments instead of the pure signal moments. We note that the expected model noise can be incorporated into the endogenous operator since F + E[M] ∈ Hθ . The expectation in the objective gap shall be taken with respect to the noise as well. ■ G. EXPERIMENTAL ADDENDUM 1) Synthetic: As a measure of nonlinearity, we only accepted problem instances in which the true reward function’s first degree Taylor approximation’s remainder empirically averaged above 10% of the maximal reward in absolute value. Concerning IGP-UCB, the underlying RKHS is defined by a squared exponential kernel k and the hyperparameters B, λ, R, δ (in the notation of [3]) ⊤ were set to (E[rx (b)])x∈X (k(x, x′ )) √x,x′ ∈X (E[rx (b)])x∈X , max( max rx (b) − Eb∼B [rx (b)]), λ, 0.1, respectively. x∈X b∈[0,1]N

The Taylor centers of SSEM-UCB-JO were placed at the means of the instantaneous rewards. The RKHSs were comprised of the elemental component functions allowed during the DAG generation. Their norm was surrogated with the norm-equivalent ℓ1 -norm. For runtime considerations, we implemented the differentiation numerically, infrequently retrained the adjacency graph only if it inadequately described new sampling data (similar to the change point detection in [23] and [24]), and set minuscule thresholds for the predicted component function coefficients and derivative values. For each super arm, an optimistic reward projection can be calculated more efficiently by considering each addend of the Taylor approximation separately, while the theoretical regret bound remains applicable. This way, each confidence score can be set to its absolutely lower or higher bound depending on the sign of the accompanying partial derivative value and whether any of the confidence intervals surround 0. We also record that scaling down Ct or decreasing the external UCB’s weight over time often led to faster convergence within the given time horizon. In a second trial instance, we alter the setup described in Section V-A by replacing the scaled periodic x 7→ sin(2πx) components with quadratic and hyperbolic tangent functions. We also include SSEM-UCB-Norm (Algorithm 3 in Appendix B2) in our analysis. Unlike SSEM-UCB-JO, it is better suited to automatic differentiation in its implementation but only entails two confidence scores per super arm to identify by optimization. A third session additionally sees the instantaneous rewards drawn from arcsine distributions with means in range [1, 1.5] and standard deviations in [0.2, 0.6]. The number of base arms is raised to 15 while the edge density is decreased to 0.15 and function coefficients are generated from [0.05, 0.15]. Comparison of the average results over 3 runs displayed in Figure 4 and Figure 5 with those shown in Section V-A indicates that SSEM-UCB-JO favors decentralized signal distributions with increased variances that intensify the nonlinearity

19

TABLE II S TATION ACRONYMS

Fig. 4. Policies’ Mean Regret on Synthetic Normal Data.

Abbreviation EISN GOTH LEIN NEUD NORD ERFU BEBR GÖTT WEIM BADH NORT GOSL KASH JENW JENP SAAL KASW JENG HALB FULD

Train Stop Eisenach Gotha Leinefelde Neudietendorf Nordhausen Erfurt Hbf Bebra Göttingen Weimar Bad Hersfeld Northeim (Han) Goslar Kassel Hbf Jena West Jena Paradies Saalfeld (Saale) Kassel-Wilhelmshöhe Jena-Göschwitz Halberstadt Fulda

Fig. 5. Policies’ Mean Regret on Synthetic Arcsine Data. Fig. 6. Super Arm Selection over Time.

of the system, while SSEM-UCB-Norm specifically addresses normally distributed signals. A dominant external exploration of valid options accommodates structural uncertainty but can incur undue initial regrets. The balance can be tuned towards earlier internal exploration according to the causal simplicity of the problem. Another trade-off lies in computation time, as low degrees of SSEM and SGB are much less demanding than higher ordered SSEM and METS. The expense is amplified by greater connectivity in the underlying graph, albeit pressing the qualitative advantage over algorithms that falsely rely on reward separability or linearity. SSEM displays robustness outside the requirements set for the theoretical study, as the hyperbolic tangent kernel is not positive-definite. 2) Real World: For broader exploration, we chose p the alterP native external UCB Et (x) = x[i]4Ct [i] ln(t)/s (see i∈[N ]

Section III-C). As in the synthetic case, we directly computed the partial derivative values around the means numerically and tuned the regularization parameter adaptively. To further increase the computational speed of the SSEM-UCB-JO implementation, we approximated the inverse (I − F̂t )−1 at each

time step t ≤ 342 with a residual neural network NΘ fit to the loss obtained from the L2 -norm of the structural equations ∥NΘ (z) − F̂t (NΘ (z)) − z∥2 on random samples z ∈ RN . Prior normalization of the signals can stabilize the training process against outliers. Scaling the average reward signal also regulates the frequency of the exploration-exploitation tradeoff. Optimal coefficients on the circular graph can be found efficiently through, e.g., ADMM. The super arm exploration of SSEM-UCB-JO during the first 150 time steps is visualized in Figure 6. Although the perpetual graph learning deployed by SSEMUCB-JO does make it agnostic to daily changes in the unshuffled data, in the future, incorporating dedicated dynamic methods [18], [23], [24] would be of particular interest to trace relations on different time frames. Larger-scale experiments could draw from further sources of periodically disclosed German rail network statistics.

Record · ID 271827 · SHA-256 200563135e4c20d4
Retrieved via Conceptio — every document is proof-bundled with source, license, and retrieval metadata.