{'fraction_non_alphanumeric': 0.05808198941277853, 'fraction_numerical': 0.013024744552505232, 'mean_word_length': 3.717853409222906, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 9, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Hitting times for second-order random walks

A second-order random walk on a graph or network is a random walk where transition probabilities depend not only on the present node but also on the previous one. A notable example is the non-backtracking random walk, where the walker is not allowed to revisit a node in one step. Second-order random walks can model physical diffusion phenomena in a more realistic way than traditional random walks and have been very successfully used in various network mining and machine learning settings. However, numerous questions are still open for this type of stochastic processes. In this work we extend well-known results concerning mean hitting and return times of standard random walks to the second-order case. In particular, we provide simple formulas that allow us to compute these numbers by solving suitable systems of linear equations. Moreover, by introducing the"pullback"first-order stochastic process of a second-order random walk, we provide second-order versions of the renowned Kac's and random target lemmas.


Introduction
Random walks on graphs are a key concept in network theory used to model various dynamics such as user navigation and epidemic spreading [29], to quantify node centrality [30,31], or to reveal communities and core-periphery structures [9,10,33,37]. The standard random walk on a graph considers a hypothetical walker that starting from node i moves to a new node of the graph choosing it uniformly at random among the neighbors of i. This defines a memoryless stochastic process which corresponds to a Markov chain on the nodes. However, the process described by a standard graph random walk is often localized and may fail to capture the underlying physical model as well as to fully exploit the global graph structure. For this reason, there is a growing interest in recent years towards so-called higherorder stochastic processes whose evolution may depend on past states, including stochastic processes modeled by transition probability tensors [5,12] as well as non-local diffusion mappings [7,11] and non-backtracking random walks [13].
The latter is one of the best known examples of second-order stochastic processes, where node transitions depend on both the current and the previous state spaces as the process is not allowed to backtrack, i.e. to immediately revisit the previous space in the next step. Non-backtracking random walks are often more appropriate than classical random walks to model real-world diffusive phenomena on networks where a message-passing or disease-spreading analogy is relevant. In particular, non-backtracking walks can avoid unwanted localization effects on the leading eigenvector of the adjacency matrix of a graph [17,20,24] and their mixing rate, in many cases, is faster than the one of classical random walks [27,6]. Furthermore, many computations with non-backtracking walks have negligible overheads compared to standard approaches [1,14].
The non-backtracking paradigm is at the basis of several second-order stochastic processes appearing in recent graph-theoretic, network mining and machine learning literature. For example, the mixing time of a clique-based version of the non-backtracking random walk on regular graphs is studied in [6]. Similarity measures based on second-order random walks are used in [38] to measure node proximity and improve classical walk-based link-prediction and clustering algorithms. Eigenvalues and eigenvectors of non-backtracking matrices have been shown to provide remarkable performances in community detection and in various centrality measurements [1,20,35]. The graph-embedding generated by the popular node2vec algorithm [15] is based on a flexible second-order random walk depending on two parameters: one used to control the likelihood of backtracking (i.e. immediately revisiting a node in the walk), the other used to bias the navigation towards or away from the immediately preceding nodes.
In this work we are concerned with hitting and return times of second-order random walks on graphs with possibly countably many vertices. Roughly speaking, a hitting time is the first time at which a stochastic process reaches a given subset of its state space, while a return time is the first time the process gets back to a given starting state. Hitting and return times are well understood for Markov chains, which include as a particular case the stochastic process that describes a standard random walk on a graph. In this case, the state space is given by the set of nodes and explicit formulas are available to compute mean hitting and return times by means of suitable systems linear equations which, for finite and connected graphs, have unique solutions. Moreover, these quantities have numerous applications to, for example, quantify node distances and similarities. For this reason, they are an important tool to evaluate network cohesion and node centrality and are at the basis of popular node classification and link prediction algorithms [22,28].
In this paper, we extend several results concerning mean hitting and return times of standard first-order stochastic processes to the second-order setting. In particular, we show that also mean hitting times of second-order random walks can be obtained from a set of linear equations, which has potential applications in higher-order classification and link prediction methods [3,4,36]. Moreover, we prove that mean return times of second-order random walks coincide with standard return times of suitably defined random walks in the original graph. In particular, mean return times for non-backtracking random walks coincide with the usual mean return times in all finite undirected graphs.
The paper is organized as follows. In Section 2 and 3 we present the notations used in this work and some basic results about random walks and second-order random walks, respectively. In Section 4 we introduce mean hitting times and return times for second-order random walks and we present first results concerning possible ways to compute them numerically via the solution of systems of linear equations. Our main results are in Section 5. There, we analyze second-order mean hitting and return times under what we call "equilibrium assumption", corresponding to the condition that the second-order random walk is at the steady state. In particular, we derive an algebraic relationship between the hitting time matrix for a second-order random walk and the one corresponding to a standard random walk on the directed line graph of the original graph. Moreover, we show second-order analogues of the Kac's and the random target lemmas for standard Markov chains. Finally, in Section 6, we illustrate some of the results of this paper on some real-world examples.

Notations and basic results
An oriented graph (or network) is a pair G = (V, E), where V is the vertex (or node) set and E is the set of oriented edges, i.e., ordered pairs of nodes. We say that (i, j) ∈ E is an edge going from i to j. From here onward, letters i, j, k denote generic nodes of G, while e, f are used to denote generic edges of G. Furthermore, • if e = (i, j) ∈ E then we set S(e) = i and T(e) = j, the source and terminal node of the edge e; • for i ∈ V we denote by I i = {e ∈ E : T(e) = i} and O i = {e ∈ E : S(e) = i} the in-neighborhood and out-neighborhood of i, respectively. The in-degree and the out-degree of node i are d − i = |O i | and d + i = |I i |, respectively.
In the present work, G can have a countably infinite number of nodes. However, we will assume that G is locally finite, that is, V is possibly infinite but each node has finite in-and out-degree.
A graph G = (V, E) can be completely described by its (possibly infinite dimensional) adjacency matrix A = (A ij ), defined as If the adjacency matrix is symmetric then G is called undirected. Hence, in an undirected graph each edge (i, j) has its own reciprocal edge (j, i). In that case we have d − i = d + i and this common value is the degree of node i, denoted as d i .
If for all i, j ∈ V there exists a walk from i to j, then the graph G is called strongly connected.
Finally, the symbols P and E denote the probability and the expectation, respectively. All vectors used in the text are column vectors, with possibly infinitely many entries. In particular, the symbol 1 denotes the column vector with all entries equal to 1, the size of which should be clear from the context.

Random walks
In this section, we recall from [18,32] some basic concepts and results in Markov chain theory. A random walk on G is a sequence of nodes v 0 , v 1 , . . . , v k , . . . where v k+1 is chosen at random between the out-neighbours of v k , according to specified probabilities p i,j : The last condition implies that d + i ≥ 1 for every i ∈ V , which will be assumed henceforth. A standard choice for p i,j is p i,j = 1/d + i , for j ∈ O i . In this case the random walk is called uniform. Our subsequent discussion is not restricted to this case.
A random walk can be viewed as a Markov chain {Y n , n ∈ N} with state space V and transition matrix P with entries P ij = p i,j . Clearly, P is a row-stochastic matrix. The random variable Y n represents the state of the walker at time n.
For a fixed subset S ⊂ V we call hitting time to S the random variable T S defined as T S = min{n ∈ N : Y n ∈ S}. In other words, T S is the first time the random walker reaches the subset S, starting from a given initial state Y 0 . Based on T S , we can define the hitting probability ϕ i→S from i to S as: i.e., the probability that, starting from node i, the random walker ever hits S. When ϕ i→S = 1, it is interesting to compute the average time taken by the random walker to reach S, starting from i. This quantity is given by and is called the mean hitting time from i to S. It is well-known that the quantities ϕ i→S and τ i→S can be obtained from the solution of certain linear systems, as shown by the next two theorems, whose proofs can be found in, e.g., [32]. Theorem 2.1. For a graph G = (V, E) and for a subset S ⊂ V , the hitting probabilities {ϕ i→S : i ∈ V } are the minimal non-negative solution of the following linear system: In the previous theorem, minimality of the solution means that if {x i : i ∈ V } is another non-negative solution of (2) then ϕ i→S ≤ x i for all i ∈ V . Theorem 2.2. For a graph G = (V, E) and for a subset S ⊂ V , suppose that ϕ i→S = 1 for all i ∈ V . Then the mean hitting times {τ i→S : i ∈ V } are the minimal non-negative solution of the following linear system: The case of finite graphs is the most usual in applications. In this case, a well-known result shows that for strongly connected graphs (3) has a unique solution, see e.g., Chap. VI of [18]. Precisely: Theorem 2.3. Let G = (V, E) be finite and strongly connected. For any given subset S ⊂ V , all the hitting probabilities to S are equal to 1 and the linear system (3) has a unique solution.
It is interesting to note that Equation (3) can be written in compact matrix-vector form as follows. Let t S and r S be the vectors entrywise defined as follows: where the numbers τ + i→S are defined as Then, the vector t S solves the singular equation In particular, among the infinitely many solutions x of this equation, the vector t S is the one characterized by the condition t S i = 0 for i ∈ S. The quantities (5) are called mean return times and represent the average number of steps required to return to S after leaving it from i ∈ S. Indeed, if we consider the random variable T + S = min{n ≥ 1 : Y n ∈ S}, then it holds τ + i→S = E(T + S |Y 0 = i). The special case when S is a singleton deserves a special attention, and we adopt the simplified notation This quantity is called mean return time to i. When G is finite and strongly connected, by means of the Perron-Frobenius theorem (see, e.g., [25]) we can conclude that the random walk {Y n , n ∈ N} admits a unique positive invariant distribution, i.e., a vector π > 0 such that π T 1 = 1 and π T P = π T . More generally, when V is countable, the same conclusion holds true if and only if the associated Markov chain is irreducible and positive recurrent, that is, τ i is finite for at least one i ∈ V , see e.g., Thm. 21.12 of [21] or Thm. 1.7.7 of [32]. The latter property also implies that all hitting probabilities (1) are equal to 1 and all hitting times are finite.
With the help of the invariant distribution we can obtain an explicit formula for the mean return time to a set S ⊂ V , This number quantifies the average number of steps taken by a random walker initially placed in S to hit again S, given that its probability of being in node i is proportional to π i . Indeed, the number π i / j∈S π j gives the conditional probability that the random walker is in i given that it is in S. The average return time τ S can be computed by means of Kac's lemma, see e.g., [16] or [21, Lemma 21.13]: Lemma 2.4. Let π be the stationary density of an irreducible Markov chain on the state space V . Then for any S ⊂ V , In particular, τ i = 1/π i for every i ∈ V .
If G is finite then the mean hitting times fulfill the so-called random target lemma, see e.g., Lemma 10.1 of [21]. Lemma 2.5. Let P be irreducible, and let π be its stationary density. Then, there exists a constant κ such that j∈V π j τ i→j = κ, The number κ = κ(P ) appearing in the preceding Lemma is the renowned Kemeny's constant of P . This constant quantifies the expected number of steps to get from node i to a node j, selected randomly according to the stationary density.
Introducing the matrix T = (τ i→j ), whose ij-th entry is the mean hitting time from i to j, the claim of Lemma 2.5 can be compactly expressed via the identity T π = κ1. Since the j-th column of T fulfills equation (6) with S = {j}, using Lemma 2.4 it is not difficult to verify that this matrix solves the equation where Diag(π) is the diagonal matrix with diagonal entries given by the entries of the vector π. Actually, the matrix T is the unique solution of (7) with zero diagonal entries.
We conclude this review on random walks with a technical lemma of independent interest and that will be useful for later. This result shows that, apart of a constant term, the vector t S defined in (4) can be expressed by a convex linear combination of the vectors {t i , i ∈ S}. Hence, mean access times to a subset can be computed from mean access times to the single elements of the set. Lemma 2.6. For any fixed subset S ⊂ V there exist positive coefficients β and {α i }, for i ∈ S, such that i∈S α i = 1, and In particular, Proof. Consider the vector a = (α i ) i∈V with entries By (6) we have Hence i α i = i∈S π i τ + i→S = 1. Moreover, from (7) we derive This proves that the difference T a − t S belongs to the kernel of I − P , that is, i∈S α i t i − t S = β1 for some β ∈ R. Considering the k-th entry of this identity for k ∈ S, we have and the proof is complete.

Second-order random walks
where the probability of choosing v k+1 depends not only on v k but also on v k−1 . More precisely, we describe this random walk by means of a stochastic process {X n , n ∈ N}, where X n ∈ V represents the node where the walker is located at time n. Hence, there exists constant probabilities p i,j,k such that, for n ≥ 0 We additionally assume that there are specific probabilities p i,j for the first transition, that is, A second-order random walk is not a Markov chain because it does not fulfill the Markov condition, as we need to remember the previous step in order to take the next step. However, it can be turned into a Markov chain by changing the state space from the vertices of the graph to the directed edges. More and only if T(e) = S(f ). Then, the second-order random walk can be seen as a (first-order) walk on G, that is, a sequence of directed edges e 0 , e 1 , . . . , e k , . . . where v 0 = S(e 0 ) and v i = S(e i ) = T(e i−1 ), for n ≥ 1. However, we must keep in mind that a walk of length m in G corresponds to a walk of length m − 1 in G.
For n ≥ 0 consider the random variable W n ∈ E defined as On the basis of (9) and (10), the sequence {W n } is a Markov chain with initial distribution P(W 0 = (i, j)) = P(X 0 = i)p i,j and transition matrix P given by Note that we allow the case P ef = 0 for some pair (e, f ) ∈ E, which may occur in some special circumstances, notably the non-backtracking random walk on G, see Example 3.2 below. Example 3.1. The uniform random walk on G is the Markov chain governed by the transition matrix P , otherwise.
As P ef do not depend on S(e), the corresponding stochastic process on G is the uniform random walk . A random walk on the Hashimoto graph is obviously related to a non-backtracking random walk on G. We can look at a random walk on H as a weighted random walk on the directed line graph G = ( V , E) by applying a null weight on the edges (e, f ) ∈ E such that T(f ) = S(e). The corresponding transition matrix P is defined as follows. If e = (i, j) and A is the adjacency matrix of G, then A directed edge e = (i, j) ∈ E such that O j = (j, i) is called dangling. The stochasticity condition f ∈E P ef = 1 implies that no edge in G is dangling, otherwise the out-degree of the corresponding node in H would be zero. Hence, the non-backtracking random walk on G is well defined if (and only if) d + i ≥ 1 for i ∈ V and there are no dangling edges. Finally, we recall from e.g., [19] that, if G is undirected then P is bi-stochastic. Example 3.3. Let P (1) and P (0) be the stochastic matrices defined in the Example 3.1 and Example 3.2, respectively. For any α ∈ (0, 1) the matrix P (α) = α P (1) + (1 − α) P (0) is stochastic (bi-stochastic, if G is undirected). The corresponding Markov chain describes random walks in G where backtracking steps are not completely eliminated, but rather the probability of performing one such step is downweighted by a factor depending on α. Thus, the corresponding second-order process can be called a backtrackdownweighted random walk. Combinatorial aspects of this kind of walks are addressed in, e.g., [2]. Remark 3.4. The Markov chain {W n } defined by (11) may not be irreducible even if G is strongly connected. For example, if G is an undirected graph consisting of a cycle, then the Hashimoto graph of G is not connected, and the matrix P built as in Example 3.2 is reducible. However, the matrix P (α) defined in Example 3.3 is irreducible in the sole assumption that G is strongly connected.

Second-order mean hitting times
In this section, we introduce the concept of mean hitting time for the second-order random walk {X n , n ∈ N} defined by (9) and (10) and we give results analogous to Theorems 2.1 and 2.2 for it. For each node k ∈ V , let T k be the random variable counting the first time the process arrives at node k, The aim of this section is to compute the corresponding mean second-order hitting time τ i→k = E( T k |X 0 = i), i.e., how many steps are required by a second-order random walk on average to reach k, starting from i. To this goal, we first compute the mean second-order hitting times to a node k, given both the first and the second node, and then we use these values to compute τ i→k by means of the formula The above formula shows that, unlike mean hitting times for Markov chains, second-order mean hitting times depend on the choice of the first transition probabilities, as defined in (10). In the sequel, we will see that there is a quite natural choice for these probabilities. On the other hand, the finiteness of the numbers τ i,j→k , which is essential for the finiteness of second-order hitting times, is not influenced by the values of the probabilities p ij . Indeed, in an infinite network it can happen that P( T k = +∞|X 1 = j, X 0 = i) > 0 for some (i, j) ∈ E and k ∈ V . For example, this happens when P is the transition matrix of the uniform random walk on G and this chain is transient. In this case, τ i,j→k = +∞ because τ i,j→k is the expectation of a random variable that is equal to +∞ with positive probability. For this reason, we first compute ϕ i,j→k = P( T k < +∞|X 1 = j, X 0 = i), the probability that, starting with X 0 = i and X 1 = j, the random walker reaches k in finite time. For these quantities, the following result holds: Theorem 4.1. For all k ∈ V , {ϕ i,j→k : (i, j) ∈ E} is the minimal non-negative solution of the following linear system: Proof. We first show that {ϕ i,j→k : (i, j) ∈ E} solves the linear system (15). Let us consider an edge (i, j) ∈ E. If i = k then X 0 = k, so T k = 0; if j = k then X 1 = k and T k = 1. In both cases it immediately follows that ϕ i,j→k = P( T k < +∞|X 1 = j, X 0 = i) = 1. In the other cases, we have that is, (15). Now let {x i,j : (i, j) ∈ E} be another non-negative solution of (15) for a fixed k. We show that x i,j ≥ ϕ i,j→k for all (i, j) ∈ E. Indeed, if i = k or j = k then x i,j = ϕ i,j→k and the claim follows. Otherwise, for i, j = k we have x i,j = p i,j,k + =k p i,j, x j,l . Let 0 = i and 1 = j. For anyn ≥ 2 we have where we used the inequality x t−1 t ≥ 0 and the empty sum (i.e., the square bracket in the n = 1 case) must be considered equal to 1. Observe that for n = 1, 2, 3, . . .
Hence we have x i,j ≥ P( T k ≤n + 1|X 1 = j, X 0 = i). Finally, and the claim follows.
On the basis of the preceding theorem we can compute the probability that the second-order process {X n } hits a node k ∈ V in a finite number of steps.
Corollary 4.2. For any fixed k ∈ V define ϕ i→k = P( T k < +∞|X 0 = i). Then, the vector (ϕ i→k ) i∈V is the vector (ϕ i,j→k ) (i,j)∈E premultiplied by the matrix M ∈ R |V |×|E| given by Proof. By elementary considerations, we have In a similar way, we can compute the second-order mean hitting times, as shown in the next theorem: Theorem 4.3. In the previous notation, if ϕ i,j→k = 1 for each (i, j) ∈ E, then { τ i,j→k : (i, j) ∈ E} is the minimal non-negative solution of the following linear system: Proof. We first show that { τ i,j→k : (i, j) ∈ E} is a solution of the linear system (17). Let (i, j) ∈ E. If i = k then P( T k = 0|X 1 = j, X 0 = i) = 1 and τ i,j→k = 0. If j = k then P( T k = 1|X 1 = j, X 0 = i) = 1 and τ i,j→k = 1. Otherwise, if both i = k and j = k then, by using the definition of expectation we obtain Now we observe that P( T k = 2|X 1 = j, X 0 = i) = P(X 2 = k|X 1 = j, X 0 = i) = p i,j,k . Then τ i,j→k is given by Finally, since m≥2 mP( T k = m|X 1 = , X 0 = j) = τ j, →k and m≥2 P( T k = m|X 1 = , X 0 = j) = 1 for = k, we have where in the last passages we have used the fact that ∈V p i,j, = 1. Now let {y i,j : (i, j) ∈ E} be another non-negative solution of (17). We show that y i,j ≥ τ i,j→k for all (i, j) ∈ E. If i = k or j = k then y i,j = τ i,j→k and the claim follows.
For i, j = k we have y i,j = 1 + p i,j,k + =k p i,j, y j,l . Now, let 0 = i and 1 = j. For anyn ≥ 2 we have where the empty sum (i.e., the square bracket when n = 1) must be considered equal to 1 and the last inequality holds since y t−1 , t ≥ 0. Rearranging terms, Finally, passing to the limitn → ∞ we obtain and the proof is complete.
Analogously to the argument of Corollary 4.2, we obtain the following characterization for the secondorder hitting times directly from Equation (14). Corollary 4.4. For any fixed k ∈ V , the vector ( τ i→k ) i∈V coincides with the vector ( τ i,j→k ) (i,j)∈E pre-multiplied by the matrix M in (16).

An alternative approach
The linear system in Theorem 4.3 allows us to compute second-order mean hitting times by means of the transition probabilities of the stochastic process {X n , n ∈ N}. The same quantities can be obtained using the fact that to any second-order random walk {X n } on G, corresponds a suitable random walk {W n } in G, as in (11). The two walks have the same transition probability, but the second is shorter than the first by one step. By correcting this discrepancy, we can obtain mean hitting times for {X n } from those for {W n }.
To this end, for any fixed subset S ⊂ E, consider the random variables T S = min{n ≥ 0 : W n ∈ S}, T + S = min{n ≥ 1 : W n ∈ S}, and the corresponding expectations These quantities can be computed by means of the techniques recalled in Section 2.1, since {W n } is a classical Markov chain. The theorem below shows how mean hitting (18) for {W n } relate to the solution of the linear system in (17). It follows that the problem of computing second-order mean hitting times in G can be expressed as the problem of computing traditional mean hitting times in G.
Let us consider a node k ∈ V . By Theorem 2.2, {τ e→I k : e ∈ E} is the minimal non-negative solution of the following linear system: Theorem 4.5. In the previous notation, it holds Proof. We show that each non-negative solution of (19) corresponds to a non-negative solution of (17) and vice-versa. Let {z i,j : (i, j) ∈ E} be a non-negative solution of (19). For all (i, j) ∈ E, we define y i,j as Then {y i,j : (i, j) ∈ E} is a non-negative solution of (17). In fact, if j = k but i = k then y i,j = 1. On the other hand, when i = k and j = k we have In the previous passages we used the fact that P is stochastic and P (i,j)( ,m) = 0 if = j. Conversely, if {y i,j : (i, j) ∈ E} is a non-negative solution of (17), then y i,j ≥ 1 if i = k. Thus, if we put, for all (i, j) ∈ E, then {z i,j : (i, j) ∈ E} is a non-negative solution of (19). The thesis follows by considering the minimal non-negative solution of both systems.

Second-order mean return times
In this section, we formulate the concept of return time for second-order random walks and we give an explicit formula for their computation in arbitrary networks, whenever the numbers τ i,j→k can be computed from (17). For a set S ⊂ V , consider the random variable T + S = min{n ≥ 1 : X n ∈ S}.
We are interested in computing the second-order mean return time to S, namely, τ S = E( T + S |X 0 ∈ S). This quantity represents the average number of steps that are required to reach again S, starting from a node in S and moving through the vertices of G according to a second-order random walk. For simplicity, in this section we limit ourselves to the single-node case S = {i}. We will go back to the general case in Section 5.1. Firstly, we observe that Introducing the auxiliary notation τ i,j = E( T + i |X 1 = j, X 0 = i), we obtain The intuition behind this formula is that the random walker chooses at first an out-neighbour j of i at random and then follows a second-order random walk that ends with node i. In fact, τ i,j can be easily computed by using the vector { τ i,j→k : (i, j) ∈ E} obtained from (17). Indeed, it holds Placing the formula above into (20), yields Apparently, this formula does not admit a simple matrix-vector form as in Corollary 4.4. In the sequel, we propose to adopt a specific choice for the numbers p i,j that allows us to obtain a result for any S ⊂ V analogous to Lemma 2.4.

The pullback of a second-order random walk
From now on, we suppose that the Markov chain associated to P is ergodic, so that there exists a (unique) vector π > 0 with π T P = π T . In this case, we can consider the additional assumption that the random walk {W n } on G is at equilibrium, i.e., P(W n = e) = π e , for n ≥ 0. Indeed, if the Markov chain {W n } is at the stationary state, then, owing to the correspondence between the events X n+1 = i and W n ∈ I i for n ≥ 0, the analysis of the process {X n } is simplified and reveals new properties.
For every e ∈ E, define The numbers λ e appearing in (21) have a simple probabilistic interpretation, namely λ e = P(W n = e|W n ∈ I T(e) ) is the conditional probability that the random walker in G is in e, given that they are also in I T(e) . Moreover, let L ∈ R |V |×|E| be the "lifting" matrix defined as follows: Note that the scalars λ e fulfill the condition e∈I i λ e = 1, i = 1, . . . , n.
This condition implies that L is stochastic, L1 = 1. If p is a probability vector on the nodes of G, then p T = p T L is a probability vector defined on the edges, and fulfills the identity p i = e∈I i p e . Furthermore, let R ∈ R |E|×|V | be the "restriction" matrix Given any probability distribution p on the edges of G, the product p T = p T R defines a probability distribution on the nodes such that p i = e∈I i p e . We note in passing that the product LR is the identity in R n . Theorem 5.1. Let P = L P R ∈ R n×n , Then, P is an irreducible stochastic matrix. Assuming that the Markov chain {W n } is at equilibrium with stationary density π, then also {X n } is at equilibrium with stationary density π T = π T R and, for n ≥ 1, P ij = P(X n+1 = j|X n = i). Moreover, the initial transition probabilities (10) are and the vector π T = π T R is the (unique) invariant vector of the Markov chain associated to P .
Proof. First, note that the matrix P is non-negative and stochastic. To prove that P is irreducible, let i, j ∈ V be arbitrary. Owing to the irreducibility of P , for any e ∈ I i and f ∈ I j there is a path in G from e to f . Let e = e 0 , e 1 , . . . , e m = f be such path. Then, the sequence T(e 0 ), T(e 1 ), . . . , T(e m ) is a path in G from i to j.
Thus, owing to the correspondence between X n+1 = i and W n ∈ I i , we have P(X n+1 = j|X n = i) = P(W n ∈ I j |W n−1 ∈ I i ) for any n ≥ 1. Moreover, consider the vector π T = π T R. In detail, To prove that π is an invariant vector of P , first observe that π T RL = π. Indeed, owing to (21), for any fixed e ∈ E we obtain ( π T RL) e = λ e f :T(f )=T(e) π f = π e .
Moreover, π T 1 = π T R1 = π T 1 = 1, eventually yielding To complete the proof it is sufficient to observe that, if {W n } is at equilibrium then also {X n } is at equilibrium and π (i,j) = P(W 0 = (i, j)) = P(X 1 = j, X 0 = i) Thus, p i,j = π (i,j) / f ∈I i π f . Finally, as f ∈I i π f = f ∈O i π f , we obtain (25), which concludes the proof.
We call the matrix P in (24) the pullback of P on G. Note that these two matrices correspond to two different stochastic processes, as briefly pointed out in Remark 5.2 below. Moreover, as we may expect, we remark that the pullback operation is not injective, as different stochastic processes on G may have the same pullback in G. This is shown by the Example 5.3 below, where the classic and non-backtracking processes are shown to have the same pullback. Remark 5.2. As shown in Theorem 5.1, the entries of the pullback matrix P correspond to first-order transition probabilities of the stochastic process {X n }. However, this does not mean that {X n } is a Markov chain, which would be true if (and only if) the sequence {X n } were to satisfy the Markov condition P(X n+1 = k|X n = j) = P(X n+1 = k|X n = j, X n−1 = i), for every (i, j), (j, k) ∈ E, at least under the considered assumptions. On the other hand, Equation (24) yields that is, P ij is a convex combinations of the numbers p h,i,j . Consequently, the Markov condition holds true if and only if the second-order transition probabilities (9) are independent of the oldest state, in which case {X n } is trivially a Markov chain.
Actually, the construction of P from P is reminiscent of a so-called lumping of {W n }, see e.g., [34] and [18,. Lumping methods are known to allow the construction of Markov chains with a reduced number of states from certain finite Markov chains with a larger states set. The corresponding transition matrices are related by an identity similar to the one defining P from P in Theorem 5.1. However, the process {X n } is not a lumping of {W n } since it is not a Markov chain, even when at equilibrium. This is the reason why second-order mean hitting times differ ostensibly from mean hitting times corresponding to P , as shown in various numerical experiments in Section 6.
Example 5.3. Let G be undirected and (strongly) connected, with adjacency matrix A, and let G be its line graph. Consider the random walk on G defined as in Example 3.1, with transition matrix P . It is not hard to verify that the pullback of P on G coincides with the transition matrix of the standard random walk on G. The same conclusion is true also for the matrix P considered in Example 3.2. Indeed, if A ij = 0 then because the condition S(f ) = T(e) is always false. On the other hand, if A ij = 1 then By linearity, also the pullback of the backtrack-downweighted random walk in Example 3.3 coincides with the standard random walk on G.

Hitting and return times at equilibrium
Due to Theorem 5.1, it is convenient to assume that the initial probabilities of the second-order stochastic process under consideration (10) are equal to λ (i,j) . Hence, from here onward we consider the "equilibrium assumption" P(X 0 = i) = π i together with (25). With this assumption, the presentation of the mean hitting and return times introduced in the previous section can be carried out with greater clarity and simplicity. To begin with, we derive an algebraic expression for the second-order hitting times matrix T = ( τ i→j ) of a very different nature from that at the basis of Equation (14). where L is as in (22) and R is as in (23).
Proof. Introduce the matrix C ∈ R |E|×|V | such that C ei = τ e→I i is defined as in Section 4.1. From the identity it follows that T = LC. Using the technical Lemma 2.6, for every e ∈ E and i ∈ V there exist coefficients α e and β i such that τ e→I i = f ∈I i α f τ e→f − β i . The latter identity can be expressed in matrix terms as Finally, as L1 = 1, we conclude.
The previous theorem yields a formula for the matrix T that is amenable to numerical computations, once the matrix T has been computed by standard methods via, e.g., (7). The next result proves that secondorder mean return times defined in Section 4.2 coincide with the usual return times of the random walk associated with the pullback of P on G.
Theorem 5.5. Let π be the stationary density of the pullback of P on G. For every S ⊂ V it holds τ S = 1/ i∈S π i . In particular, τ i = 1/π i for i ∈ V .
Proof. Let S = ∪ i∈S I i . Then, τ S can be reformulated as the mean return time in S of the chain {W n }. From Lemma 2.4 and (26), we get and we have the claim.
Example 5.6. Let G = (V, E) be finite and undirected. For α ∈ [0, 1], the stationary density π of the matrix P (α) in Example 3.3 is a constant vector, π e = 1/|E|, as the transition matrix is bi-stochastic. Consequently, the second-order mean return times coincides with the mean return times of the uniform random walk on G, namely, τ i = ( j d j )/d i , independently of α and the topology of G. Moreover, in this example Equation (25) becomes λ (i,j) = 1/d i being G undirected, i.e., the first step in the second-order process coincides with one step of a standard uniform random walk on G.
Theorem 5.5 yields a second-order analogous to Kac's lemma, see Lemma 2.4. For a family of networks that exhibit some specific symmetries, Lemma 2.5 can also be extended to second-order random walks, as shown in the sequel.
The coefficients α e appearing in the proof of Theorem 5.4 come from Lemma 2.6 and, according to (8), can be expressed as α e = π e τ + e→I i where i = T(e) and τ + e→I i represents the mean return time to I i of {W n }, after leaving it from e ∈ I i , as defined in (18). Suppose that these return times depend on i but not on e ∈ I i . This assumption is clearly verified if for every e, f ∈ I i there is an automorphism of G that exchanges e and f , i.e., the nodes e and f in G can be exchanged without altering the topology of G. The latter condition is satisfied if, for example, {W n } is a uniform random walk on a bipartite regular graph, that is, an undirected, bipartite graph whose adjacency matrix has the form in which the diagonal blocks can have different sizes.
For notational simplicity, denote by ρ i = τ + e→I i that common value. Then, from Theorem 5.1 we have Thus ρ i = 1/π i . Moreover, for e ∈ I i we obtain α e = π e ρ i = π e /π i . Hence, in the notation of Theorem 5.4, we have that Diag(a)Rπ = π, and thus where b = (β i ) i∈V is the coefficient vector defined in Theorem 5.4. By Lemma 2.5, we have T π = κ1 where κ is the Kemeny's constant of T , and we eventually obtain: Combining all the observations above, we obtain the following second-order equivalent of the random target lemma: before the removal after the removal nodes edges diameter nodes edges diameter Guppy  99  726  6  98  725  5  Dolphins  62  159  8  53  150  7  Householder93  104  211  7  73  180  5   Table 1: Dimension of the networks used in our experiments before and after the removal of leaf nodes. Corollary 5.7. Let π be the stationary density given in Theorem 5.1. In the equilibrium assumption discussed at the beginning of Section 5.1 and using the previous notation, we have that, if for every i ∈ V the return time τ + e→I i depends on i but not on e ∈ I i , then j∈V π j τ i→j = κ, 6 Numerical experiments on real-world networks In this section, we present numerical results on the computation of second-order mean hitting times on some real-world networks. We consider three undirected and unweighted networks: Guppy [8], that represents the social interactions of a population of free-ranging guppies (Poecilia reticulata), Dolphins [23], that represents the social structure of a group of dolphins, and Householder93 [26], that represents a collaboration network of a group of mathematicians. Since these graphs contain dangling edges, we remove in sequence all the nodes of degree 1 until each vertex has degree greater or equal to 2. Table  1 shows the dimension of these networks before and after removing such leaf nodes. Our experiments illustrate various statistics on the second-order hitting times τ i→j for non-backtracking random walks in these graphs. Note that hitting times computed with the pullback transition matrix coincide with the hitting times τ i→j for the standard random walk, see Example 5.6.
In Figure 1 we scatter plot the average non-backtracking mean hitting times m j = i∈V τ i→j /|V | against the average traditional mean hitting times m j = i∈V τ i→j /|V |. The dotted red line is the bisector of the first quadrant. We observe that, in all these networks, the ratio between m j and m j is always smaller than 1. Intuitively, this means that, on average, non-backtracking random walks are shorter than classical random walks.   i→j denote second-order mean hitting times associated with the transition matrix P (α) . To better illustrate the effect of penalizing backtracking steps, we measure the mean hitting times ratio The curves in Figure 3 illustrate the behavior of the average value (magenta line), along with the maximum and the minimum value (black lines) of these ratios. The values for α = 0 correspond to the results shown in Figure 1.

Conclusions
Second-order random walks, which include non-backtracking random walks as a particular case, provide stochastic models for navigation and diffusion processes on networks that can offer improved modeling capacity with respect to classical random-walks, and thus allow us to better capture real-world dynamics on networks. In this paper, we propose a very general framework to properly define mean hitting times and mean return times for second-order random walks for directed graphs with finite or countably many nodes. Our approach is based on the correspondence between the second-order random walk and a standard random walk obtained by first lifting the second-random walk to the directed line graph and then pulling back onto the original graph. The "pullback" Markov chain obtained this way has the same stationary density as the initial second-order random process and allow us to transfer several well-known results for standard random walks to the second-order case. For example, we prove that second-order mean return times coincide with the mean return times of the corresponding pullback Markov chain and we show how to extend the Kac's and random target lemmas to second-order random walks. These results prompt us to further investigations concerning computational issues as well as theoretical properties of second-order random walks. For instance, it is certainly interesting to devise more efficient numerical methods to compute second-order hitting times than the one based on Corollary 4.4 and Theorem 5.4. Furthermore, numerical experiments show that non-backtracking hitting times are consistently smaller than their corresponding classical, backtracking versions. In fact, reducing the probability of backtracking results in a reduction of the mean hitting times. However, a formal explanation of this fact, at the best of our knowledge, is still missing.