Christopher Sebastian

Hidalgo Calva

Alejandro P Riascos

Optimal exploration of random walks with local bias on networks

We propose local-biased random walks on general networks where a Markovian walker is defined by different types of biases in each node to establish transitions to its neighbors depending on their degrees. For this ergodic dynamics, we explore the capacity of the random walker to visit all the nodes characterized by a global mean first passage time. This quantity is calculated using eigenvalues and eigenvectors of the transition matrix that defines the dynamics. In the first part, we illustrate how our framework leads to optimal exploration for small-size graphs through the analysis of all the possible bias configurations. In the second part, we study the most favorable configurations in each node by using simulated annealing. This heuristic algorithm allows obtaining approximate solutions of the optimal bias in different types of networks. The results show how the local bias can optimize the exploration of the network in comparison with the unbiased random walk. The methods implemented in this research are general and open the doors to a broad spectrum of tools applicable to different random walk strategies and dynamical processes on networks.

I. INTRODUCTION

Markovian random walks on networks are found in different contexts and have been studied intensively over the last decades, from a purely theoretical perspective [START_REF] Hughes | Random Walks and Random Environments[END_REF][START_REF] Masuda | [END_REF][START_REF] Michelitsch | Fractional Dynamics on Networks and Lattices[END_REF][START_REF] Riascos | [END_REF] to diverse applications such as database exploration [START_REF] Leskovec | Mining of Massive Datasets[END_REF][START_REF] Blanchard | Random walks and diffusions on graphs and databases. An introduction[END_REF], web page ranking [7,8], epidemic spreading [9], encounter networks [10], human mobility [11], among many others. These dynamics have been implemented in diverse fields as processes that are able to efficiently reach hidden targets or to simply explore a particular region of space [START_REF] Riascos | [END_REF].

In particular, degree-biased random walks describe a dynamics where a walker in a particular network node uses the degree information of neighbors to choose randomly a new node to visit [12][13][14]; in this case, a tunable parameter β controls the bias in the complete network. This process is of interest in theory and applications to understand how β affects the exploration of a network in different topologies [START_REF] Riascos | [END_REF][12][13][14]. The spectrum of applications includes traffic dynamics and routing protocols [13], the study of simultaneous random walkers [15][16][17][18], extreme events [19], epidemic spreading [20], just to mention a few examples. Recently, degree-biased random walks have been generalized to include multiple biases [21], in potential-driven random walks [22], to study the influence of damage and aging in complex systems [23,24].

In different studies, degree-biased random walks are defined with a global parameter β; in this manner, its value is the same for all the nodes in the network. However, this parameter could be local, i.e., with values in each node. This type of dynamics is illustrated in Fig. 1. In this example, the transition probability from a node i to one of its neighbors depends on a parameter β i . Hence, the transition probability matrix that defines the random walker is determined by a vector β = (β 1 , β 2 , . . . , β N ), being N the size of the network. In connected undirected networks, finite values of β i define a dynamics capable of reaching any node from any initial condition.

In this paper, we explore degree-biased random walk dynamics with local bias defined by a vector β. In Sec. II, we present a summary of the formalism to obtain mean-first passage times and its weighted average over all the nodes to describe the transport of ergodic random walkers on networks. In Sec. III, we introduce the formal definitions of random walks with a local bias and the motivation for this stochastic process. We analyze general properties and the optimal dynamics in graphs with small sizes (N = 4, 5) through the exploration of all the possible bias configurations. In Sec. IV, we analyze networks with different sizes and topologies using simulated annealing. This heuristic method allows us to efficiently identify particular values of the local bias that are most favorable (or close to an optimal) to reach all FIG. 1. A random walker with local bias illustrated as different transition probabilities between the nodes in a network. The possible movements from a node i to its neighbors are given by transition probabilities defined by β = (β1, β2, . . . , βN ), where N is the size of the network and βs a parameter associated to each node s. The walker can modify each element of β trying to discover a strategy that optimizes its capacity to reach all the nodes of the network.

the nodes of the network. We explore this algorithm in different types of networks. The framework introduced in this research is general and can be extended to other random walk strategies with local parameters to define the dynamics. This approach may have applications in routing processes, to address common problems in traffic planning in public transportation, among other cases where it is necessary to tune parameters in the parts of a complex system.

II. RANDOM WALKS ON NETWORKS

A. Random walks and mean first passage times

In this section, we discuss the mean first passage times (MFPTs) and a global time to characterize random walks on networks. We consider the dynamics on connected undirected networks with N nodes i = 1, . . . , N described by an adjacency matrix A with elements A ij = A ji = 1 if there is an edge between the nodes i and j and A ij = 0 otherwise; in particular, A ii = 0 to avoid edges connecting a node with itself. The degree of a node i is given by k i = N ℓ=1 A iℓ , this is the number of neighbors of i. The occupation probability to find a Markovian random walker in a node j at time t, starting from i at t = 0, is given by P ij (t) and obeys the master equation [START_REF] Hughes | Random Walks and Random Environments[END_REF][START_REF] Riascos | [END_REF]25]

P ij (t + 1) = N ℓ=1 P iℓ (t)w ℓ→j , (1)
where we consider discrete times t = 0, 1, . . .. In Eq. ( 1), w ℓ→j is the transition probability to randomly hop in one step from node ℓ to j. These probabilities define a N × N matrix W with elements W ij = w i→j . All the information of the dynamics is defined by W. In the following, we study the dynamics with transition matrices W for which the random walker can reach in finite time any node starting from any initial condition; in this manner, W defines an ergodic process.

Using the probability P ij (t) and the stationary distribution

P ∞ j ≡ lim T →∞ 1 T T t ′ =0 P ij (t ′ ), the moment R (0) ij is defined as [25] R (0) ij ≡ ∞ t=0 {P ij (t) -P ∞ j }. (2)
In terms of these quantities, the analytical form of the MFPT T ij that gives the average number of steps to start in i and reach for the first time the node j is [25]

T ij = 1 P ∞ j R (0) jj -R (0) ij + δ ij . (3)
In order to calculate T ij it is necessary to find P ij (t).

In the following, we use Dirac's notation, then

P ij (t) = i| W t |j , (4)
where {|m } N m=1 represents the canonical base of R N . For ergodic random walks on connected undirected networks, the matrix W can be diagonalized. For right eigenvectors of W we have W |φ i = λ i |φ i for i = 1, . . . , N , where the eigenvalues λ i can take complex values and are sorted in the form λ 1 = 1 and 1 ≥ |λ l | ≥ 0 for l = 2, 3, . . . , N . With this information, we define a matrix Z with elements Z ij = i|φ j [START_REF] Riascos | [END_REF]. The matrix Z is invertible, and a new set of vectors φi is obtained by means of (Z -1 ) ij = φi |j , then

δ ij = (Z -1 Z) ij = N l=1
φi |l l|φ j = φi |φ j and I = ZZ -1 = N l=1 |φ l φl , where I is the N × N identity matrix.

By using the spectral form of W, Eq. ( 4) takes the form

P ij (t) = N l=1 λ t l i|φ l φl |j . (5)
From Eq. ( 5), the stationary probability distribution is P ∞ j = i|φ 1 φ1 |j , where the result i|φ 1 = constant makes P ∞ j independent of the initial condition. Now, as a consequence of the definition of

R (0) ij [4] R (0) ij = N l=2 1 1 -λ l i|φ l φl |j . (6)
Therefore, for i = j in Eq. (3), the MFPT T ij is [4]

T ij = N l=2 1 1 -λ l j|φ l φl |j -i|φ l φl |j j|φ 1 φ1 |j , (7)
whereas T ii = ( i|φ 1 φ1 |i ) -1 is the mean first return time to the node i. In addition, we have the time τ j independent of the initial condition [START_REF] Riascos | [END_REF]26]

τ j ≡ N l=2 1 1 -λ l j|φ l φl |j j|φ 1 φ1 |j . ( 8
)
Using this definition, we have

N i=1 P ∞ i T ij = 1 + τ j , (9)
where we use the relation T ij characterizes the capacity of the dynamics to reach the node j from any initial condition i = j (see Ref. [25] for additional details). Hence, this value allows defining a global time T to quantify the exploration of the whole network [24,26]

T ≡ 1 N N j=1 τ j = 1 N N j=1 i =j P ∞ i T ij . (10)
In this relation, we see that T is a global time that gives the weighted average of the number of steps to reach any node of the network. In the following, we use the global time T to quantify the capacity of the random walker to explore the network using the eigenvalues and eigenvectors of the transition matrix [START_REF] Michelitsch | Fractional Dynamics on Networks and Lattices[END_REF][START_REF] Riascos | [END_REF]. Other alternatives for a global characterization of the dynamics may include mixing times of Markov chains [START_REF] Blanchard | Random walks and diffusions on graphs and databases. An introduction[END_REF] or the cover time defined as the time a random walker requires to visit every node in the network at least once [27].

B. Degree-biased random walks

In this section, we discuss degree-biased random walks. For this case, the random walker hops with transition probabilities w i→j depending on the degrees of the neighbors of the node i. Degree-biased random walks are defined by [14]

w i→j = A ij k β j N l=1 A il k β l , (11)
where β is a real parameter. In Eq. ( 11), β > 0 describes the bias to hop to neighbor nodes with a higher degree, whereas for β < 0 this behavior is inverted and, the walker tends to hop to nodes less connected. When β = 0, the normal random walk with

w i→j = A ij /k i is recovered.
The random walk in Eq. ( 11) is also defined in terms of a symmetric matrix of weights Ω with elements Ω ij = A ij (k i k j ) β with transition probabilities [START_REF] Riascos | [END_REF]

w i→j = Ω ij N l=1 Ω il = Ω ij S i .
Here S i = N l=1 Ω il and represents the total weight of the node i (see Ref. [START_REF] Riascos | [END_REF] for a review of different random walks defined using a symmetric matrix of weights). In terms of this formalism, in connected undirected networks, degree biased random walks are ergodic for β finite and the stationary distribution is given by [4]

P ∞ i = S i N l=1 S l = N l=1 (k i k l ) β A il N l,m=1 (k l k m ) β A lm . ( 12
)
Degree biased random walks have been studied extensively in the literature in different contexts as varied as routing processes [13], chemical reactions [28], extreme events [19,29], multiple random walks on networks [15,18], among others [14,30,31]. Additionally, mean We show the average global time T obtained numerically using Eq. ( 10) for degree-biased random walks defined by the transition probabilities in Eq. ( 11) with β = -2.5, -2.25, . . . , 2.25, 2.5 (continuous lines are used as a guide). The values with dots represent the results for the discrete values β = -1, 0, 1.

field approximations have been explored for diverse cases [14,28,32].

In Fig. 2 we show the global time T in Eq. ( 10) for degree-biased random walkers defined by Eq. ( 11) as a function of β. We examine the results for the dynamics on three random structures with N = 100: an Erdős-Rényi network [33] with probability p = ln(N )/N to define the lines, a Watts-Strogatz network with rewiring probability p = 0.05 [34] and a scale-free Barabási-Albert network, generated with the preferential attachment rule where each newly introduced node connects to m previous nodes, we choose m = 1 [35]. The results are obtained numerically for β ∈ [-2.5, 2.5], our findings illustrate how the parameter β modifies the global exploration of the network in different topologies. In particular, in Fig. 2(a) we see that an optimal exploration is obtained for β = 0, whereas in Fig. 2(b) the random walker may benefit using the bias β = -1 that gives preference to the pass to less connected nodes. On the contrary, for the complex network in Fig. 2(c), the best choice is to use a small bias towards highly connected nodes, for example, implementing a random walk with β ∈ (0, 1).

In addition, regarding the values of T to characterize the dynamics on heterogeneous networks, for |β| ≫ 1 the dy-namics with bias requires a huge number of steps to reach any node since the random walker gets trapped with revisits to the highly connected nodes (for β ≫ 1) or in nodes with reduced connectivity (for β ≪ -1) yielding an optimum for finite values β.

III. RANDOM WALKS WITH LOCAL BIAS

A. General definition

In the case presented in Eq. ( 11), the parameter β defines the global bias of the random walker. However, another alternative is to choose this parameter at each node; in this manner, by using a vector β = (β 1 , β 2 , . . . , β N ) with the local bias β i ∈ R at node i, we define

w i→j ( β ) = A ij k βi j N l=1 A il k βi l . (13)
Then, the vector β and the network structure defines the random walk with local bias, a strategy that allows to control at each node the transitions from a node to more or less connected neighbors. Note that, when β = (β, β, . . . , β), we recover the random walk with global bias discussed in Sec. II B.

In particular, for finite values of β i , the random walker with transition matrix W( β) with elements w i→j ( β ) is ergodic in connected undirected networks, i.e. maintains the capacity to reach any node from any initial condition.

In the general case, W( β) is not defined in terms of a symmetric matrix of weights and we cannot deduce explicitly the stationary distribution. However, the results in Eqs. ( 7)-( 10) are still valid for the analysis of the dynamics with local bias by using the respective eigenvalues and eigenvectors of W( β). In particular, the eigenvalues λ i can take complex values, with λ 1 = 1 unique and 0 ≤ |λ l | ≤ 1 for l = 2, 3, . . . , N , in agreement with the Perron-Frobenius theorem (see Refs. [24,36] for two examples of ergodic dynamics with bias and complex eigenvalues).

Once defined the dynamics with local bias, the question arises as to whether it is possible to find configurations of the β i or the vector β to optimize the capacity of the random walker to reach any node of the network. If these optimal strategies exist, are they unique? How do we find the best values? How are the optimal β related to the network topology?

In the following, we explore those questions. For simplicity, we limit the analysis to a case where the values β i ∈ {-1, 0, 1} representing three types of bias. However, the approach implemented is general and can be applied to other sets of parameters. The capacity to explore the network is evaluated with the global time T ( β) defined in Eq. ( 10) and obtained numerically with the eigenvalues and eigenvectors of W( β). Our main goal is to determine if there exists a reduction of T ( β) in comparison with the dynamics without bias recovered when β = (0, 0, . . . , 0).

B. Dynamics on small graphs

In this section, we analyze random walks with local bias in graphs with N = 4 and N = 5 nodes. Due to their relative simplicity, we can explore the whole configuration space of β and find those that minimize the value of the global time T . For instance, the domain of β for a network with N = 4 nodes has 3 4 = 81 different configurations if β i ∈ {-1, 0, 1}, whereas for N = 5, 3 5 = 243 configurations are possible. In these cases, the resulting random walk dynamics and their performance to explore the network can be evaluated efficiently using Eqs. ( 8) and (10). However, as the number of nodes in the graph increases, computing T for each configuration becomes a difficult task. Our findings for small graphs are reported in Figs. 3 and4.

In Fig. 3, we present the analysis for all the nonisomorphic connected graphs with N = 4. In Fig. 3(a), we show the values of the global time T (sorted in decreasing order) for each configuration of β for a particular network presented as inset. The global time T unbiased for the unbiased case is presented with a dashed line, three optimal configurations with T optimal minimizing T are shown with blue dots whereas the unbiased dynamics is presented with a black dot. Here, it is worth noticing that exploring the whole domain of β, it is possible to find configurations more effective to reach all the nodes of the network than the unbiased random walk, these cases have values T < T unbiased , i.e., with values below the dashed line. Another important finding is that optimal strategies for the local bias are not always unique. In Fig. 3(b), we illustrate the multiplicity of configurations with minimal T for all the non-isomorphic connected networks with N = 4 (the 6 networks G 1 , G 2 , . . . , G 6 were obtained from [37]), where we represent the graphs with different types of nodes according to the values β for the optimal configurations with T optimal . For instance, in our example in inset (a) we have the graph G 3 . In this case, three optimal configurations with β = (β 1 , β 2 , β 3 , β 4 ) are obtained with

β 1 ∈ {-1, 0, 1}, β 2 = -1, β 3 = β 4 = 0.
In addition, in Fig. 3(b) we can see that there exists nodes unaltered by the bias. In these nodes (represented as pie charts) β i ∈ {-1, 0, 1}. These particular nodes fulfill the following condition valid for connected undirected networks:

If all the neighbors of the node i have the same degree k, then

w i→j ( β ) = A ij k βi N l=1 A il k βi = A ij k i . (14)
FIG. 3.

Random walks with local bias on graphs with N = 4. (a) Global times T for the 81 configurations of β with βi ∈ {-1, 0, 1} for the particular graph presented as inset. The results are sorted in decreasing order and the dashed line represents the value T unbiased obtained for β = (0, 0, 0, 0). (b) Non-isomorphic connected networks with N = 4. In this case, the nodes represent the multiple values that βi can take for the optimal configuration (with minimal T ), nodes shown as a pie chart represent βi ∈ {-1, 0, 1}, in other nodes, the bias βi is unique with its value included in the respective node.

In this manner, the pass from i to its neighbor j is independent of β i . This condition is observed in two particular nodes in G 5 in Fig. 3(b). In addition, if a node i has only one neighbor,

k i = 1, then w i→j ( β ) = A ij 1 βi N l=1 A il 1 βi = A ij = 1.
This particular case is observed in graphs G 1 , G 2 , G 3 in Fig. 3(b). On the other hand, in regular graphs all the nodes have the same degree k fulfilling the condition (14). As a consequence the random walk is not affected by the vector β. In Fig. 3(b), this is the case in G 4 (ring with k = 2) and G 6 (fully connected graph with k = 3).

complement the exploration of optimal configurations with the analysis of all the non-isomorphic connected networks with N = 5 in Fig. 4, the graphs G 1 , G 2 , . . . , G 21 were obtained from [37]. In Fig. 4(a), we present all the networks with different nodes representing the possible values of the bias β i in the optimal configurations. As a consequences of the condition ( 14) fulfilled for all the nodes, the dynamics on graphs G 1 , G 7 , G 14 and G 21 are not modified with the introduction of local bias. In other cases, optimal configurations occur with different biases and multiplicities represented with diverse types of nodes. In this respect, a pie chart denotes multiple values of β i for the bias minimizing the global time T . Here, it is worth noticing that only G 16 has a unique configuration for the optimal bias. On the other hand, in Fig. 4(b) we present a circular plot with the global times T optimal and T unbiased . Finally, in Fig. 4(c), we present the relative difference between T optimal and T unbiased as a percentage given by %∆ ≡ T unbiased -T optimal T unbiased × 100%.

We can see that in some cases, the improvement in the exploration of the network with the introduction of local bias is significant. For example, in G 6 we have %∆ = 9.82% whereas in G 13 , the optimal value is %∆ = 14.74% better.

To summarize, we have proved that it is possible to find a better strategy than the unbiased case by adding a local bias to the random walker strategy, in some cases the optimal configurations are non unique. In the next section we expand our results to more complex and larger structures.

IV. OPTIMAL RANDOM WALKS

In small-size networks, we saw that it is possible to explore all the values of the vector β and identify configurations with the optimal local bias. However, if at each node i, β i can take M values, the exploration of the complete domain β requires the evaluation of M N configurations. This exponential dependence with the size N makes the exploration of optimal configurations that minimize T a difficult task for N ≫ 1. In this section, we apply techniques for this minimization using random search algorithms developed to explore the configuration space efficiently to seek optimal ones. These type of algorithms have been applied on different branches of science like systems design [38], operational research [39], FIG. 4. Random walks with local bias on graphs with N = 5. (a) Multiplicity of β for the optimal configurations (with minimal T = T optimal ) for the G1, G2, . . . , G21 non-isomorphic connected networks with N = 5 nodes. (b) Values of the optimal time T optimal and for the unbiased random walk T unbiased obtained for β = (0, 0, 0, 0, 0). (c) Difference between T optimal and T unbiased expressed as the percentage %∆ defined in Eq. ( 15).

computational biology [START_REF] Cedersund | Uncertainty in Biology[END_REF], spin glasses and quantum algorithms [START_REF] Bapst | [END_REF]42], machine learning [START_REF] Goodfellow | Deep Learning[END_REF], among many others.

In the following part, we implement the simulated annealing (SA) algorithm, a Monte Carlo algorithm based on the physical process in which a solid is slowly cooled to reach a minimum energy configuration [START_REF] Van Laarhoven | Simulated annealing: Theory and applications[END_REF][START_REF] Berg | [END_REF][START_REF] Hansmann | The generalized-ensemble approach for protein folding simulations[END_REF]. The basic idea of SA is to use random moves to reach minimum values of a cost function, such as the global time T as a function of β, that physically plays the role of the system's energy with a transition probability given by the Boltzmann distribution law. At each move, we accept a new optimal configuration proposal β ′ with probability

ρ = min 1, exp - T ( β ′ ) -T ( β) T . ( 16
)
Here, T is a parameter that emulates the temperature and it is reduced at each Monte Carlo iteration in such a way that the probability of accepting a new configuration that does not favor the minimization of the cost function decreases, the values T ( β ′ ) and T ( β) are obtained numerically using Eqs. ( 8) and (10). The principal characteristic of SA is that it provides a form to escape from local minimums by occasionally accepting a configuration with a higher cost function. The algorithm has been widely used in a variety of problems such as clustering [START_REF] Brown | [END_REF], deep learning [48], quantum computing [49], chaotic systems [50], just to mention a few examples.

The convergence of the algorithm mostly depends on the type of temperature decrease as well as the way the configurations are modified [51]. In our implementation of the algorithm, the temperature decreases logarithmically, that is

T ℓ = T 0 ln(1 + ℓ) , (17)
where T ℓ is the temperature at iteration ℓ and T 0 is the initial temperature, we choose the initial configuration β with uniformly distributed random entries β i ∈ {-1, 0, 1}. Moreover, the variation of a given configuration is performed according to a trigger parameter α. Given an actual configuration β, for each iteration we generate a random number p ∈ (0, 1) uniformly distributed; then, if p < α, only two entries are modified randomly and for p ≥ α, half of its entries change randomly. This second alternative allows the algorithm to continue exploring the whole configuration space and not to get stuck in a local minimum when the temperature is quite small. Note that, when α is close to one, the modifications of β lie around the current configuration except for rare occasions where the change is abrupt. In addition, based on the results in Sec. III B, we only modify the entries of β assigned to nodes with degree k > 1 or where their nearest neighbors have different degrees. In the left panels, we present the values T for 100 realizations with 500 iterations of the algorithm implemented with α = 0.85 and T0 = 1. We show as an inset the results of the entropies found using Eq. ( 18) for the optimal configurations, the results are sorted in increasing order. In the right panels, we depict each network analyzed with their nodes representing the entropy value encoded in the color bar.

In this manner, we center the implementation of the SA algorithm to produce modifications in the global time T .

In Fig. 5 we illustrate the implementation of the SA algorithm to find optimal configurations of the local bias β reducing the global time T . We explore two networks: the Karate club network [52] with N = 34 presented in Fig. 5(a), and the Dolphins network [53] with N = 62 in Fig. 5(b). The results are obtained considering 100 Monte Carlo realizations with 500 iterations of the SA algorithm with α = 0.85 and T 0 = 1. Here, it is worth mentioning that in each realization, a different optimal value of T is obtained. This is due to the fact that SA is a heuristic method to get closer to a global minimum and obtain a similar result in each realization is unlikely if the configuration space is large where different minimums present multiplicities in the distribution of values β i as we illustrate in the analysis of small graphs in Figs. 34.

Moreover, the uncertainty of the multiple configurations can be quantified associating to each node i an entropy defined by

S i = - m∈{-1,0,1} f (i) m ln(f (i) m ), (18)
where f

(i)
m is the relative frequency of the values m that the i-th entry of β adopts in the final configurations obtained with the SA algorithm. This compact representation of our results allows us to determine how well defined is the local bias in a particular node. On the left panels in Fig. 5, we present as an inset the values of S i found for the final configurations produced with the SA algorithm. On the right, we show the respective networks with the nodes colored according to their entropy. In particular, nodes with S i = 0 are defined by a unique bias β i whereas S i = ln [START_REF] Michelitsch | Fractional Dynamics on Networks and Lattices[END_REF] indicates that β i can take with equal probability the values β i ∈ {-1, 0, 1}.

The analysis made in Fig. 5 can be implemented to other types of networks. In Table I, we expand our results to diverse graphs with different nature, topology, and size. In the first columns, we present a general description of the network features such as the size N , number of edges |E| = 1 2 l,m A lm and average degree k = 1 N N i=1 k i . We also include T unbiased for the unbiased case and the ensemble average T optimal for the optimal configurations found with the SA algorithm for 100 realizations with 500 iterations using α = 0.85 and T 0 = 1. The standard deviation of the results for all the realizations is denoted as σ T optimal . We also include the average of the entropy for all the nodes S = 1 N N i=1 S i and the respective standard deviation σ S for the values obtained in the nodes. In the last column, we present the relative difference %∆ defined in Eq. ( 15) to compare the unbiased case with the value T optimal . In Fig. 6 we present the networks analyzed, the color of each node i represents the normalized value S i /S max , with S max = ln(3).

The networks analyzed in Table I include two small size networks with 9 and 10 nodes, the Karate club and the Dolphins network analyzed in Fig. 5, three synthetic random networks generated with Erdős-Rényi model [33], the Watts-Strogatz algorithm [34] generated from a ring with nearest-neighbor and next-nearest-neighbor links and a rewiring probability of p = 0.05 and a complex network of the Barabási-Albert type, generated with the preferential attachment rule [35].

In addition, in different urban transportation networks, the main function of the infrastructure is to communicate efficiently to all the nodes. In some cases, models with random walkers can help us to characterize the mobility in a particular transportation mode [11,18,56] or to understand the effect of damage [START_REF] Eraso-Hernandez | Evolution of transport under cumulative damage in metro systems[END_REF]. Due to possible applications in the context of human mobility and urban transport, we analyze the metro networks in New York, Paris, and London. The networks were obtained from Refs. [START_REF] Derrible | Wolfram demonstrations project: World metro networks[END_REF][START_REF] Derrible | [END_REF], in this case, the nodes represent stations where users can change between lines or the end stations of a line (see details in Ref. [START_REF] Derrible | [END_REF]).

The results for all the networks explored show that the SA algorithm allows obtaining local bias configurations with T < T unbiased , the values for the average over realizations T optimal present small variations quantified with the standard deviation σ T optimal . Although the final results for T are similar, the local bias configurations can be very varied, something that is observed with the entropy values S averaged between all nodes, with standard deviations σ S that are high showing that there may exist nodes with well-defined local biases with their respective S i ≈ 0 or nodes where local bias does not affect dynamics with S i ≈ ln(3). Regardless of this multiplicity, it is always possible to improve the walker without bias, something that is revealed in the relative differences measured with %∆.

In this manner, the optimization of T obtained with the modification of a random walk strategy with local parameters is important because maintains the same adjacency matrix A without the introduction of new edges or rewiring of nodes, a fact that, in some systems may signify a huge increase of the operational costs; for example, the introduction of new lines in metro systems. All the results in this section reveal that beneficial improvements in the bias can be obtained through the implementation of the SA algorithm while the resulting average entropy gives us an idea of the multiplicity that these strategies can have.

V. CONCLUSIONS

In this research, we introduce a random walk strategy with a local bias. In this case, transition probabilities between nodes are defined in terms of values β i in each node, a vector β = (β 1 , β 2 , . . . , β N ) contains all the information of the biases. We explore the capacity of this random walker to visit the nodes of the network by using a global time T expressed in terms of eigenvalues and eigenvectors of the transition matrix defining the random walker.

Through the analysis of connected non-isomorphic graphs with N = 4, 5 nodes, we explore all the possible configurations of local bias for β i ∈ {-1, 0, 1}. The results reveal the existence of multiple configurations minimizing T . In different graphs, optimal bias improves the capacity of the walker to explore the network in comparison with the unbiased dynamics (standard random walk). In other cases, the introduction of local bias does not alter the dynamics. We deduce rules to identify nodes where the local bias does not change the elements of the transition matrix.

As the identification of optimal configurations is a problem with multiple minimal, we implement a simulated annealing algorithm to explore configurations close to an optimal value of T for diverse types of real and synthetic networks with different sizes and topologies. This approach allows deducing configurations for which local biases improve the dynamics in comparison with the unbiased random walker. The multiplicity of the local bias in these configurations is quantified in terms of an entropy obtained from the statistical analysis of the results produced by multiple realizations of the algorithm.

The framework introduced in this research is general and provides a tool for the exploration of different random walk dynamics defined in terms of parameters taking values at each node. The understanding of stochastic processes with this characteristic and its optimization can be useful in different contexts, for example in routing processes as well as the planning of routes in urban transportation systems.

  for l = 2, 3, . . . N , this is a consequence of the relation φ1 |φ ℓ = δ 1ℓ . In this manner τ j = i =j P ∞ i

FIG. 2 .

 2 FIG.2. Biased random walks on networks with N = 100 nodes: (a) Erdős-Rényi, (b) Watts-Strogatz, (c) Barabási-Albert network (details in the main text). We show the average global time T obtained numerically using Eq. (10) for degree-biased random walks defined by the transition probabilities in Eq. (11) with β = -2.5, -2.25, . . . , 2.25, 2.5 (continuous lines are used as a guide). The values with dots represent the results for the discrete values β = -1, 0, 1.

FIG. 5 .

 5 FIG. 5. Monte Carlo simulations of the SA algorithm to optimize the transport with a local bias. (a) The Karate club network, (b) the Dolphins network.In the left panels, we present the values T for 100 realizations with 500 iterations of the algorithm implemented with α = 0.85 and T0 = 1. We show as an inset the results of the entropies found using Eq. (18) for the optimal configurations, the results are sorted in increasing order. In the right panels, we depict each network analyzed with their nodes representing the entropy value encoded in the color bar.

TABLE I .

 I Characterization of the optimal local bias obtained through the SA algorithm in different network topologies. For each network we present the number nodes N , the total number of edges |E |, the average degree k, the global time T unbiased and the values T optimal , σT optimal , S, σS, %∆ characterizing the final results obtained with 100 Monte Carlo realizations of the SA algorithm implemented with 500 iterations, α = 0.85 and T0 = 1 (see details in the main text). FIG. 6. Networks analyzed in Table I. The nodes represent the value Si/Smax, with Smax = ln(3).

	Network	N |E |	k T unbiased T optimal σT optimal	S	σS %∆
	Network 9	9 13 2.889 11.3864 10.522	0.018 0.668 0.456 7.59
	Network 10	10 21 4.2 13.5057 11.4096	0.02 0.561 0.5 15.52
	Karate club [52]	34 78 4.588 63.916	58.554	0.488 0.794 0.334 8.39
	Dolphins [53]	62 159 5.129 155.231 129.479	2.273 0.756 0.325 16.59
	Erdős-Rényi [33]	50 69 2.76 143.055 131.97	1.827 0.875 0.253 7.75
	Watts-Strogatz [34]	50 100 4.0 137.887 132.909	0.085 0.656 0.408 3.61
	Barabási-Albert [35] 50 97 3.88 88.823	84.661	0.53 0.932 0.15 4.69
	Metro New York [54, 55] 77 109 2.831 280.02 260.369	4.894 0.685 0.374 7.02
	Metro Paris [54, 55]	78 125 3.205 239.196 227.816	2.479 0.714 0.404 4.76
	Metro London [54, 55] 83 121 2.916 329.138 311.495	4.709 0.721 0.386 5.36

ACKNOWLEDGMENTS

CSHC acknowledges support from CONACYT México. APR acknowledges financial support from by PAPIIT-UNAM grant No. IN116220.