Add Bernoulli random graph

src/Graphs.jl

@@ -295,6 +295,8 @@ export
     kronecker,
     dorogovtsev_mendes,
     random_orientation_dag,
+    bernoulli_graph,
+    rho_correlated_bernoulli_graphs,
 
     # community
     modularity,

src/SimpleGraphs/SimpleGraphs.jl

@@ -74,6 +74,8 @@ export AbstractSimpleGraph,
     static_scale_free,
     kronecker,
     random_orientation_dag,
+    bernoulli_graph,
+    rho_correlated_bernoulli_graphs,
     # generators
     complete_graph,
     star_graph,

src/SimpleGraphs/generators/randgraphs.jl

@@ -1422,3 +1422,54 @@ function random_orientation_dag(
     end
     return g2
 end
+
+"""
+    bernoulli_graph(Λ; rng = nothing, nodes_type = Int64)
+
+Given the symmetric matrix ``\\Lambda \\in [0,1]^{n \\times n}``, return a Bernoulli graph with ``n`` vertices. Each edge ``(i,j)`` exists with probability ``\\Lambda[i,j]``.
+"""
+function bernoulli_graph(
+    Λ::AbstractMatrix{<:AbstractFloat};
+    rng::Union{Nothing,AbstractRNG}=nothing,
+    nodes_type::Type=Int64,
+)
+    size(Λ)[1] != size(Λ)[2] &&
+        throw(ArgumentError("The probability matrix must be a square matrix"))
+    !issymmetric(Λ) && throw(ArgumentError("Λ must be symmetric"))
+    n = size(Λ)[1]
+    g = SimpleGraph{nodes_type}(n)
+    for j in 1:n
+        for i in (j + 1):n
+            if rand(rng) <= Λ[i, j]
+                add_edge!(g, i, j)
+            end
+        end
+    end
+    return g
+end
+
+"""
+rho_correlated_bernoulli_graphs(Λ, ρ; rng = nothing, nodes_type = Int64)
+
+Given the symmetric matrix ``\\Lambda \\in [0,1]^{n \\times n}`` and a real number ``\\rho \\in [0,1]`` return a Tuple with two ``\\rho``-correlated Bernoulli graphs with ``n`` vertices. It means that, calling ``A`` and ``B`` the adjacency matrix of the outputed graphs, for ``i,j`` such that ``i``!=``j`` the Pearson correlation coefficient for ``A_{i,j}`` and ``B_{i,j}`` is ``\\rho``.
+"""
+function rho_correlated_bernoulli_graphs(
+    Λ::AbstractMatrix{<:AbstractFloat},
+    ρ::Float64;
+    rng::Union{Nothing,AbstractRNG}=nothing,
+    nodes_type::Type=Int64,
+)
+    (ρ < 0.0 || ρ > 1.0) && throw(ArgumentError("ρ must be in [0,1]"))
+    n = size(Λ)[1]
+    g2 = SimpleGraph{nodes_type}(n)
+    g1 = bernoulli_graph(Λ; rng, nodes_type=eltype(g2))
+    g1_adj = Int.(adjacency_matrix(g1))
+    for j in 1:n
+        for i in (j + 1):n
+            if rand(rng) <= ((1 - ρ) * Λ[i, j] + ρ * g1_adj[i, j])
+                add_edge!(g2, i, j)
+            end
+        end
+    end
+    return g1, g2
+end

test/runtests.jl

@@ -21,7 +21,7 @@ const testdir = dirname(@__FILE__)
 end
 
 @testset verbose = true "Code formatting (JuliaFormatter.jl)" begin
-    @test format(Graphs; verbose=false, overwrite=false, ignore=["vf2.jl"])  # TODO: remove ignore kwarg once the file is formatted correctly
+    @test format(Graphs; verbose=false, overwrite=false, ignore="vf2.jl")  # TODO: remove ignore kwarg once the file is formatted correctly
 end
 
 @testset verbose = true "Doctests (Documenter.jl)" begin

test/simplegraphs/generators/randgraphs.jl

@@ -384,4 +384,33 @@
         @test μ1 - sv1 <= 0.3 * 5 <= μ1 + sv1 # since the stdev of μ1 is around sv1/sqrt(N), this should rarely fail
         @test μ2 - sv2 <= 0.7 * 3 <= μ2 + sv2
     end
+
+    @testset "bernoulli graphs" begin
+        n = 50
+        Λ = sparse(rand(n, n))
+        @test_throws ArgumentError("Λ must be symmetric") bernoulli_graph(Λ; rng)
+        @test_throws ArgumentError("The probability matrix must be a square matrix") bernoulli_graph(
+            sparse(rand(2, 3)); rng
+        )
+        Λ = Symmetric(Λ)
+        @test_throws ArgumentError("ρ must be in [0,1]") rho_correlated_bernoulli_graphs(
+            Λ, -1.0; rng
+        )
+        ρ = 1.0 # isomorphism case
+        g1, g2 = rho_correlated_bernoulli_graphs(Λ, ρ; rng)
+        g1_adj = adjacency_matrix(g1)
+        g2_adj = adjacency_matrix(g2)
+        @test g1_adj == g2_adj
+        @test diag(g1_adj) == diag(g2_adj) == zeros(n)
+        ρ = 0.5 # non isomorphism case
+        g3, g4 = rho_correlated_bernoulli_graphs(Λ, ρ; rng)
+        g3_adj = adjacency_matrix(g3)
+        g4_adj = adjacency_matrix(g4)
+        @test g3_adj != g4_adj
+        @test diag(g3_adj) == diag(g4_adj) == zeros(n)
+        g5 = bernoulli_graph(Λ; rng)
+        for g in testgraphs(g5)
+            @test nv(g) == n
+        end
+    end
 end