edit

Counting motifs

In this use case, we will count the number of motifs in a food web. Specifically, we will count how many times there is a linear food chain (A→B→C) between three species.

using EcologicalNetwork

# Load a network
N = stony();

# List of motifs
m = unipartitemotifs();
Dict{Symbol,EcologicalNetwork.UnipartiteNetwork} with 13 entries:
  :D4 => EcologicalNetwork.UnipartiteNetwork(Bool[false true false; false false…
  :S1 => EcologicalNetwork.UnipartiteNetwork(Bool[false true false; false false…
  :D1 => EcologicalNetwork.UnipartiteNetwork(Bool[false true true; false false …
  :D6 => EcologicalNetwork.UnipartiteNetwork(Bool[false true true; true false t…
  :D3 => EcologicalNetwork.UnipartiteNetwork(Bool[false false true; false false…
  :D5 => EcologicalNetwork.UnipartiteNetwork(Bool[false true false; false false…
  :S4 => EcologicalNetwork.UnipartiteNetwork(Bool[false true false; false false…
  :S3 => EcologicalNetwork.UnipartiteNetwork(Bool[false true false; false false…
  :D8 => EcologicalNetwork.UnipartiteNetwork(Bool[false true true; true false f…
  :S2 => EcologicalNetwork.UnipartiteNetwork(Bool[false true true; false false …
  :D2 => EcologicalNetwork.UnipartiteNetwork(Bool[false true true; false false …
  :D7 => EcologicalNetwork.UnipartiteNetwork(Bool[false true true; true false f…
  :S5 => EcologicalNetwork.UnipartiteNetwork(Bool[false true true; false false …

The m object has 13 different motifs, named as in Stouffer et al. (2007). The function unipartitemotifs will generate them when needed.

The function to count motifs is called motif, and returns a count: how many triplets of species are in a given conformation. For example:

s1 = motif(N, m[:S1])
1035.0

We may be interested in knowing whether this motif is over or under-represented in the empirical network, compared to a random expectation. To determine this, we will shuffle interactions around in a way that preserves the number of interactionsand the degree distribution of all species, using swaps. We will create 100 replicated networks to test.

permutations = swaps(N, 100, constraint=:degree)

ms1 = map(x -> motif(x, m[:S1]), permutations)
100-element Array{Float64,1}:
 522.0
 569.0
 636.0
 578.0
 556.0
 554.0
 576.0
 573.0
 597.0
 564.0
   ⋮
 583.0
 570.0
 606.0
 516.0
 528.0
 599.0
 575.0
 629.0
 538.0