Add newman_watts_strogatz graph generator (#198)

In a later paper by Newman and Watts they generate networks that are
slightly different from the original small-world networks by Watts and
Strogatz. Instead of rewiring edges, they keep every original edge and
add new long-range connections along the network.

Author: blwh

src/Graphs.jl

@@ -278,6 +278,7 @@ export
     erdos_renyi,
     expected_degree_graph,
     watts_strogatz,
+    newman_watts_strogatz,
     random_regular_graph,
     random_regular_digraph,
     random_configuration_model,

src/SimpleGraphs/SimpleGraphs.jl

@@ -57,6 +57,7 @@ export AbstractSimpleGraph,
     erdos_renyi,
     expected_degree_graph,
     watts_strogatz,
+    newman_watts_strogatz,
     random_regular_graph,
     random_regular_digraph,
     random_configuration_model,

src/SimpleGraphs/generators/randgraphs.jl

@@ -276,11 +276,11 @@ end
     watts_strogatz(n, k, β)
 
 Return a [Watts-Strogatz](https://en.wikipedia.org/wiki/Watts_and_Strogatz_model)
-small world random graph with `n` vertices, each with expected degree `k`
+small world random graph with `n` vertices, each with expected degree `k` 
 (or `k - 1` if `k` is odd). Edges are randomized per the model based on probability `β`.
 
 The algorithm proceeds as follows. First, a perfect 1-lattice is constructed,
-where each vertex has exactly `div(k, 2)` neighbors on each side (i.e., `k` or
+where each vertex has exacly `div(k, 2)` neighbors on each side (i.e., `k` or
 `k - 1` in total). Then the following steps are repeated for a hop length `i` of
 `1` through `div(k, 2)`.
 
@@ -296,6 +296,7 @@ For `β = 0`, the graph will remain a 1-lattice, and for `β = 1`, all edges wil
 be rewired randomly.
 
 ### Optional Arguments
+- `remove_edges=true`: if false, do not remove the original edges.
 - `is_directed=false`: if true, return a directed graph.
 - `rng=nothing`: set the Random Number Generator.
 - `seed=nothing`: set the RNG seed.
@@ -317,7 +318,8 @@ function watts_strogatz(
     n::Integer,
     k::Integer,
     β::Real;
-    is_directed=false,
+    is_directed::Bool=false,
+    remove_edges::Bool=true,
     rng::Union{Nothing,AbstractRNG}=nothing,
     seed::Union{Nothing,Integer}=nothing,
 )
@@ -363,14 +365,66 @@ function watts_strogatz(
             d == s && continue          # Self-loops prohibited
             d == t && break             # Rewired to original target
             if add_edge!(g, s, d)       # Was this valid (i.e., unconnected)?
-                rem_edge!(g, s, t)      # True rewiring: Delete original edge
-                break                   # We found a valid target
+                remove_edges && rem_edge!(g, s, t)     # True rewiring: Delete original edge
+                break                                   # We found a valid target
             end
         end
     end
     return g
 end
 
+"""
+    newman_watts_strogatz(n, k, β)
+
+Return a Newman-Watts-Strogatz small world random graph with `n` vertices, each
+with expected degree `k(1 + β)` (or `(k - 1)(1 + β)` if `k` is odd). Edges are
+randomized per the model based on probability `β`.
+
+The algorithm proceeds as follows. First, a perfect 1-lattice is constructed,
+where each vertex has exacly `div(k, 2)` neighbors on each side (i.e., `k` or
+`k - 1` in total). Then the following steps are repeated for a hop length `i` of
+`1` through `div(k, 2)`.
+
+1. Consider each vertex `s` in turn, along with the edge to its `i`th nearest
+   neighbor `t`, in a clockwise sense.
+
+2. Generate a uniformly random number `r`. If `r < β`, `s` is connected to some
+   vertex `d`, chosen uniformly at random from the entire graph, excluding `s` and
+   its neighbors. (Note that `t` is a valid candidate.)
+
+For `β = 0`, the graph will remain a 1-lattice, and for `β = 1`, all edges will
+be rewired randomly.
+
+Note: In the original paper by Newman and Watts, self-loops and double edges were
+allowed, which is not the case here. However, for large enough networks and small
+enough `p` and `k`, this should not deviate much from the original model.
+
+### Optional Arguments
+- `is_directed=false`: if true, return a directed graph.
+- `rng=nothing`: set the Random Number Generator.
+
+## Examples
+```jldoctest
+julia> newman_watts_strogatz(10, 4, 0.3)
+{10, 26} undirected simple Int64 graph
+
+julia> newman_watts_strogatz(Int8(10), 4, 0.8, is_directed=true, seed=123)
+{10, 36} directed simple Int8 graph
+```
+
+### References
+- Scaling and percolation in the small-world network model, M. E. J. Newman, Duncan J. Watts. [https://doi.org/10.1103/PhysRevE.60.7332](https://doi.org/10.1103/PhysRevE.60.7332)
+"""
+function newman_watts_strogatz(
+    n::Integer,
+    k::Integer,
+    β::Real;
+    is_directed::Bool=false,
+    rng::Union{Nothing,AbstractRNG}=nothing,
+)
+    return watts_strogatz(n, k, β; is_directed=is_directed, remove_edges=false, rng=rng)
+end
+
 function _suitable(edges::Set{SimpleEdge{T}}, potential_edges::Dict{T,T}) where {T<:Integer}
     isempty(potential_edges) && return true
     list = keys(potential_edges)

test/simplegraphs/generators/randgraphs.jl

@@ -86,6 +86,33 @@
         @test is_directed(ws) == true
     end
 
+    @testset "Newman-Watts-Strogatz" begin
+        # Each iteration creates, on average, n*k/2*(1+β) edges (if the network is large enough).
+        # As such, we need to average to check if the function works.
+        nsamp = 100
+        n = 1000
+        β = 0.2
+        k = 4
+        mean_num_edges = n * k / 2 * (1 + β)
+        nes = 0
+        for i in 1:nsamp
+            nws = newman_watts_strogatz(n, k, β; rng=rng)
+            nes += ne(nws)
+            @test nv(nws) == n
+            @test is_directed(nws) == false
+        end
+        @test abs(sum(nes) / nsamp - mean_num_edges) / mean_num_edges < 0.01
+
+        nes = 0
+        for i in 1:nsamp
+            nws = newman_watts_strogatz(n, k, β; is_directed=true, rng=rng)
+            nes += ne(nws)
+            @test nv(nws) == n
+            @test is_directed(nws) == true
+        end
+        @test abs(sum(nes) / nsamp - mean_num_edges) / mean_num_edges < 0.01
+    end
+
     @testset "Barabasi-Albert" begin
         ba = barabasi_albert(10, 2; rng=rng)
         @test nv(ba) == 10