Directed random geometric graphs: structural and spectral properties

In this work we analyze structural and spectral properties of a model of directed random geometric graphs: given n vertices uniformly and independently distributed on the unit square, a directed edge is set between two vertices if their distance is smaller than the connection radius ℓ , which is randomly drawn from a Pareto distribution. This Pareto distribution is characterized by the power-law decay α and the lower bound of its support ℓ0 ; thus the graphs depend on three parameters G(n,α,ℓ0) . By increasing ℓ0 , for fixed (n,α) , the model transits from isolated vertices ( ℓ0≈0 ) to complete graphs ( ℓ0=2 ). We first propose a phenomenological expression for the average degree ⟨k(G)⟩ which works well for α > 3, when k is a self-averaging quantity. Then we numerically demonstrate that 〈Vx(G)〉≈n[1−exp(−〈k〉] , for all α, where Vx(G) is the number of nonisolated vertices of G. Finally, we explore the spectral properties of G(n,α,ℓ0) by the use of adjacency matrices represented by diluted random matrix ensembles; a non-Hermitian and a Hermitian one. We find that ⟨k⟩ is a good scaling parameter of spectral and eigenvector properties of G mainly for large α.


Introduction
Graphs and networks are nowadays well accepted tools to model complex systems. Moreover, if the system under study is embedded in a geometric space (like road networks, power networks or networks of wireless communication devices) the best suited representation belongs to the so-called spatial networks [1]. In this respect, probably the most studied model of spatial networks is the random geometric graph (RGG) model [2,3] introduced by Gilbert in 1959 [4].
In Gilbert's original version a RGG consists of n vertices uniformly and independently distributed on the infinite plane, where two vertices u and v are connected by an edge if their Euclidean distance is less or equal than the connection radius ℓ [4,5]. Moreover, since Gilbert's pioneering work, several variants of RGGs have been proposed to account for the properties of realistic complex systems: Commonly, RGGs are considered inside bounded domains such as the unit square [2,6], the unit circle, see e.g. [7], or d-dimensional cubes and spheres [8]. Also, recently, RGGs over more intricate bounding domains (i.e. elongated rectangles [9][10][11] or annuli [12,13]) have been studied. RGGs with non-uniform vertex distributions have been explored too [14][15][16]; while quasirandom vertex distributions have also been considered in [17]. Additionally, motivated by wireless ad hoc networks, there is a generalization of RGGs, known as soft RGGs [18], where the connection between vertices obeys a probabilistic function; see also [12,19].
It is interesting though that even when most of the real-world spatial networks that could be modeled by RGGs (i.e. transportation networks like cargo ship networks, infrastructure networks like power networks or road networks, communication networks, etc) and the processes associated to them (i.e. disease spreading, urban sprawl, commuter behavior, etc) are intrinsically directed, models of directed RGGs have been scarcely studied. For exceptions see [20,21] where models of directed RGGs have been used to study wireless networks and word association networks, respectively. Therefore, to contribute to fill this gap, in this work we numerically explore structural and spectral properties of a model of directed RGGs.

Random graph model
In this work, we consider the following model of directed RGGs (dRGGs): given n vertices uniformly and independently distributed on the unit square we randomly assign the connection radius ℓ u to the vertex u, where ℓ u is taken from a Pareto distribution, Here, α > 0 is the power-law decay of the Pareto distribution, ℓ 0 ∈ (0, √ 2] is the lower bound of the support of the distribution, and is a normalization factor such that´√ 2 ℓ0 f(ℓ)dℓ = 1. In figure 1 we show some examples of probability distribution functions f(ℓ) for different combinations of α and ℓ 0 . Then, if the vertex v lies within the disc of radius ℓ u , centered at u, there is a directed edge from u to v, u → v. Note that since ℓ u ̸ = ℓ v , the existence of u → v does not imply u ← v, so the graph is directed. Therefore, this model of dRGGs depends on three parameters: G(n, α, ℓ 0 ).
It is fair to mention that the random graph model G(n, α, ℓ 0 ) we study here is inspired by that introduced in [21], with three important differences: (a) we consider vertices distributed on the unit square instead of on a two dimensional torus. Indeed, in [21] the graph model is more general since vertices are considered on a d-dimensional torus. (b) We consider α > 0 instead of α > d + 1, so here we can approach the case of a flat f(ℓ) when α → 0. (c) We consider ℓ 0 as a free parameter in the range (0, √ 2], while in [21] it is set to ℓ 0 = (log n/πn) 1/2 (for the two dimensional torus), such that the resulting graph is almost surely connected. So here we can explore the limit ℓ 0 → 0 producing graphs with mostly isolated vertices. Thus, the model of dRGGs defined above shows a transition by increasing ℓ 0 , for fixed (n, α), from isolated vertices when ℓ 0 ≈ 0 to complete graphs for ℓ 0 = √ 2. However, large values of α also enhance the existence of isolated vertices; see for example the black distribution in figure 1(c) which is strongly peaked at ℓ ≈ 0.
In section 4 we explore structural and spectral properties of dRGGs. In particular, we analyze the spectral properties of the model by the use of a recently introduced approach under which the adjacency matrices of random graphs are represented by diluted random matrix theory (RMT) ensembles 1 . In this respect, Erdös-Rényi graphs and RGGs have been represented by diluted versions of the Gaussian orthogonal ensemble (GOE) [11,24], multilayer and multiplex networks have been studied through diluted-banded random matrix ensembles [25], while bipartite and mutualistic networks were analyzed by the use of block-diluted random matrix ensembles [26]. Here, within the same approach, we use two adjacency matrix representations that we define as follows.

Randomly weighted non-Hermitian adjacency matrix
We define the elements of the non-Hermitian adjacency matrix A nH of the dRGGs as where ϵ uv are independent random variables drawn from a normal distribution with zero mean and unit variance; ϵ uv ̸ = ϵ vu since G(n, α, ℓ 0 ) is directed. According to this definition, diagonal random matrices (known in RMT as the Poisson ensemble (PE) [27]) are obtained for ℓ 0 → 0, when the graphs consist mostly of disconnected vertices; whereas the real Ginibre ensemble (RGE) (i.e. full real and non-symmetric random matrices [28]) is recovered for ℓ 0 = √ 2, when the graphs are complete. Therefore, we expect to observe a transition from the PE to the RGE for increasing ℓ 0 . Thus, the ensemble defined by A nH corresponds to a diluted RGE.

Randomly weighted Hermitian adjacency matrix
Additionally, given the standard binary adjacency matrix B extracted from the directed graph G(n, α, ℓ 0 ), we consider the following Hermitian adjacency matrix A H : Again, ϵ uv are independent random variables taken from a normal distribution with zero mean and variance one. It is clear that In fact, A H , with the name of magnetic adjacency matrix, was already used in [29] to represent randomly-weighted Erdös-Rényi directed graphs; while the unweighted version of A H was introduced in [30]. The main advantage of A H over A nH is that the former has a real spectrum, thus some spectral measures from RMT can be computed easier.
Therefore, the magnetic random matrix ensemble (MRME) defined by A H transits for increasing ℓ 0 from the PE, when ℓ 0 ≈ 0, to the GOE (i.e. full real and symmetric random matrices [27]), when ℓ 0 = √ 2. However, note that the MRME can not be considered as a diluted GOE since for any ℓ 0 < √ 2, A H is a complex Hermitian matrix.

Structural measures
In addition to the standard quantities used to characterize the structure of random graphs and networks, such as the statistical properties of the vertex degree k and the average number of nonisolated vertices V x , here we also use topological indices.
Given a simple graph G = (V(G), E(G)) with V(G) and E(G) the vertex set and edge set, respectively, the Randić connectivity index is defined as [31] while the Harmonic index is given by [32] Here uv denotes the edge connecting vertices u and v and k u the degree of vertex u: where k in u and k out v are the in-degree and out-degree of vertex u. From definitions (5) and (6), when ℓ 0 ≈ 0 (i.e. when the model of dRGGs produces mostly isolated vertices), R(G) ≈ 0 and H(G) ≈ 0; while for ℓ 0 = √ 2 (i.e. when the graph model produces complete graphs), R(G) = N/2 and H(G) = N/2. That is, when applying R(G) and H(G) to ensembles of random graphs characterized by the parameter set (n, α, ℓ 0 ), we expect to observe a transition from ⟨R(G)⟩ ≈ 0 and ⟨H(G)⟩ ≈ 0 to ⟨R(G)⟩ = N/2 and ⟨H(G)⟩ = N/2 when increasing ℓ 0 from ℓ 0 ≈ 0 to √ 2. We stress that the random weights we impose to the adjacency matrices A nH and A H , as defined respectively in equations (3) and (4), do not play any role in the computation of the degree-based topological indices.
In fact, within a statistical approach to random graphs, it has been recently shown that the average values of R(G) and H(G), normalized to the order of the graph n, scale with the average degree ⟨k⟩; see e.g. [33][34][35]. That is, ⟨R(G)⟩ /n and ⟨H(G)⟩ /n are functions of ⟨k⟩ only. Moreover, it was also found that ⟨R(G)⟩ and ⟨H(G)⟩ are highly correlated with the Shannon entropy of the eigenvectors of the adjacency matrix of RGGs [36]. This is a notable result because it puts forward the application of degree-based topological indices beyond mathematical chemistry (it is relevant to add that random graphs have also been studied by means of eigenvalue-based topological indices such as the Estrada index, the Laplacian Estrada index, and Rodriguez-Velazquez indices, see e.g. [37,38]). Therefore, the study of R(G) and H(G) on dRGGs seems to be pertinent.

Spectral measures
We use standard RMT measures to characterize spectral and eigenvector properties of the non-Hermitian A nH and Hermitian A H adjacency matrices defined, respectively, in equations (3) and (4).
Regarding spectral properties, for A nH which has complex spectra, we use the ratio between nearest-and next-to-nearest-neighbor spacings r C , with the ith ratio defined as [39] where λ NN i and λ NNN i are, respectively, the nearest and the next-to-nearest neighbors of λ i in C. While, given the real and ordered spectrum of A H , λ 1 > λ 2 > . . . > λ n , we compute the ratio between consecutive level spacings r R , where the ith ratio is given by [40] .
Note that r C can also be computed for real spectra, so it works for both A nH and A H . Regarding eigenvector properties, given the normalized eigenvectors Ψ i (i.e. ∑ n m=1 |Ψ i m | 2 = 1) of A nH or A H , we compute the inverse participation ratios [41] ] −1 (10) and the Shannon entropies [42] Both IPR and S measure the extension of eigenvectors on a given basis, so in the case of random graphs they can be used to prove percolation transitions; see e.g. [24,36]. From definitions (8)-(11), when ℓ 0 ≈ 0 (i.e. when A nH and A H are mostly diagonal), when A nH and A H are full random matrices), all spectral and eigenvector measures are expected to acquire their corresponding full random matrix values.

Structural properties
We start by exploring some statistical properties of the vertex degree k of the model of dRGGs. In figure 2 we present examples of degree distributions P(k) for different graph sizes n (different columns correspond to different values of n) and two values of α: α = 1.5, panels (a-d), and α = 7, panels (i-l). We note that, for fixed α, the support of P(k) grows with n while keeping the shape; indeed the shape of P(k) evolves with α. Specifically, while for small α, P(k) develops a two-peak structure (one broad peak for k → 1 and one sharp peak at k = n − 1, see e.g. panels (a-d) for α = 1.5), P(k) is unimodal for large α (see e.g. panels (i-l) for α = 7). Moreover, we observe a smooth transition in the shape of P(k) from the two-peak structure to the Notice, though, that in figures 3(a)-(e) we are plotting distributions of the normalized degree k = (k − ⟨k⟩) /σ(k), where σ 2 (k) is the variance of k, to show that the shape of P(k) depends on α only.
Also, in figure 2 we plot the average degree ⟨k⟩ as a function of ℓ 0 for dRGGs characterized by α = 1.5, panels (e)-(h), and α = 7, panels (m)-(p). In particular, we noticed that the curves ⟨k⟩ vs. ℓ 0 for α = 7 can be well described by the expression for ⟨k(n, ℓ)⟩, derived in [9], for standard RGGs when replacing ℓ by ℓ 0 : see the magenta dashed lines in figures 2(m)-(p). This fact makes sense to us because, for large α, f(ℓ) is strongly peaked at ℓ 0 so most radii get the value of ℓ 0 . Moreover, taking into account that we are dealing with a distribution of radii ℓ, we found that describes the average degree of dRGGs with α > 3 even better than equation (12). In equation (13), C i are fitting constants and ⟨ℓ⟩ is the average connection radius that we computed directly from equation (1) as where η is given in equation (2). An important consequence of equation (13)    the numerical data for α > 3, while for smaller α we see good correspondence for large ℓ 0 only. Notice form table 1 that the value of the fitting constant C 1 is approximately equal to π for all α; therefore, the prediction for ⟨k(n, ℓ)⟩ for standard RGGs from [9] with ℓ → ⟨ℓ⟩ approximately coincides with equation (13) to the leading order.
To complete the statistical inspection of k, in figures 3(k)-(o), we show that the variance of k grows with ℓ 0 , approaches a maximum, and then decreases to zero for ℓ → √ 2; where the maximum is higher the larger the value of α is. Also, in figures 3(p)-(t), we plot the ratio σ 2 (k)/⟨k⟩ 2 as a function of ℓ 0 . It is interesting to note that for small α, σ 2 (k)/⟨k⟩ 2 remains constant for increasing graph size (that is, the curves in panels (p, q) fall one on top of the other for different n), meaning that k is not a self-averaging quantity. On the contrary, for large α and small ℓ 0 , the ratio σ 2 (k)/⟨k⟩ 2 decreases with n (see the left part of the curves in panels (s, t)); so, under these conditions, k is a self-averaging quantity. Now, we recall that increasing ℓ 0 from ℓ 0 ≈ 0 to ℓ 0 = √ 2 drives the model of dRGGs from mostly isolated vertices to complete graphs. This transition, as a function of ℓ 0 , could be well characterized by computing the average number of nonisolated vertices ⟨V x (G)⟩ as well as the average values of the Randić and the Harmonic indices, see e.g. [36]. Certainly, we expect ⟨V x (G)⟩ ≈ 0 when ℓ 0 ≈ 0 and ⟨V x (G)⟩ = n for ℓ 0 =  dRGGs characterized by α = 3. We show curves for four different graph sizes n in each panel. As can be clearly seen, the three quantities ⟨X(G)⟩ (where X represents V x , R or H) undergo a smooth transition (in log scale) from ⟨X(G)⟩ ≈ 0 for small ℓ 0 to the corresponding maxima max[X(G)] at ℓ 0 = √ 2. Moreover, the normalized curves ⟨X(G)⟩/max[X(G)] vs. ℓ 0 have a very similar functional form but they are displaced to the left on the ℓ 0 -axis for increasing n; see figures 4(b), (e) and (h). In fact it was shown, for standard RGGs [36] and also for non-uniform RGGs [14], that these three normalized quantities scale with ⟨k⟩; that is, the curves ⟨X(G)⟩/max[X(G)] vs. ⟨k⟩ fall one on top of the other for different graph sizes. Therefore, we validate this scaling for dRGGs in figures 4(c), (f) and (i) where we plot ⟨V x (G)⟩/n, ⟨R(G)⟩/(n/2) and ⟨H(G)⟩/(n/2) vs. ⟨k⟩ in log-log scale (main panels) as well as in semilog scale (insets). While we observe a very good scaling for the average number of nonisolated vertices, the scaling of both the Randić and the Harmonic indices is good for ⟨k⟩ < 1 only.
Furthermore, in figure 5 we explore the scaling of ⟨X(G)⟩/max[X(G)] with ⟨k⟩ for dRGGs characterized by values of α ranging from 1.5 to 7. We do confirm that ⟨V x (G)⟩/n shows a perfect scaling with ⟨k⟩ for all the values of α examined in this work; see figures 5(a)-(e). In contrast, the scaling of ⟨R(G)⟩/n and ⟨H(G)⟩/n with ⟨k⟩ is reasonably good for α > 3 only; indeed for α < 3 these two indices do not scale in any range of ⟨k⟩, see figures 5(f), (g), (k) and (l). We surmise that the lack of scaling of ⟨R(G)⟩/n and ⟨H(G)⟩/n with ⟨k⟩, for small α, may be related to the non self-averaging property of k for those values of α.
We also recall that an expression for V x (G) is known for standard RGGs [43]:  (16). Each column corresponds to a different value of the power-law decay α. All averages were computed from 10 6 /n random graphs.
which was recently generalized for both RGGs and non-uniform RGGs as [14] ⟨V Then, given the perfect scaling of ⟨V x (G)⟩/n with ⟨k⟩ (meaning that ⟨V x (G)⟩/n is a function of ⟨k⟩ only), we test equation (16) on dRGGs by plotting it in figures 5(a)-(e); see the magenta dahsed lines. Indeed, we do confirm that equation (16) describes perfectly well our numerical data for all α. This is a remarkable result because it stresses the universality of equation (16) on RGGs models.

Spectral properties
We now explore some spectral properties of dRGGs. We perform exact numerical diagonalization to obtain the eigenvalues λ i and eigenvectors Ψ i (i = 1, . . . , n) of large ensembles of randomly-weighted adjacency matrices A nH and A H characterized by the parameter set (n, α, ℓ 0 ). In figures 6(a), (d) and (g) we present, respectively, the average ratio between consecutive eigenvalue spacings ⟨r C (G)⟩, the average inverse participation ratio ⟨IPR(G)⟩, and the average Shannon entropy ⟨S(G)⟩ as a function of ℓ 0 for dRGGs characterized by α = 3. Note that here the dRGGs are represented by the diluted RGE; i.e. by the non-Hermitian adjacency matrices A nH of equation (3). In each panel we report results for graphs of four different sizes n. As for the structural measures (see figures 4(a), (d) and (g) as a reference), the three spectral measures ⟨X(G)⟩ (where X represents r C , IPR or S) undergo a smooth transition (in log scale) from ⟨X(G)⟩ ≈ X PE for ℓ 0 ≈ 0 to the corresponding maxima max[X(G)] = X RGE at ℓ 0 = √ 2. Here, (r C ) PE ≈ 0.5 [29], (r C ) RGE ≈ 0.737 [29], IPR PE = 1, IPR RGE ≈ n/2.04, S PE = 0, and S RGE ≈ ln(n/1.56) [29]. We note that since we did not find an expression for IPR RGE in the literature we computed it numerically for the graph sizes used in this work.
To ease the analysis of these average quantities we conveniently normalize them as In contrast, we note that the scaling of ⟨IPR(G)⟩ is not good, see figure 6(f). However, note that in figure 6 we are reporting results for dRGGs characterized by α = 3 only. Then, in figure 7 we explore the scaling of ⟨X(G)⟩ with ⟨k⟩ for dRGGs characterized by several values of α. As for the topological indices reported in figure 5, here we observe better scaling of spectral and eigenvector properties the larger the value of α is, even for ⟨IPR(G)⟩; see figures 7(f)-(j).
Recall that the results in figures 6 and 7 are for dRGGs represented by the non-Hermitian adjacency matrices A nH of equation (3). The question now is whether the Hermitian adjacency matrices A H of equation (4) provide equivalent results. Then in figure 8 we plot the average normalized ratios ⟨r C (G)⟩ and ⟨r R (G)⟩ (recall that r C can be computed for real spectra too) for dRGGs represented by the MRME and  characterized by several values of the power-law decay α. In figure 8 we have used the following normalizations: We just want to add that since the MRME reproduces the GOE only when the corresponding graphs are complete, the values for spectral and eigenvector measures have different values for ℓ → √ 2 and exactly at ℓ = √ 2, as shown in figure 8 where we can see an abrupt decay of the curves at ℓ = √ 2. Specifically, while (r C ) MRME(ℓ→ √ 2) ≈ 0.6175 and (r R ) MRME(ℓ→ √ 2) ≈ 0.5995 [36], (r C ) MRME(ℓ= √ 2) = (r C ) GOE ≈ 0.5688 [36] and (r R ) MRME(ℓ= √ 2) = (r R ) GOE ≈ 0.5307 [40]. From figure 8 it is clear that both ratios, ⟨r C (G)⟩ and ⟨r R (G)⟩, scale better with ⟨k⟩ the larger the value of α is. However, even for small α (see figure 8(a), (f) corresponding to α = 1.5) the scaling is reasonably good. This is in contrast with ⟨r C (G)⟩ when computed on A nH , see figure 7(a), where the curves of ⟨r C (G)⟩ for different graph sizes n are clearly displaced on the ⟨k⟩-axis. We have also verified (see the appendix) that the eigenvector properties of the MRME scale with ⟨k⟩ better than those computed for the diluted RGE.

Discussion and conclusions
We have analyzed structural and spectral properties of a model of dRGGs. The graphs of this model depend on three parameters G(n, α, ℓ 0 ): the number of vertices n uniformly and independently distributed on the unit square, the power-law decay α > 0 of a Pareto distribution (which characterizes the random connection radii ℓ used to define directed edges among the vertices), and the lower bound of the support ℓ ∈ [ℓ 0 , √ 2] of the Pareto distribution. For a fixed pair (n, α), by increasing ℓ 0 this model of dRGGs transits from isolated vertices, when ℓ 0 ≈ 0, to complete graphs, for ℓ 0 = √ 2; so we computed structural and spectral properties of dRGGs along this transition.
First, we explored some structural properties of the model. We discovered that for α > 3 the degree k(G) (defined as the sum of the in-degree and out-degree, see equation (7)) is a self-averaging quantity, while it is not for α < 3. Then, we were able to propose a phenomenological expression for the average degree ⟨k(G)⟩, see equation (13) which indeed fails when k is a non-self-averaging quantity. Moreover, we demonstrated that the average number of nonisolated vertices of G, normalized to the graph size, ⟨V x (G)⟩/n, perfectly scales with ⟨k⟩, see figures 5(a)-(e); following a quite simple relation, see equation (16). In fact, equation (16) allow us to write down an explicit expression for ⟨V x (G)⟩ in terms of the model parameters (n, α > 3, ℓ 0 ) by the use of equations (2), (13) and (14): We also computed topological indices (in particular the Randić index R(G) and the Harmonic index H(G)) on dRGGs. In fact, we found that ⟨k⟩ scales properly ⟨R(G)⟩/(n/2) and ⟨H(G)⟩/(n/2) when k is a self-averaging quantity. Second, we explored some spectral properties of the model by means of standard random matrix theory measures such as the ratio between nearest-and next-to-nearest-neighbor spacings r C , the ratio between consecutive level spacings r R , the inverse participation ratios of eigenvectors IPR, and the Shannon entropies of eigenvectors S. We computed these measures on randomly-weighted adjacency matrices of G(n, α, ℓ 0 ) represented by a non-Hermitian and a Hermitian random matrix ensemble; see equations (3) and (4) for A nH and A H , respectively. While the ensemble defined by A nH corresponds to a diluted RGE, we named the ensemble defined by A H as the MRME. Even though, both random matrix ensembles provide similar results, we found that the scaling of the spectral and eigenvector properties of G with ⟨k⟩ is observed in a wider range of α for the MRME; see the appendix. This puts forward a deeper study of Hermitian adjacency matrices to represent directed graphs.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.   In addition, it is useful to look at the values of the exponent γ of equation (A.1) obtained from the fittings of the curves ⟨k⟩ * vs. n, which are reported in tables A1-A3 for ⟨r C (G)⟩, ⟨IPR(G)⟩ and ⟨S(G)⟩, respectively. From table A1 we can clearly see, for all the values of α reported here, that the values of γ for the MRME are smaller than those for the diluted RGE. Moreover, since γ ≈ 0 for the MRME, we can claim that ⟨k⟩ is indeed the scaling parameter of ⟨r C (G)⟩ for this ensemble of adjacency matrices, as was already observed in figure 8. The situation is quite different for ⟨IPR(G)⟩ and ⟨S(G)⟩, see tables A2 and A3: there the values of γ for the MRME are smaller than those for the diluted RGE only when α ⩾ 3. This means that only for α ⩾ 3 (when k is a self-averaging quantity) the eigenvector measures scale better with ⟨k⟩ for the MRME than for the diluted RGE. Anyway, the scaling of the eigenvector measures of both random adjacency matrix ensembles seems to be acceptable for α ⩾ 3 only, where the curves ⟨X(G)⟩ vs. ξ are clearly invariant (i.e. they fall one on top of the other for different n). So, a deeper study of the statistical properties of k, in particular the self-averaging property, seems to be necessary.