Note that this code has been tested with Julia 1.5.1 and SimpleHypergraphs.jl version v0.1.13
.
using Pkg
#pkg"add SimpleHypergraphs"
using SimpleHypergraphs
The GraphPlot
package will make it possible to vizualize classic graph representations of a hypergraph such as bisection or 2-section representation.
# pkg"add GraphPlot"
The second option is to use for visualisation the Python's HyperNetX. In order for this integration to work you need to install PyCall.jl
and Conda.jl
along with appropiate Python modules.
#pkg"add PyCall Conda"
#using PyCall
#using Conda
#Conda.runconda(`install matplotlib --yes`)
#Conda.runconda(`install networkx --yes`)
#run(`$(PyCall.python) -m pip install hypernetx==1.2.5`)`
# since the output of this command is long it is not included in this notebook
Include the library package with using.
using SimpleHypergraphs
Usually SimpleHypergraphs is used together with standard LightGraphs.jl library
import Graphs
We start by creating an empty hyperaph with 5 vertices and 4 hyperedges
h = Hypergraph{Float64}(5,4)
5×4 Hypergraph{Float64, Nothing, Nothing, Dict{Int64, Float64}}: nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing nothing
Now we add the connection weights
h[1:3,1] .= 1.5
h[3,4] = 2.5
h[2,3] = 3.5
h[4,3:4] .= 4.5
h[5,4] = 5.5
h[5,2] = 6.5
h
5×4 Hypergraph{Float64, Nothing, Nothing, Dict{Int64, Float64}}: 1.5 nothing nothing nothing 1.5 nothing 3.5 nothing 1.5 nothing nothing 2.5 nothing nothing 4.5 4.5 nothing 6.5 nothing 5.5
#basic operations over the hg h
@assert add_vertex!(h) == 6
@assert add_hyperedge!(h) == 5
h[5,5] = 1.2
h[6,5] = 1.3
h
6×5 Hypergraph{Float64, Nothing, Nothing, Dict{Int64, Float64}}: 1.5 nothing nothing nothing nothing 1.5 nothing 3.5 nothing nothing 1.5 nothing nothing 2.5 nothing nothing nothing 4.5 4.5 nothing nothing 6.5 nothing 5.5 1.2 nothing nothing nothing nothing 1.3
To visualize a given hypergraph h
, the user needs to specify two mandatory parameters:
h
to drawh
GraphBased
represents each hyperedge he
with a fake vertex fv
to which each vertex v ∈ he
is connected.HyperNetX
renders an Euler diagram of the hypergraph where vertices are black dots and hyper edges are convex shapes containing the vertices belonging to the edge set.GraphBased
visualization¶with_node_labels=true
, but node_labels
is not specified, vertex ids will be used as their label.SimpleHypergraphs.draw(h,
GraphBased;
width=500,
height=500,
radius=10, #same radius for each node
node_color = "yellow", #same color for each node
node_stroke="orange", #same stroke for each node
stroke_width=2, #same stroke-width value for each node
node_opacity=0.5, #same opacity for each node
with_node_labels=true, #wheter displaying or not node labels
with_node_metadata_hover=true,
)
SimpleHypergraphs.draw(
h,
GraphBased;
width=500,
height=500,
radius=10, #same radius for each node
node_color = "yellow", #same color for each node
node_colors = ["yellow", "yellow", "yellow", "blue", "red", "red", "blue"],
node_stroke = "orange", #same stroke for each node
node_strokes = ["orange", "orange", "orange", "orange", "black", "black", "black"],
stroke_width=2, #same stroke-width value for each node
node_opacity=0.5, #same opacity for each node
with_node_labels=true, #whether displaying or not node labels
node_labels=["A","B","C","D","E","F","G"],
with_node_metadata_hover=true,
)
with_node_weight=true
, each vertex weight within the hyperedges it belongs to will be displayed.SimpleHypergraphs.draw(
h,
GraphBased;
width=500,
height=500,
radius=10, #same radius for each node
node_color = "yellow", #same color for each node
node_stroke="orange", #same stroke for each node
stroke_width=2, #same stroke-width value for each node
node_opacity=0.5, #same opacity for each node
with_node_labels=true, #whether displaying or not node labels
node_labels=["A","B","C","D","E","F","G"],
with_node_metadata_hover=true,
with_node_weight=true
)
draw(
h,
GraphBased;
width=500,
height=500,
radius=10, #same radius for each node
node_color = "yellow", #same color for each node
node_stroke="orange", #same stroke for each node
stroke_width=2, #same stroke-width value for each node
node_opacity=0.5, #same opacity for each node
with_node_labels=true, #whether displaying or not node labels
with_node_metadata_hover=true,
with_node_weight=true, #whether displaying vertices metadata on mouse hover
he_colors=["green", "blue", "red", "yellow","black"], #hyperedges colors
with_he_labels=true, #whether displaying or not hyperedge labels
he_labels=["a","b","c","d"], #hyperedges labels
with_he_metadata_hover=true #whether displaying hyperedges metadata on mouse hover
)
SimpleHypergraphs integates the Python library HyperNetX to let the user visualize a hypergraph h
exploiting an Euler-diagram visualization. For more details, please refer to the library HyperNetX.
draw(h, HyperNetX; width=5, height=5, no_border=true)
There are many options for Hypergraph
plotting. Type ?draw
to see them all.
?draw # press Ctrl+Enter to see documentation for `draw`
search: draw RoundNearestTiesAway dirname isdirpath deg2rad @deprecate
function draw(
h::Hypergraph,
type::Type{GraphBased};
element::Union{String, Int}=get_next_div_id(),
width::Int=500,
height::Int=500,
radius::Real=10,
node_radii::Union{AbstractVector{<:Real}, Nothing}=nothing,
node_color::String="#999",
node_colors::Union{AbstractVector{String}, Nothing}=nothing,
node_stroke::Union{String, Nothing} = nothing,
node_strokes::Union{AbstractVector{String}, Nothing}=nothing,
stroke_width::Real=0,
stroke_widths::Union{AbstractVector{<:Real}, Nothing}=nothing,
node_opacity::Real=1,
node_opacities::Union{AbstractVector{<:Real}, Nothing}=nothing,
stroke_opacity::Real=1,
stroke_opacities::Union{AbstractVector{<:Real}, Nothing}=nothing,
with_node_labels::Bool=false,
node_labels::Union{AbstractVector{String}, Nothing}=nothing,
with_node_metadata_hover::Bool=false,
with_node_weight::Bool=false,
he_colors::Union{AbstractVector{String}, Nothing}=nothing,
with_he_labels::Bool=false,
he_labels::Union{AbstractVector{String}, Nothing}=nothing,
with_he_metadata_hover::Bool=false
)
Draw a hypergraph h
in a web-based environment (e.g. Jupyter Notebook), using a js script based on the library (D3)[https://d3js.org/]. Each hyperedge he
is represented by a fake vertex fv
to which each vertex v ∈ he
is connected.
Arguments
h
: the hypergraph to drawtype
: how the hypergraph will be drawn. If type=GraphBased
, each hyperedgewill be represented as a vertex (see above)
width
: width of the figureheight
: height of the figureradius
: same default radius for each vertex (represented as a circle)node_radii
: distinct radius values for each vertexnode_color
: same default color for each vertexnode_colors
: distinct node colors for each vertexnode_stroke
: same default stroke for each vertexnode_strokes
: distinct node strokes for each vertexstroke_width
: same default stroke-width for each vertexstroke_widths
: distinct stroke-width values for each vertexnode_opacity
: same default opacity for each vertexnode_opacities
: distinct node-opacity values for each vertexstroke_opacity
: same default stroke-opacity for each vertexstroke_opacities
: distinct stroke-opacity values for each vertexwith_node_labels
: whether displaying node labelsnode_labels
: node labels to be shownwith_node_metadata_hover
: whether displaying node metadata when hovering each vertexwith_node_weight
: whether displaying node weights within each hyperedgehe_colors
: distinct hyperedge colors for each hyperedgewith_he_labels
: whether displaying hyoeredges labelswith_he_metadata_hover
: whether displaying hyperedge metadata when hovering each hyperedgedraw(
h::Hypergraph,
type::Type{HyperNetX};
width::Int=10,
height::Int=10,
node_labels::Union{Dict{Int, String}, Nothing}=nothing,
edge_labels::Union{Dict{Int, String}, Nothing}=nothing,
collapse_nodes::Bool=false,
collapse_edges::Bool=false,
pos::Union{Dict{Int,Pair{Int,Int}}, Nothing}=nothing,
with_color::Bool=true,
with_node_counts::Bool=false,
with_edge_counts::Bool=false,
layout::PyObject=nx.spring_layout,
layout_kwargs::Dict=Dict{String, Any}(),
ax::Union{PyObject, Nothing}=nothing,
no_border::Bool=false,
edges_kwargs::Dict=Dict{String, Any}(),
nodes_kwargs::Dict=Dict{String, Any}(),
edge_labels_kwargs::Dict=Dict{String, Any}(),
node_labels_kwargs::Dict=Dict{String, Any}(),
with_edge_labels::Bool=true,
with_node_labels::Bool=true,
label_alpha::Float64=.35
)
Draw a hypergraph h
as an Euler diagram, using the library HyperNetX.
Arguments
h
: the hypergraph to drawtype
: how the hypergraph will be drawn. If type=HyperNetX
, the hypergraph will be represented as a Euler Diagramwidth
: width of the figureheight
: height of the figurenode_labels
: node labels to be shownedge_labels
: edge labels to be showncollapse_nodes
: draws the hypergraph gotten by identifying nodes contained by the same edges (from HyperNetX)collapse_edges
: draws the hypergraph gotten by identifying edges containing the same nodes (from HyperNetX)no_border
: indicates wheter the figure should have a borderFor more details about the other parameters, please refer to the library HyperNetX.
The type BipartiteView
represents a non-materialized view of a bipartite representation hypergraph h
. Note this is a view - changes to the original hypergraph will be automatically reflected in the view.
The bipartite view of a hypergraph is suitable for processing with the LightGraphs.jl
package.
Several LightGraphs methods are provided for the compability.
b = BipartiteView(h)
{11, 11} undirected simple Int64 graph
The BipartiteView
provide LightGraphs.jl compability.
supertype(typeof(b))
Graphs.SimpleGraphs.AbstractSimpleGraph{Int64}
We add here a edge to a parent Hypergraph of a bisection view. Note that this change will be reflected in the bipartite view
add_vertex!(h)
7
This graph can be plotted using LightGraphs
tools.
using GraphPlot
using Graphs
nodes, hyperedges = size(h)
nodes_membership = fill(1, nodes)
hyperedges_membership = fill(2, hyperedges)
membership = vcat(nodes_membership, hyperedges_membership)
nodecolor = ["lightseagreen", "orange"]
#membership color
nodefillc = nodecolor[membership]
gplot(b, nodefillc=nodefillc, nodelabel=1:nv(b), layout=circular_layout)
The functionality of LightGraphs
can be used directly on a bipartite view of a hypergraph.
Graphs.a_star(b, 1, 3)
2-element Vector{Graphs.SimpleGraphs.SimpleEdge{Int64}}: Edge 1 => 8 Edge 8 => 3
#number of vertices
nv(b)
12
#number of edges
ne(b)
11
#neighbors
sort(collect(outneighbors(b,5)))
3-element Vector{Int64}: 9 11 12
#neighbors
sort(collect(inneighbors(b,9)))
1-element Vector{Int64}: 5
#shortest path - it does not consider the nodes associated with a hyperedge
shortest_path(b,1,4)
3-element Vector{Int64}: 1 2 4
Represents a two section view of a hypergraph h
. Note this is a view - changes to the original hypergraph will be automatically reflected in the view.
The bipartite view of a hypergraph is suitable for processing with the LightGraphs.jl
package.
Several LightGraphs methods are provided for the compability.
Note that the view will only work correctly for hypergraphs not having overlapping hyperedges. To check whether a graph has overlapping edges try has_overlapping_hedges(h) - for such graph you need to fully materialize it rather than use a view. This can be achieved via the get_twosection_adjacency_mx(h) method.
# This condition is required for an unmaterialized `TwoSectionView` representation of a hypergraph to make sense
@assert SimpleHypergraphs.has_overlapping_hedges(h) == false
t = TwoSectionView(h)
{7, 8} undirected simple Int64 graph
gplot(t, nodelabel=1:nv(t))
#number of vertices
nv(t)
7
#number of edges
ne(t)
8
#neighbors
sort(collect(outneighbors(t,5)))
3-element Vector{Int64}: 3 4 6
#neighbors
sort(collect(inneighbors(t,1)))
2-element Vector{Int64}: 2 3
#shortest path
shortest_path(t,1,5)
3-element Vector{Int64}: 1 3 5
Let us consider the following hypergraph
h = Hypergraph{Float64}(8,7)
h[1:3,1] .= 1.5
h[3,4] = 2.5
h[2,3] = 3.5
h[4,3:4] .= 4.5
h[5,4] = 5.5
h[5,2] = 6.5
h[5,5] = 5.5
h[5,6] = 6.5
h[6,7] = 5.5
h[7,7] = 6.5
h[8,7] = 6.5
h[8,6] = 6.5
h
8×7 Hypergraph{Float64, Nothing, Nothing, Dict{Int64, Float64}}: 1.5 nothing nothing nothing nothing nothing nothing 1.5 nothing 3.5 nothing nothing nothing nothing 1.5 nothing nothing 2.5 nothing nothing nothing nothing nothing 4.5 4.5 nothing nothing nothing nothing 6.5 nothing 5.5 5.5 6.5 nothing nothing nothing nothing nothing nothing nothing 5.5 nothing nothing nothing nothing nothing nothing 6.5 nothing nothing nothing nothing nothing 6.5 6.5
Let us search for communities in the hypergraph h
best_comm = findcommunities(h, CFModularityCNMLike(100))
display(best_comm.bm)
display(best_comm.bp)
0.3193650793650794
2-element Vector{Set{Int64}}: Set([4, 2, 3, 1]) Set([5, 6, 7, 8])
And now we visualize them in 2-section view
t = TwoSectionView(h)
function get_color(i, bp)
color = ["red","green","blue","yellow"]
for j in 1:length(bp)
if i in bp[j]
return color[j]
end
end
return "black"
end
gplot(t, nodelabel=1:nv(t), nodefillc=get_color.(1:nv(t), Ref(best_comm.bp) ))