Add plotting examples and additional documentation

ragibson · ragibson · commit a26a8ea52559 · 2021-03-13T15:45:08.000-05:00
diff --git a/TODO.md b/TODO.md
@@ -1,6 +1,5 @@
 ## TODO items
 
-- [ ] Add plotting information (with examples) to documentation
 - [ ] Add forthcoming citations/acknowledgements to README
 
 ## Potential Future Work
diff --git a/docs/example_figures/figure_plot_2d_domains_with_estimates.png b/docs/example_figures/figure_plot_2d_domains_with_estimates.png
diff --git a/docs/example_figures/figure_plot_2d_domains_with_num_communities.png b/docs/example_figures/figure_plot_2d_domains_with_num_communities.png
diff --git a/docs/example_figures/figure_plot_estimates.png b/docs/example_figures/figure_plot_estimates.png
diff --git a/docs/example_figures/figure_plot_multiplex_community.png b/docs/example_figures/figure_plot_multiplex_community.png
diff --git a/docs/index.rst b/docs/index.rst
@@ -18,6 +18,7 @@ communities when desired.
    source/basic_multilayer_example
    source/usage
    source/multilayer_usage
+   source/plotting_examples
    :maxdepth: 2
    :caption: Contents:
 
diff --git a/docs/source/multilayer_usage.rst b/docs/source/multilayer_usage.rst
@@ -67,4 +67,4 @@ These functions provide utilities related to the parameter estimation of `Newman
 modularitypruning.plotting
 --------------------------
 
-Documentation and examples forthcoming.
+See :doc:`plotting_examples`.
diff --git a/docs/source/plotting_examples.rst b/docs/source/plotting_examples.rst
@@ -0,0 +1,147 @@
+Plotting Examples
+=================
+
+We have included a few functions to generate plots of CHAMP domains and resolution parameter estimates in the
+singlelayer and multilayer cases. We encourage you to make your own application-specific plotting code as well (the
+source code in modularitypruning.plotting is very straightforward and may provide inspiration).
+
+Plotting Single-Layer CHAMP Domains and Resolution Parameter Estimates
+----------------------------------------------------------------------
+
+.. autofunction:: modularitypruning.plotting.plot_estimates
+
+An example on the Karate Club network is as follows.
+
+.. code-block:: python
+
+  import igraph as ig
+  import matplotlib.pyplot as plt
+  from modularitypruning.champ_utilities import CHAMP_2D
+  from modularitypruning.parameter_estimation_utilities import ranges_to_gamma_estimates
+  from modularitypruning.louvain_utilities import repeated_parallel_louvain_from_gammas
+  from modularitypruning.plotting import plot_estimates
+  import numpy as np
+
+  # get Karate Club graph in igraph
+  G = ig.Graph.Famous("Zachary")
+
+  # run louvain 100K times on this graph from gamma=0 to gamma=2 (takes ~2-3 seconds)
+  partitions = repeated_parallel_louvain_from_gammas(G, np.linspace(0, 2, 10 ** 5))
+
+  # run CHAMP to obtain the dominant partitions along with their regions of optimality
+  ranges = CHAMP_2D(G, partitions, gamma_0=0.0, gamma_f=2.0)
+
+  # append gamma estimate for each dominant partition onto the CHAMP domains
+  gamma_estimates = ranges_to_gamma_estimates(G, ranges)
+
+  # plot gamma estimates and domains of optimality
+  plt.rc('text', usetex=True)
+  plt.rc('font', family='serif')
+  plot_estimates(gamma_estimates)
+  plt.title(r"Karate Club CHAMP Domains of Optimality and $\gamma$ Estimates", fontsize=14)
+  plt.xlabel(r"$\gamma$", fontsize=14)
+  plt.ylabel("Number of communities", fontsize=14)
+  plt.show()
+
+which generates
+
+.. image:: ../example_figures/figure_plot_estimates.png
+
+Plotting Multi-Layer CHAMP domains and Resolution Parameter Estimates
+---------------------------------------------------------------------
+
+.. autofunction:: modularitypruning.plotting.plot_2d_domains_with_estimates
+
+An example on a realization of multilayer synthetic network in :doc:`../source/basic_multilayer_example` is as follows.
+
+.. code-block:: python
+
+  from modularitypruning.champ_utilities import CHAMP_3D
+  from modularitypruning.louvain_utilities import (repeated_parallel_louvain_from_gammas_omegas,
+                                                   check_multilayer_louvain_capabilities)
+  from modularitypruning.parameter_estimation_utilities import domains_to_gamma_omega_estimates
+  from modularitypruning.plotting import plot_2d_domains_with_estimates
+  import matplotlib.pyplot as plt
+  import numpy as np
+
+  # fail if version of louvain does not support multilayer optimization
+  check_multilayer_louvain_capabilities()
+
+  # run louvain on a uniform grid (10K samples) of gamma and omega (takes ~3 seconds)
+  gamma_range = (0.5, 1.5)
+  omega_range = (0, 2)
+  parts = repeated_parallel_louvain_from_gammas_omegas(G_intralayer, G_interlayer, layer_vec,
+                                                       gammas=np.linspace(*gamma_range, 100),
+                                                       omegas=np.linspace(*omega_range, 100))
+
+  # run CHAMP to obtain the dominant partitions along with their regions of optimality
+  domains = CHAMP_3D(G_intralayer, G_interlayer, layer_vec, parts,
+                     gamma_0=gamma_range[0], gamma_f=gamma_range[1],
+                     omega_0=omega_range[0], omega_f=omega_range[1])
+
+  # append resolution parameter estimates for each dominant partition onto the CHAMP domains
+  domains_with_estimates = domains_to_gamma_omega_estimates(G_intralayer, G_interlayer, layer_vec,
+                                                            domains, model='temporal')
+
+  # plot resolution parameter estimates and domains of optimality
+  plt.rc('text', usetex=True)
+  plt.rc('font', family='serif')
+  plot_2d_domains_with_estimates(domains_with_estimates, xlim=omega_range, ylim=gamma_range)
+  plt.title(r"CHAMP Domains and ($\omega$, $\gamma$) Estimates", fontsize=16)
+  plt.xlabel(r"$\omega$", fontsize=20)
+  plt.ylabel(r"$\gamma$", fontsize=20)
+  plt.gca().tick_params(axis='both', labelsize=12)
+  plt.tight_layout()
+  plt.show()
+
+which generates
+
+.. image:: ../example_figures/figure_plot_2d_domains_with_estimates.png
+
+There are other similar functions
+  * :meth:`~modularitypruning.plotting.plot_2d_domains` essentially plots the above without arrows showing resolution
+    parameter estimates
+  * :meth:`~modularitypruning.plotting.plot_2d_domains_with_ami` plots CHAMP domains, colored by the AMI between each
+    partition and a ground truth
+  * :meth:`~modularitypruning.plotting.plot_2d_domains_with_num_communities` plots CHAMP domains, colored by number
+    of communities
+
+For example, the last function generates
+
+.. image:: ../example_figures/figure_plot_2d_domains_with_num_communities.png
+
+Plotting Multiplex Community Structure of a Single Partition
+------------------------------------------------------------
+
+.. autofunction:: modularitypruning.plotting.plot_multiplex_community
+
+An example from the Lazega Law Firm network is as follows
+
+.. code-block:: python
+
+  from modularitypruning.plotting import plot_multiplex_community
+  import matplotlib.pyplot as plt
+  import numpy as np
+
+  num_layers = 3
+  layer_vec = [i // 71 for i in range(num_layers * 71)]
+  membership = [1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1,
+                1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
+                2, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1,
+                1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
+                2, 2, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0,
+                0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,
+                2, 2, 2]
+
+  plt.rc('text', usetex=True)
+  plt.rc('font', family='serif')
+  ax = plot_multiplex_community(np.array(membership), np.array(layer_vec))
+  ax.set_xticks(np.linspace(0, num_layers, 2 * num_layers + 1))
+  ax.set_xticklabels(["", "Advice", "", "Coworker", "", "Friend", ""], fontsize=14)
+  plt.title(f"Multiplex Communities", fontsize=14)
+  plt.ylabel("Node ID", fontsize=14)
+  plt.show()
+
+which generates
+
+.. image:: ../example_figures/figure_plot_multiplex_community.png
diff --git a/docs/source/usage.rst b/docs/source/usage.rst
@@ -64,4 +64,4 @@ These functions provide utilities related to the parameter estimation of `Newman
 modularitypruning.plotting
 --------------------------
 
-Documentation and examples forthcoming.
+See :doc:`plotting_examples`.
diff --git a/setup.py b/setup.py
@@ -9,7 +9,7 @@
 
 setup(
     name='modularitypruning',
-    version='1.1.1',
+    version='1.1.2',
     package_dir={'modularitypruning': 'utilities'},
     packages=['modularitypruning'],
     url='https://github.com/ragibson/ModularityPruning',
diff --git a/utilities/parameter_estimation_utilities.py b/utilities/parameter_estimation_utilities.py
@@ -522,6 +522,11 @@ def prune_to_multilayer_stable_partitions(G_intralayer, G_interlayer, layer_vec,
 
     See **[CITATION FORTHCOMING]** for more details.
 
+    NOTE: This method truncates omega estimates to ``omega_end - 1e-3`` in order to properly identify stable partitions
+    with infinite interlayer coupling estimates (e.g. when all membership labels persist across layers). If
+    ``omega_end`` is set too low, such partitions may be incorrectly identified as stable. As such, you should be
+    somewhat wary of the returned partitions with zero community structure differences across layers.
+
     :param G_intralayer: intralayer graph of interest
     :type G_intralayer: igraph.Graph
     :param G_interlayer: interlayer graph of interest
diff --git a/utilities/plotting.py b/utilities/plotting.py
@@ -15,7 +15,11 @@ def plot_adjacency(adj):
 
 
 def plot_estimates(gamma_estimates):
-    """Plot partition dominance ranges with gamma estimates."""
+    """Plot partition dominance ranges with gamma estimates.
+
+    :param gamma_estimates: gamma estimates as returned from
+        :meth:`~modularitypruning.plotting.ranges_to_gamma_estimates`
+    """
 
     ax = plt.gca()
     num_parts = len(gamma_estimates)
@@ -68,7 +72,12 @@ def plot_estimates(gamma_estimates):
 def plot_2d_domains(domains, xlim, ylim, flip_axes=False, use_current_axes=False):
     """Plot partition dominance ranges in the (gamma, omega) plane, using the domains from CHAMP_3D.
 
-    Limits output to xlim and ylim dimensions. Note that the plotting here has x=gamma and y=omega."""
+    Limits output to xlim and ylim dimensions. Note that the plotting here has x=gamma and y=omega.
+
+    :param domains: CHAMP domains as returned from :meth:`~modularitypruning.champ_utilities.CHAMP_3D`.
+    :param xlim: plotting x limits
+    :param ylim: plotting y limits
+    """
     if use_current_axes:
         ax = plt.gca()
     else:
@@ -117,7 +126,13 @@ def plot_2d_domains_with_estimates(domains_with_estimates, xlim, ylim, plot_esti
     """Plot partition dominance ranges in the (gamma, omega) plane, using the domains from CHAMP_3D with their gamma
     and omega estimates overlaid.
 
-    Limits output to xlim and ylim dimensions. Note that the plotting here has x=omega and y=gamma."""
+    Limits output to xlim and ylim dimensions. Note that the plotting here has x=omega and y=gamma.
+
+    :param domains_with_estimates: CHAMP domains and resolution parameter estimates as returned from
+        :meth:`~modularitypruning.parameter_estimation_utilities.domains_to_gamma_omega_estimates`.
+    :param xlim: plotting x limits
+    :param ylim: plotting y limits
+    """
     assert flip_axes
     fig, ax = plt.subplots()
     patches = []
@@ -134,6 +149,8 @@ def plot_2d_domains_with_estimates(domains_with_estimates, xlim, ylim, plot_esti
 
         if gamma_est is not None and omega_est is not None:
             width = xlim[1] - xlim[0]
+
+            # prune arrows too far outside plotting range to avoid excess visual noise
             if centroid_x < xlim[0] - 0.05 * width or centroid_x > xlim[1] + 0.05 * width:
                 continue
 
@@ -222,7 +239,13 @@ def plot_2d_domains_with_ami(domains_with_estimates, ground_truth, xlim, ylim, f
 
 
 def plot_multiplex_community(membership, layer_vec):
-    """Plot multiplex community membership"""
+    """Plot a visualization of multiplex community membership
+
+    :param membership: partition membership vector
+    :type membership: tuple[int]
+    :param layer_vec: list of each vertex's layer membership
+    :type layer_vec: list[int]
+    """
     T = max(layer_vec) + 1
     N = len(membership) // T
     communities_per_layer = np.zeros((N, T))
