Tackling Stackelberg Network Interdiction against a Boundedly Rational Adversary

1 Introduction

Network interdiction is a well-studied topic in Artificial Intelligence. There are many practical problems [Smith and Song, 2020], such as in cyber systems, that can be modeled as a network interdiction problem. In literature, many variations in models of network interdiction exist, and consequentially, a variety of techniques have been used for solving different types of these problems. Our work focuses on a particular type in which there is a set of critical nodes to protect within a larger network. We employ a popular network interdiction model [Fulkerson and Harding, 1977; Israeli and Wood, 2002], where the interdictor (defender)

uses a randomized allocation of limited defense resources for the critical nodes. The adversary traverses the graphs starting from an origin and reaching a destination. There is an interaction with the defender only if the adversary crosses any critical node. The interaction is modeled using a leader-follower (Stackelberg) game where the defender first randomly allocates resources and then the adversary chooses its path accordingly.

Motivated by the fact that human adversaries in real-world security domains often act non-optimally [Tambe, 2011], we model the adversary behavior in our game setting using a dynamic Quantal Response model (an instance of well-known dynamic discrete choice (DDC) models [Rust, 1987; Aguirregabiria and Mira, 2010]). While many real world security applications have benefited from bounded rational Quantal Response model in single shot game settings [Tambe, 2011; Yang et al., 2012; Fang et al., 2016; Bose et al., 2022], to the best of our knowledge, existing works in sequential network interdiction unrealistically assume perfectly rational adversaries and make use of the linearity to utilize linear programming techniques to tackle the problem [Smith et al., 2009; Smith and Song, 2020]. We are the first to explore the DDC model of bounded rational adversaries in the network interdiction setting and formulate the defender's problem as a nonlinear optimization, leading to the requirement of solving the network interdiction problems via nonlinear optimization techniques.

While there is a closed form of the DDC adversary choice probabilities in our game setting, which is mathematically interesting in itself, the closed form presents computational challenges as the naive computation of any such probability involves reasoning about exponentially many paths from origin to destination and is a non-convex problem. This presents challenges beyond those observed in the single shot setting with quantal responding adversary [Fang et al., 2016; Mai and Sinha, 2022]. Thus, we address the challenge of solving such complex non-convex optimization problem for the defender with two different novel approximation algorithms.

First, we introduce an MILP-based method, named LiSD (Linearization via Sampling and Discretization). The solution of LiSD is a bounded approximation for the interdiction problem. LiSD is the result of an innovative combination of path-sampling with piece-wise linear approximation (PL) techniques. Path sampling tackles the computational challenge of exponentially many paths while PL provides a nearoptimal defender strategy solution with a guaranteed bound.

Second, we propose an efficient Dynamic Programming method, named DynP. Essentially, DynP provides a compact and tractable formulations of the defender utility function and the optimization objective's gradient even though these terms involve summing over exponentially many paths. This is accomplished by exploiting recursive relationships among adversary utility-related terms across different paths that involves in the defender's optimal strategy computation. By employing dynamic programming, we can follow a gradient descent approach that is computationally efficient at each step to optimize the defender strategy. Furthermore, while DynP is computationally efficient, it does not guarantee global optimality due to the non-convexity of the defender problem. We thus identify a special case in which the adversary can visit only one (any one) critical node and show that the optimization is unimodal in that case, implying that this problem can be solved optimally in a tractable manner using gradient descent. We further identify specific conditions under which the solution to the restricted problem provides approximation guarantees for the original unrestricted one.

Notation: Boldface characters represent matrices or vectors or sets, and ai denotes the i-th element of a if a is indexable. We use [m], to denote the set ¶1↪ ▷ ▷ ▷ ↪m♦.

2 RelatedWork

Dynamic discrete choice (DDC) models. From the seminal work of [Rust, 1987], DDC models have been widely studied and used to analyze sequential looking-forward choice behaviors and have various applications, e.g., on fertility and child mortality [Wolpin, 1984], on job matching and occupational choice [Miller, 1984], on bus engine replacement [Rust, 1987], and on route choice analysis [Fosgerau et al., 2013; Mai et al., 2015]. Among existing DDC models, the logit-based DDC has been popular due to its closed-form formulation [Rust, 1987]. This model can be viewed as a dynamic version of the well-known multinomial logit (or Quantal Response) model [McFadden, 1981; Train, 2003]. In transportation modeling, logit-based DDC was utilized to develop models to predict people's boundedly rational path-choice behavior [Fosgerau et al., 2013; Mai et al., 2015]. As highlighted in [Zimmermann and Frejinger, 2020], such a model presents synergies with the stochastic shortest path problem [Bertsekas and Tsitsiklis, 1991].

Network interdiction. Our work is a boundedly rational version of the well-studied shortest path interdiction problem [Fulkerson and Harding, 1977; Israeli and Wood, 2002]. Existing work only consider perfectly rational adversaries [Smith et al., 2009; Smith and Song, 2020]. The shortest path and other network interdiction problems with perfectly rational adversaries are generally NP-hard and have strong connections with the areas of bi-level optimization [Dempe et al., 2015] and robust optimization [Ben-Tal and Nemirovski, 2002]. We refer the readers to [Smith and Song, 2020] for a comprehensive review. Our work explores the DDC framework to model bounded rational adversaries, resulting in a significantly more challenging defender problem as it involves complex nonlinear optimization. Besides, there are other variant models where the problem data is not perfectly known to players [Cormican et al., 1998], or where the players repeatedly make their actions [Sefair and Smith, 2016], or where online learning is involved [Borrero et al., 2016].

Network security games and others. Our work also relates to static Stackelberg security game models with Quantal Response adversaries [Yang et al., 2011; Yang et al., 2012; Haghtalab et al., 2016; Mai and Sinha, 2022; ˇCern`y et al., 2021; Milec et al., 2020; Bose et al., 2023b]. In dynamic models named as network security games [Jain et al., 2011], the set-up is different from our work as in this work the rational adversary aims to reach a target and stop, whereas in our work the boundedly rational adversary can attack multiple targets. Other related works along this line only consider zero-sum network security game setting [Xue et al., 2021; Xue et al., 2022]. A Quantal Response type relaxation for network security game was also studied, where the focus in on smart predict and optimize [Wang et al., 2020], however, the optimize part is done using standard non-linear solver such sequential quadratic program with no guarantees.

There are other related game models where players act in a graph-based environment, including pursuit-evasion and security patrol games [Zhang et al., 2019; Basilico et al., 2009; Basilico et al., 2017]. However, these works do not consider the attacker's bounded rationality. Additionally, their strategy spaces and problem settings are characterized differently which involve aspects of real-time information or alarm signals., etc.

3 ProblemFormulation

3.1 Stackelberg Network Interdiction Games

Our network interdiction problem is a leader-follower game with a single adversary. The game is played on a network (graph) (S↪ A) where S is a set of nodes S = ¶1↪ 2↪ ▷ ▷ ▷ ↪ ♣S♣♦, and A is a set of arcs. We formulate the problem as a twoplayer network interdiction game . The follower (adversary) takes a path through this network, which is sampled from a distribution as described below. The origin so 2 S is a given starting node. In our problem, we also assume the existence of a sink (or destination) node sd 2 S that the adversary ultimately reaches. Let L be the set of critical nodes (i.e., subset of nodes in the network) that the defender can interfere or alter. From the leader's (defender's) viewpoint, the aim is to assign M resources to nodes s 2 L; each such assignment is a defender pure strategy. Further, nodes and resources are of certain types such that nodes of a given type can only be protected by resources of that same kind. Let there be K types of nodes. Let the number of resources of each type k be Mk, hence Pk2[K] Mk = M. Also, let ¶Lk♦k2[K] be a partition of the set of nodes L by the types of the nodes.

A mixed strategy is a randomized allocation resulting in a coverage vector x = ¶xs↪ s 2 L↪ Ps2Lk xs  Mk↪ 8k 2 [K]♦ where xs is the marginal probability of covering node s, which then impacts the adversary's path choice probabilities. Given a node s 2 S, if the adversary crosses this node, then the defender gets a node-specific reward rl (s↪ xs). The Stackelberg equilibrium can be computed by solving the following problem [Yang et al., 2012; Mai and Sinha, 2022]:

where Rl (⌧♣x) = Ps2L\⌧ rl (s↪ xs) is the defender's accumulated reward on path ⌧ and P f(⌧♣x) is the probability the attacker follows the path ⌧ (of which computation is discussed in the behavior modeling part). Here, [Lx ↪U x ] represent the required lower bound and upper bound on the coverage probability for each node in the critical set L.

3.2 Boundedly Rational Adversary Behavior

We model the adversary's bounded rational behavior using the dynamic discrete choice framework (and specifically the logit-based recursive path choice model [Fosgerau et al., 2013]). A known property in this setting is that the bounded rational adversary chooses a policy that is equivalent to a static multinomial logit (MNL) discrete choice model over all possible paths [Fosgerau et al., 2013].

Concretely, let U(⌧ ♣s0↪ x)= Ps2⌧ v(s; x) be the deterministic long-term utility of the adversary when starting in s0; if s0=so, then we simply write U(⌧ ♣x). Here, v(s; x) is the adversary's utility associated with node s when the defender's strategy is x. Given x, the probability the adversary follows a path ⌧ can be computed as follows [Fosgerau et al., 2013]:

given ⌦ is the set of all possible paths and θ is the parameter which governs the follower's rationality. Thus, we can view the logit-based dynamic discrete choice formulation as a soft version of the shortest weighted path problem from the source so to destination sd. Given the adversary behavior model, the adversary's expected utility can be computed as an expectation over all paths, as follows:

Our Prop. 1 shows that the adversary's expected utility approaches the best accumulated utility (smallest path weight) as θ tends to zero (we drop the fixed strategy x for simplicity).

As a result, limθ!0 E f = U(⌧ ⇤ ). 1

1 All proofs, if not presented, are included in the appendix.

4 Common Binary Search Framework

Overall, (OPT) is computationally challenging since the objective not only involves an exponential number of paths in the network but also is non-convex. To address this computational challenge, we propose two new different algorithms which share the common underlying binary search framework. The purpose is to reduce the original fractional (OPT) to a simpler non-fractional problem. These algorithms then differ in applying different efficient techniques to solve each binary search step. We elaborate them in subsequent sections.

Essentially, we write the objective of (OPT) as follows:

F l (x) has a fractional non-convex form. A typical way to simplify this structure is to use the Dinkelbach transform and a binary search algorithm [Dinkelbach, 1967] to convert the original problem into a sequence of simpler ones. We use binary search to write (OPT) equivalently as: maxλ nλ ∣∣∣ 9x s.t. F l (x) ≥ λo which is equivalent to finding a maximum value of λ 2 R such that the following subproblem:

has a non-negative optimal objective value. Overall, (3) is still non-convex, but no longer fractional. In addition, the set ⌦ of all feasible paths can be huge and may not be enumerable. Therefore, we propose two different algorithms (as elaborated next) to tackle these challenges in solving (3).

5 Linearization via Sampling and Discretizing

We describe our first near-optimal method, LiSD, which involves exploring path-sampling with piece-wise linear approximation (PL) techniques to approximate (3) by a MILP. Path sampling tackles the computational challenge of exponentially many paths while PL provides a near-optimal defender strategy solution with a guaranteed bound for (3).

5.1 Sample Average Approximation

We first approximate the sum over ⌦ via sample average approximation. That is, we select a feasible solution x0 to create a fixed distribution over paths in ⌦. By dividing the objective of (3) by P⌧2⌦ exp(U(⌧;x0)◁θ), which is a constant, we aim to maximize the following objective function:

where e U(⌧ ♣x) = U(⌧ ;x) θ - U(⌧;x0) θ ↪ and D(x0) is the distribution over paths ⌧ with probabilities P f (⌧ ♣x0) (Eq. 2).

We now can approximate the objective function g(x↪ λ) by sample average approximation. Specifically, let ⌧1↪ ▷ ▷ ▷ ↪ ⌧N be N samples from D(x0), we approximate g(x↪ λ) by:

Essentially, the approximation bGN (x↪ λ) converges to G(x↪ λ) almost surely as N ! 1 and the approximation errors can be bounded as shown in Proposition 2.

Proposition 2. For any given ⇠ > 0, we have:

Proposition 2 implies that bGN (x↪ λ) will converge to the true function G(x↪ λ) in probability with an exponential rate as the number of samples N increases. This is a direct result fromHoeffding's inequality [Hoeffding, 1994].

5.2 Piece-wise Linear (PL) Approximation

We now further approximate bGN (x↪ λ) by a PL function, allowing the subproblem to be solved to near-optimality via a MILP solver. First, for each ⌧n, we introduce new variables un =Rl (⌧n♣x) and vn = e U(⌧n♣x) θ . We now can re-write the

Let Ln and Un be an lower and upper bounds of vn. The PL approximation can be done by partitioning each interval [Ln↪Un] into K sub-intervals of equal size, and introducing K binary variables z1n ↪ ▷ ▷ ▷ ↪ znK such that z1n ≥ z2n ≥ ▷▷▷ ≥ zK n , to represent each interval. Intuitively, zkn = 1 implies the kth sub-interval involves in the approximation of exp(vn) and zkn = 0, otherwise. Let ∆n = (Un-Ln)◁K (i.e., the size of each interval) and δkn , k 2 [K] is the slop of function evn in the interval [Ln +∆n(k - 1)↪ Ln +∆nk]:

Each component exp(vn) can be approximated as follows:

We then can re-write the sub-problem (3) as follows:

which maximizes the piece-wise approximation of bGN (x↪ λ). The additional variable n captures the gap between vn and the binary approximation Ln + ∆n Pk2[K] zkn .

Finally, there are only some bi-linear terms left to be linearized in the objective function. We do that using McCormick inequalities. Specifically, let Lun and U u n be lower and upper bounds of un, we introduce new variables skn to present (un-λ)zkn , we can now linearize the bi-linear term (un-λ)zkn with the following additional constraints:

The above three constraints guarantee that when zkn = 1, then skn = un - λ. Conversely, when zkn = 0, then skn = 0.

By combining the above new variable skn and constraints with (MINLP), we obtain the MILP reformulation:

s.t. Constraints (6-12) are satisfied.

We further establish a performance bound for PL approxima-

Proposition 3. Assume that λ 2 ⇥ minn¶Lun ♦↪maxn¶U u n ♦⇤, let bxNK be an optimal solution to (MILP) and x⇤ be optimal for average approximation sub-problemmaxx bGN (x↪ λ), then we obtain the following inequality:

So, it is sufficient to consider λ 2 [minn¶Lun ♦↪maxn¶U u n ♦]. This allows us to state Proposition 3 below.

Theorem 1. Assume that λ 2 ⇥ minn¶Lun ♦↪maxn¶U u n ♦⇤. Let bxNK be an optimal solution to (MILP) and x⇤ be optimal for maxx G(x↪ λ), then given any ⇠ > 0, if we choose N,K such that NK  ⇠ 6B , then we have:

From an intuitive standpoint, augmenting K would diminish the approximation error of the PL approximation. Conversely, augmenting N has a dual effect: while it lessens the error arising from path sampling, it simultaneously heightens the cumulative error stemming from all the samples. In fact, to drive the bound closer to zero, Proposition 3 indicates that it's necessary that the rate of increase for K should surpass that of N. We further investigate this dual effect by looking at the quality of a solution returned from (MILP) w.r.t the original sub-problem maxx G(x↪ λ). A performance bound is provided in Theorem 1, which implies that, under the condition N  ⇠K 6B , the solution given by the PL approximation will converge in probability to a true optimal solution, with an exponential rate.

This result can be employed to establish (theoretical) estimates for N and K to achieve a desired performance.

The above estimates might shed light on how N↪K depends on the performance criteria ↵↪ β. We note that these estimates would be conservative, as in practice we may need much smaller N↪K to achieve the desired performance. A final note is that one can employ an off-the-shelf solver (e.g. CPLEX or GUROBI) to solve (MILP). Although this program would be large in size, SOTA solvers can efficiently handle very large MILPs, aided by powerful machines.

6 Dynamic Programming Based Solution

The above MILP approximation involves binary variables and would be intractable in large scenarios. We thus propose an alternative new algorithm, DynP that also follows binary search, but at each binary step, (i) it presents a non-trivial compact representation of the objective function based on the creation of a dynamic program, which handles an exponential number of paths; and (ii) it applies a gradient ascent-based method to efficiently solve the resulting compact problem.

6.1 Compact Representation

We can rearrange terms in the objective of sub-problem (3) according to critical nodes as follows:

Since g(x↪ λ) is differentiable, this maximization problem can be solved for a local maximum by a gradient-based method. One of the key challenges is the computation of g(x↪ λ), which, if done naively, would require enumerating exponentially many paths on ⌦. We next show that g(x↪ λ) has a compact form, which allows us to compute g(x↪ λ) and its gradient efficiently via dynamic programming.

For a compact representation of g(x↪ λ), we introduce the following new terms for all nodes s↪ s0 2 S:

where ⌦sd (s) is the set of all paths from s to the destination sd and ⌦(s0 ↪ s) is the set of all paths from s0 to s.

The objective g(x↪ λ) can be now re-formulated as follows:

where so is the origin. Although these new terms still involve exponentially many paths in ⌦sd (s) and ⌦(s0 ↪ s), they can be computed efficiently via dynamic programming.

Indeed, ¶Zs♦s can be computed recursively as follows:

Algorithm 1: Dynamic Programming based algorithm (DynP) to solve Maximizing g(x↪ λ)

Input: λ 2 R and an initial value of x

while not converged do

where N(s) = ¶s0 2 S♣ (s↪ s0 ) 2 A♦, is the set of possible next nodes that can be reached in one hop from node s 2 S. Let M be a matrix of size ♣S⇥S♣ with entries defined as:

Then Z = ¶Zs↪ s2S♦ is a solution to the linear system Z = MZ + b, where b is of size ♣S♣ ⇥ 1 with zero entries except bsd =1. Similarly, we can compute Ys =¶Y s s0♦s0 recursively:

Clearly, Ys is a solution to the linear system Ys = MYs +bs , where bs is of size ♣S♣ with zeros everywhere except bss = 1. Since Ys and Z are solutions to the systems Ys = MYs + bs and Z = MZ + b, respectively, 8s 2 S, the objective g(x↪ λ) can be computed via solving ♣L♣ + 1 system of linear equations. Finally, we see that all the above linear systems rely on the common matrix M. We can group them all into only one linear system. Let H be a matrix of size (♣S♣) ⇥ (♣L♣ + 1) in which the 1st to ♣L♣-th columns are vectors Ys , s 2 L and the last column is Z. Let B be a matrix of size (♣S♣)⇥(♣L♣+1) in which the 1st to ♣L♣-th columns are vectors bs , s 2 L and the last column is b. We see that H is a solution to the linear system (I -M)H = B. Thus, in general, we can solve only one linear system to obtain all Ys and Z. This way should be scalable when the size of L increases.

6.2 Gradient Computation

We aim at employing the gradient-based approach to solve the binary search step: maxx ¶g(x↪ λ)♦ (aka. Eq. 3). The core is to compute the gradient {@g(x↪λ)◁@xs . According to Eq. 14, this gradient computation requires differentiating through the matrices Z and ¶Ys ♦ (or equivalently, differentiating through the matrix H). We first present our Proposition 4:

Proposition 4. (I -M) is invertible in a cycle-free network.

Prop. 4 allows us to compute the matrix H as: H = (I - M)-1B▷ By taking the derivatives of both sides w.r.t xj, j 2 L, we obtain the following: for all j 2 L,

where J H↪j and J M↪j are the gradient matrices of H and M w.r.t xj , i.e., J H↪j is a matrix of size ♣S♣⇥(♣L♣+1) with entries J H↪j ss0 = @Hss0◁@xj , and J M↪j is a matrix of size (♣S♣ ⇥ ♣S♣)

with entries J M↪j ss0 = @Mss0◁@xj , for any j 2 L. Let Rl (x) be a matrix of size 1 ⇥ ♣L♣ with entries rl (s↪ xs) for s 2 L. We use AS↪T to denotes a sub-matrix of A which uses the rows in set S and columns in set T. If S or T is a singleton, e.g., S = ¶so♦ or T = ♣L♣ + 1, then we write it as so or ♣L♣ + 1.

As a result, we now can compute the required gradient as follows for all s 2 L where ◦ denotes Hadamard product:

We summarize the main steps to optimize g(x↪ λ) in Alg. 1.

Remark 1. Alg. 1 only guarantees a local optimum due to the non-convexity of g(x↪ λ). The complexity is determined by the matrix inversion which, in worst case, is in O(♣S♣3 ). The gradient descent loop runs O(1◁✏) to provide an additive ✏ approximation. Thus, the total complexity is O((1◁✏)♣S♣♣A♣). In practice, the gradients can be found via auto differentiation techniques, providing significantly more speed-up.

6.3 ANatural Special Case

Separation of critical resources and/or privileges is an important concept in security [Lin et al., 2023]. Following this principle, we analyze a special yet natural security design scenario where that the critical nodes L are well separated. Specifically, we assume that the cost of travelling between nodes in L is high. More formally, given a critical node s 2 L, let ∆+(s) be the set of paths that cross s and at least another critical node in L. Let β1, β2 > 0 such that:

Intuitively, β1 and β2 are expected to be small if the cost of traveling between any two critical nodes in L is large. That is, β1↪ β2 ! 0 as P⌧2⌦(s↪s0 ) exp (U(⌧;x)◁θ ) ! 0, where ⌦(s↪ s0 ) consists of all paths from s to s0 , for any s↪ s0 2 L. Surprisingly, even though DynP only finds a locally optimal solution for (OPT) due to its non-convexity, we show that assuming small β1 and β2 provides approximation guarantees for the globally optimal solution value.

For this approximation, we need mild assumptions that the utilities have a linear form: v(s; x) = wf s xs + tfs and rl (s; x) = rl (s↪ xs) = wlsxs + tls for some constants wf s ↪ tfs ↪ wls ↪ tls . We assume that wf s <0 and wls >0, i.e., more resources xs at s will lower adversary's utilities, and increase the defender's utility. This setting is intuitive for security settings [Yang et al., 2012; Mai and Sinha, 2022].

We first introduce a restricted interdiction problem that can be solved optimally in a tractable manner using our efficient gradient descent-based method. We then present an important theoretical result showing how the restricted problem's solution yields an approximate solution of with the original problem for well separated critical nodes.

Let ∆(s) be the set of paths that cross a critical node s and do not cross any other node in L. We consider the following restricted interdiction problem:

Intuitively, in this restricted problem (Approx-OPT), the adversary's path choices are restricted to a subspace of paths in the network which only cross a single critical node in L. We denote by X, the feasible set of the defender's interdiction strategies x that satisfy the constraints in (Approx-OPT).

Solution Relation with Original Problem (OPT)

We now theoretically analyze (Approx-OPT)'s solution in relation to our original problem (OPT). We prove that:

Theorem2. Let x⇤ be an approx. solution to (Approx-OPT): maxx2X e F(x) such that e F(x⇤ ) ≥ (1 - ✏)maxx e F(x) for given ✏ > 0, let  = maxx2X Ps2L ♣rl (s↪ xs)♣ be the maximal absolute reward that the defender can possibly achieve at a critical node, then we obtain:

Additionally, if x⇤ is an approx. solution with an additive error ✏ > 0, we obtain the following bound:

+minx2X e Note that maxx2X¶F l (x)♦ is the original problem (OPT) to find an optimal defender strategy. As stated previously, when the cost of traveling between any two critical nodes is high, (β1↪ β2) is close to zero, meaning the RHS of both inequalities in Theorem 2 will closely approximate the optimal solution value of (OPT).

Solving the Restricted Problem: To solve the restricted problem, we also apply binary search. It can be demonstrated that each sub-problem of the binary search can be effectively solved to optimality (or near-optimality) using a gradientbased method. Due to limited space, the details of this approach is provided in the appendix along with the sub-results that lead to the main Theorem 2.

7 Numerical Experiments

To illustrate the efficacy of our proposed algorithms, we perform experiments on synthetic data.

7.1 Experimental Settings

Data generation. We generate random graphs (cycle-free) with ♣S♣ vertices and edge probability p. We randomly choose ♣L♣ vertices (except source and destination) as the critical nodes that can be attacked. We set ♣L♣ = 0▷8 ⇥ ♣S♣. In addition, the defender weights {(wlj ↪ tlj ) ∣∣ j 2 ⇥♣L♣⇤ are generated uniformly at random from the interval [0↪ 1] and the adversary weights {(wf j ↪ tfj ) ∣∣ j 2 ⇥♣S♣⇤ are generated at random from the interval [-1↪ 0]. Moreover, we used p = 0▷8↪ θ = 2.

Baseline. We approximate the sums over exponentially many paths in Equation 3 by sampling paths from the network and run gradient descent on this expression to estimate the optimal decision variable. To sample paths for the baseline, a resource allocation x is assigned to L and the follower is initially placed at the origin so. Its next node s1 is sampled from the distribution

where s 2 ¶nodes having an edge to s0♦ and N(s0) is the set of outgoing nodes from s0. Similarly s2↪ ▷ ▷ ▷ etc. are sampled till the destination sd is reached. This sampling is repeated 1000 times per iteration and then average is taken to get the objective. Based on the gradients, the resource allocation x is updated which changes the transition probabilities and the process is repeated again until convergence. Ten different values of x were taken and the seed with the lowest loss was reported. We will compare this baseline against our near optimal MILP-based algorithm, LiSD and our dynamic programming based algorithm, DynP. To ensure fairness, all algorithms were run with the same number of epochs.

Choosing N andK for the LiSD To justify our choices of N and K for the LiSP described, we first fix N = 90 and vary K from 5 to 100. For each value of K, we generate 10 independent instances and solve them using the MILP approach. We then observe that the optimal values given by K = 20 are only 0.3% different from those given by the largest value of K (i.e., K = 100). Moreover, the optimal values given by N = 90 are only about 3% different from those from the largest value of N, i.e., N = 150. The numerical details can be found in the appendix. We therefore choose N =90 and K=20 for the LiSD approach. According to the above analyses, these choices would suffice to guarantee low practical approximation errors stemming from both path-sampling and PL approximation. We use GUROBI (a SOTA MILP solver) to solve (MILP). All our experiments were run on a 2.1 GHz CPU with 128GB RAM.

7.2 Numerical Comparison

We vary the number of nodes from 20 to 100. For each choice of number of nodes, we generate 20 independent instances and solve them by the three methods (Baseline, DynP, and LiSD). The rates of giving best objective values (over 20 instances) are reported in Table 1. LiSD consistently outperforms the other methods, by a large margin, in terms of the number of times it returns the best objective values. DynP performs better than Baseline for large-size settings, but worse than Baseline for small-size ones. Note that, among the four methods, only LiSD can guarantee near-optimal solutions. The computing times are not directly comparable, as GUROBI ran on several cores while the other gradient-based methods use only one CPU core. We however observe that, for instances of 100 nodes, the average computing times for the Baseline, DynP, and LiSD are approximately 3, 15 and 1▷8 minutes.

We further compare the objective values returned by the four methods by computing the percentage improvement of each objective w.r.t. the lowest objective value given by the four approaches. Specifically, we solve each instance to obtain 4 objective values. We then compute the percentage improvement of each objective value w.r.t the lowest value among the four values. The average percentage values are reported in Table 2 below, which show that LiSD performs the best, and DynP outperforms the Baseline.

8 Conclusion

Network interdiction game problems present a set of challenges that appear intractable to start with. In this work, we address some of these challenges by providing novel methods that efficiently solve a class of network interdiction problems with approximation guarantees. We are the first to study the dynamic Quantal Response model in the type of network interdiction studied in [Fulkerson and Harding, 1977; Israeli and Wood, 2002]. We believe this modeling and methodology contribution provides suggestions for other future research directions, such as a variant where the adversary's objective is to maximize a flow through the network or a setting where the leader would also need to make dynamic and time-dependent decisions. This online nature of the problem suggests possible future work in online learning problems such as [Bose et al., 2023a]. It is interesting to analyse scenarios where network structures arise naturally such as ride pool matching [Bose and Varakantham, 2021] and are under threat from adversaries such as competing providers.

Acknowledgments

Dr. Tien Mai and Dr. Arunesh Sinha were supported by the National Research Foundation Singapore and DSO National Laboratories under the AI Singapore Programme

• (AISG Award No: AISG2- RP-2020-017). Dr. Nguyen was supported by grant W911NF-20-1-0344 from the US Army Research Office.