Improve the performance and add options to the plot methods

Project.toml

@@ -1,7 +1,7 @@
 name = "EcologicalNetworksPlots"
 uuid = "9f7a259d-73a7-556d-a7a2-3eb122d3865b"
 authors = ["Timothée Poisot <timothee.poisot@umontreal.ca>"]
-version = "0.1.0"
+version = "0.1.1"
 
 [deps]
 EcologicalNetworks = "f03a62fe-f8ab-5b77-a061-bb599b765229"

docs/make.jl

@@ -1,4 +1,4 @@
-push!(LOAD_PATH,"../src/")
+push!(LOAD_PATH, "../src/")
 
 using Documenter, EcologicalNetworksPlots
 
@@ -19,12 +19,12 @@ const pages = [
 ]
 
 makedocs(
-    sitename = "EcologicalNetworksPlots",
-    authors = "Timothée Poisot",
-    pages = pages
-    )
+    sitename="EcologicalNetworksPlots",
+    authors="Timothée Poisot",
+    pages=pages
+)
 
 deploydocs(
-    repo = "github.com/EcoJulia/EcologicalNetworksPlots.jl.git",
-    push_preview = true
+    repo="github.com/EcoJulia/EcologicalNetworksPlots.jl.git",
+    push_preview=true
 )

docs/src/advanced/attributes.md

@@ -28,6 +28,21 @@ plot(I, Unes, aspectratio=1)
 scatter!(I, Unes, bipartite=true, nodesize=degree(Unes))
 ```
 
+Note that you can change the range of sizes for the nodes using the
+`nodesizerange` argument (a tuple), as well as the symbol for bipartite networks
+using the `bipartiteshapes` (a tuple too):
+
+```@example default
+Unes = web_of_life("M_SD_033")
+I = initial(BipartiteInitialLayout, Unes)
+position!(NestedBipartiteLayout(0.4), I, Unes)
+plot(I, Unes, aspectratio=1)
+scatter!(I, Unes, bipartite=true, nodesize=degree(Unes), nodesizerange=(1.0, 7.0), bipartiteshapes=(:square, :circle))
+```
+
+For quantitative networks, the `plot` method has a `linewidthrange` argument
+that is, similarly, a tuple with the lowest and highest widths allowed.
+
 ## Node annotations
 
 ```@example default

src/circular.jl

@@ -68,6 +68,6 @@ function position!(
     for n in species(N)
         i = 2 * R[Θ[n]] * π / S
         x, y = LA.radius * cos(i), LA.radius * sin(i)
-        L[n] = NodePosition(x, y)
+        L[n] = NodePosition(x=x, y=y)
     end
 end

src/forcedirected.jl

@@ -1,3 +1,15 @@
+"""
+    _force(dist, coeff, expdist, expcoeff)
+
+Return the force for two objects at a distance d with a movemement coefficient c
+and exponents a and b, so that 𝒻 = dᵃ×cᵇ. This works for both attraction and
+repulsion, and the different families of FD layouts are determined by the values
+of a and b.
+"""
+function _force(dist::Float64, coeff::Float64, expdist::Float64, expcoeff::Float64)::Float64
+    return (dist^expdist)*(coeff^expcoeff)
+end
+
 """
     ForceDirectedLayout
 
@@ -48,7 +60,7 @@ end
 
 TODO
 """
-ForceDirectedLayout(ka::Float64, kr::Float64; gravity::Float64=0.75) = ForceDirectedLayout((true,true), (ka,kr), (2.0, -1.0, -1.0, 2.0), gravity, 0.0, true)
+ForceDirectedLayout(ka::Float64, kr::Float64; gravity::Float64=0.75) = ForceDirectedLayout((true, true), (ka, kr), (2.0, -1.0, -1.0, 2.0), gravity, 0.0, true)
 
 """
     FruchtermanRheingold(k::Float64; gravity::Float64=0.75)
@@ -74,7 +86,7 @@ ForceAtlas2(k::Float64; gravity::Float64=0.75) = ForceDirectedLayout((true, true
 In the spring electric layout, attraction is proportional to distance, and
 repulsion to the inverse of the distance squared.
 """
-SpringElectric(k::Float64; gravity::Float64=0.75) = ForceDirectedLayout((true,true), (k, k), (1.0, 1.0, -2.0, 1.0), gravity, 0.0, true)
+SpringElectric(k::Float64; gravity::Float64=0.75) = ForceDirectedLayout((true, true), (k, k), (1.0, 1.0, -2.0, 1.0), gravity, 0.0, true)
 
 """
 Stops the movement of a node position.
@@ -87,36 +99,48 @@ end
 """
 Repel two nodes
 """
-function repel!(LA::T, n1::NodePosition, n2::NodePosition, fr) where {T <: ForceDirectedLayout}
+function repel!(LA::T, n1::NodePosition, n2::NodePosition) where {T<:ForceDirectedLayout}
+    # Distance between the points
     δx = n1.x - n2.x
     δy = n1.y - n2.y
-    Δ = sqrt(δx^2.0+δy^2.0)
-    Δ = Δ == 0.0 ? 0.0001 : Δ
+    Δ = max(1e-4, sqrt(δx^2.0 + δy^2.0))
+    # Effect of degree
+    degree_effect = LA.degree ? max(n1.degree, 1.0) * (max(n2.degree, 1.0)) : 1.0
+    # Raw movement
+    𝒻 = EcologicalNetworksPlots._force(Δ, LA.k[2], LA.exponents[3:4]...)
+    # Calculate the movement
+    movement = (degree_effect * 𝒻) / Δ
+    movement_on_x = δx * movement
+    movement_on_y = δy * movement
+    # Apply the movement
     if LA.move[1]
-        n1.vx += δx/Δ*fr(Δ)
-        n2.vx -= δx/Δ*fr(Δ)
+        n1.vx += movement_on_x
+        n2.vx -= movement_on_x
     end
     if LA.move[2]
-        n1.vy += δy/Δ*fr(Δ)
-        n2.vy -= δy/Δ*fr(Δ)
+        n1.vy += movement_on_y
+        n2.vy -= movement_on_y
     end
 end
 
 """
 Attract two connected nodes
 """
-function attract!(LA::T, n1::NodePosition, n2::NodePosition, fa) where {T <: ForceDirectedLayout}
+function attract!(LA::T, n1::NodePosition, n2::NodePosition, w; gravity=false) where {T<:ForceDirectedLayout}
     δx = n1.x - n2.x
     δy = n1.y - n2.y
-    Δ = sqrt(δx^2.0+δy^2.0)
+    Δ = sqrt(δx^2.0 + δy^2.0)
+    # Raw movement
+    𝒻 = EcologicalNetworksPlots._force(Δ, LA.k[1], LA.exponents[1:2]...)
     if !iszero(Δ)
+        μ = gravity ? ((LA.gravity * 𝒻) / Δ) : ((w^LA.δ * 𝒻) / Δ)
         if LA.move[1]
-            n1.vx -= δx/Δ*fa(Δ)
-            n2.vx += δx/Δ*fa(Δ)
+            n1.vx -= δx * μ
+            n2.vx += δx * μ
         end
         if LA.move[2]
-            n1.vy -= δy/Δ*fa(Δ)
-            n2.vy += δy/Δ*fa(Δ)
+            n1.vy -= δy * μ
+            n2.vy += δy * μ
         end
     end
 end
@@ -125,10 +149,10 @@ end
 Update the position of a node
 """
 function update!(n::NodePosition)
-    Δ = sqrt(n.vx^2.0+n.vy^2.0)
+    Δ = sqrt(n.vx^2.0 + n.vy^2.0)
     if !iszero(Δ)
-        n.x += n.vx/Δ*min(Δ, 0.01)
-        n.y += n.vy/Δ*min(Δ, 0.01)
+        n.x += n.vx / Δ * min(Δ, 0.05)
+        n.y += n.vy / Δ * min(Δ, 0.05)
     end
     stop!(n)
 end
@@ -148,39 +172,28 @@ With the maximal displacement set to 0.01, we have found that k ≈ 100 gives
 acceptable results. This will depend on the complexity of the network, and its
 connectance, as well as the degree and edge strengths distributions.
 """
-function position!(LA::ForceDirectedLayout, L::Dict{K,NodePosition}, N::T) where {T <: EcologicalNetworks.AbstractEcologicalNetwork} where {K}
-
-    degdistr = degree(N)
+function position!(LA::ForceDirectedLayout, L::Dict{K,NodePosition}, N::T) where {T<:EcologicalNetworks.AbstractEcologicalNetwork} where {K}
 
-    # Exponents and forces - the attraction and repulsion functions are
-    # (Δᵃ)×(kₐᵇ) and (Δᶜ)×(kᵣᵈ)
-    a,b,c,d = LA.exponents
-    ka, kr = LA.k
-    fa(x) = (x^a)*(ka^b)
-    fr(x) = (x^c)*(kr^d)
-
-    plotcenter = NodePosition(0.0, 0.0, 0.0, 0.0)
+    # Center point
+    plotcenter = NodePosition(x=0.0, y=0.0)
 
     for (i, s1) in enumerate(species(N))
-        attract!(LA, L[s1], plotcenter, (x) -> LA.gravity*fa(x))
+        if LA.gravity > 0.0
+            attract!(LA, L[s1], plotcenter, LA.gravity; gravity=true)
+        end
         for (j, s2) in enumerate(species(N))
             if j > i
-                if LA.degree
-                    repel!(LA, L[s1], L[s2], (x) -> (degdistr[s1]+1)*(degdistr[s2]+1)*fr(x))
-                else
-                    repel!(LA, L[s1], L[s2], fr)
-                end
+                repel!(LA, L[s1], L[s2])
             end
         end
     end
-    
+
     for int in interactions(N)
-        # We can do Bool^δ and it returns the Bool, so that's tight
-        attract!(LA, L[int.from], L[int.to], (x) -> N[int.from, int.to]^LA.δ*fa(x))
+        attract!(LA, L[int.from], L[int.to], N[int.from, int.to])
     end
 
     for s in species(N)
         update!(L[s])
     end
-    
-end
+
+end

src/initial_layouts.jl

@@ -5,28 +5,30 @@ Random disposition of nodes in a circle. This is a good starting point for any
 force-directed layout. The circle is scaled so that its radius is twice the
 square root of the network richness, which helps most layouts converge faster.
 """
-function initial(::Type{RandomInitialLayout}, N::T) where {T <: EcologicalNetworks.AbstractEcologicalNetwork}
-  L = Dict([s => NodePosition() for s in species(N)])
-  _adj = 2sqrt(richness(N))
-  for s in species(N)
-    L[s].x *= _adj
-    L[s].y *= _adj
-  end
-  return L
+function initial(::Type{RandomInitialLayout}, N::T) where {T<:EcologicalNetworks.AbstractEcologicalNetwork}
+    k = degree(N)
+    L = Dict([s => NodePosition(x=rand(), y=rand(), degree=float(k[s])) for s in species(N)])
+    _adj = 2sqrt(richness(N))
+    for s in species(N)
+        L[s].x *= _adj
+        L[s].y *= _adj
+    end
+    return L
 end
 
 """
     initial(::Type{BipartiteInitialLayout}, N::T) where {T <: EcologicalNetworks.AbstractBipartiteNetwork}
 
 Random disposition of nodes on two levels for bipartite networks.
 """
-function initial(::Type{BipartiteInitialLayout}, N::T) where {T <: EcologicalNetworks.AbstractBipartiteNetwork}
-  level = NodePosition[]
-  for (i, s) in enumerate(species(N))
-    this_level = s ∈ species(N; dims=1) ? 1.0 : 0.0
-    push!(level, NodePosition(rand(), this_level, 0.0, 0.0))
-  end
-  return Dict(zip(species(N), level))
+function initial(::Type{BipartiteInitialLayout}, N::T) where {T<:EcologicalNetworks.AbstractBipartiteNetwork}
+    level = NodePosition[]
+    k = degree(N)
+    for s in species(N)
+        this_level = s ∈ species(N; dims=1) ? 1.0 : 0.0
+        push!(level, NodePosition(x=rand(), y=this_level, degree=k[s]))
+    end
+    return Dict(zip(species(N), level))
 end
 
 """
@@ -36,13 +38,13 @@ Random disposition of nodes on trophic levels for food webs. Note that the
 continuous trophic level is used, but the layout can be modified afterwards to
 use another measure of trophic rank.
 """
-function initial(::Type{FoodwebInitialLayout}, N::T) where {T <: EcologicalNetworks.AbstractUnipartiteNetwork}
-  L = initial(RandomInitialLayout, N)
-  tl = trophic_level(N)
-  for s in species(N)
-    L[s].y = tl[s]
-  end
-  return L
+function initial(::Type{FoodwebInitialLayout}, N::T) where {T<:EcologicalNetworks.AbstractUnipartiteNetwork}
+    L = initial(RandomInitialLayout, N)
+    tl = trophic_level(N)
+    for s in species(N)
+        L[s].y = tl[s]
+    end
+    return L
 end
 
 """
@@ -51,15 +53,16 @@ end
 Random disposition of nodes on a circle. This is the starting point for
 circle-based layouts.
 """
-function initial(::Type{CircularInitialLayout}, N::T) where {T <: EcologicalNetworks.AbstractEcologicalNetwork}
-  level = NodePosition[]
-  n = richness(N)
-  for (i, s) in enumerate(species(N))
-    θ = 2i * π/n
-    x, y = cos(θ), sin(θ)
-    push!(level, NodePosition(x, y, i))
-  end
-  return Dict(zip(species(N), level))
+function initial(::Type{CircularInitialLayout}, N::T) where {T<:EcologicalNetworks.AbstractEcologicalNetwork}
+    level = NodePosition[]
+    k = degree(N)
+    n = richness(N)
+    for (i, s) in enumerate(species(N))
+        θ = 2i * π / n
+        x, y = cos(θ), sin(θ)
+        push!(level, NodePosition(x=x, y=y, r=i, degree=k[s]))
+    end
+    return Dict(zip(species(N), level))
 end
 
 """
@@ -70,13 +73,13 @@ axis is the omnivory index. Note that the *fractional* trophic level is used,
 but the layout can be modified afterwards to use the continuous levels. See the
 documentation for `UnravelledLayout` to see how.
 """
-function initial(::Type{UnravelledInitialLayout}, N::T) where {T <: EcologicalNetworks.AbstractUnipartiteNetwork}
-  layout = Dict([s => NodePosition() for s in species(N)])
-  tl = fractional_trophic_level(N)
-  oi = omnivory(N)
-  for s in species(N)
-    layout[s].x = float(oi[s])
-    layout[s].y = float(tl[s])
-  end
-  return layout
+function initial(::Type{UnravelledInitialLayout}, N::T) where {T<:EcologicalNetworks.AbstractUnipartiteNetwork}
+    layout = Dict([s => NodePosition() for s in species(N)])
+    tl = fractional_trophic_level(N)
+    oi = omnivory(N)
+    for s in species(N)
+        layout[s].x = float(oi[s])
+        layout[s].y = float(tl[s])
+    end
+    return layout
 end

src/recipes.jl

@@ -10,6 +10,9 @@ end
     nodesize=nothing,
     nodefill=nothing,
     bipartite=false,
+    nodesizerange=(2.0,8.0),
+    linewidthrange=(0.5,3.5),
+    bipartiteshapes = (:dtriangle, :utriangle)
 ) where {T<:AbstractEcologicalNetwork} where {K}
 
     # Node positions
@@ -33,7 +36,7 @@ end
                 linecolor --> :darkgrey
                 if typeof(network) <: QuantitativeNetwork
                     linewidth --> EcologicalNetworksPlots._scale_value(
-                        interaction.strength, int_range, (0.5, 3.5)
+                        interaction.strength, int_range, linewidthrange
                     )
                 end
                 if typeof(network) <: ProbabilisticNetwork
@@ -49,23 +52,19 @@ end
             if nodesize !== nothing
                 nsi_range = (minimum(values(nodesize)), maximum(values(nodesize)))
                 markersize := [
-                    EcologicalNetworksPlots._scale_value(nodesize[s], nsi_range, (2, 8)) for
+                    EcologicalNetworksPlots._scale_value(nodesize[s], nsi_range, nodesizerange) for
                     s in species(network)
                 ]
             end
 
             if nodefill !== nothing
-                nfi_range = (minimum(values(nodefill)), maximum(values(nodefill)))
-                marker_z := [
-                    EcologicalNetworksPlots._scale_value(nodefill[s], nfi_range, (0, 1)) for
-                    s in species(network)
-                ]
+                marker_z := [nodefill[s] for s in species(network)]
             end
 
             if bipartite
                 m_shape = Symbol[]
                 for (i, s) in enumerate(species(network))
-                    this_mshape = s ∈ species(network; dims=1) ? :dtriangle : :circle
+                    this_mshape = s ∈ species(network; dims=1) ? bipartiteshapes[1] : bipartiteshapes[2]
                     push!(m_shape, this_mshape)
                 end
                 marker := m_shape