Refactoring

.JuliaFormatter.toml

@@ -1 +1,4 @@
-style = "sciml"
+indent = 4
+margin = 92
+normalize_line_endings = "unix"
+style = "sciml"

.github/dependabot.yml

@@ -0,0 +1,10 @@
+version: 2
+updates:
+  - package-ecosystem: "github-actions"
+    directory: "/"
+    schedule:
+      interval: "monthly"
+    open-pull-requests-limit: 99
+    labels:
+      - "dependencies"
+      - "github-actions"

.github/workflows/CI.yml

@@ -18,7 +18,7 @@ jobs:
       fail-fast: false
       matrix:
         version:
-          - '1.8'
+          - '1.11'
           # - 'nightly'
         os:
           - ubuntu-latest
@@ -35,7 +35,7 @@ jobs:
       - uses: julia-actions/julia-buildpkg@v1
       - uses: julia-actions/julia-runtest@v1
       - uses: julia-actions/julia-processcoverage@v1
-      - uses: codecov/codecov-action@v3
+      - uses: codecov/codecov-action@v5
         with:
           files: lcov.info
           token: ${{ secrets.CODECOV_TOKEN }}

.github/workflows/benchmark_pr.yml

@@ -0,0 +1,78 @@
+name: Benchmark a pull request
+
+on:
+  pull_request_target:
+    branches:
+      - master
+
+permissions:
+  pull-requests: write
+
+jobs:
+    generate_plots:
+        runs-on: ubuntu-latest
+
+        steps:
+            - uses: actions/checkout@v4
+            - uses: julia-actions/setup-julia@v2
+              with:
+                version: "1.11"
+            - uses: julia-actions/cache@v2
+            - name: Extract Package Name from Project.toml
+              id: extract-package-name
+              run: |
+                PACKAGE_NAME=$(grep "^name" Project.toml | sed 's/^name = "\(.*\)"$/\1/')
+                echo "::set-output name=package_name::$PACKAGE_NAME"
+            - name: Build AirspeedVelocity
+              env:
+                JULIA_NUM_THREADS: 2
+              run: |
+                # Lightweight build step, as sometimes the runner runs out of memory:
+                julia -e 'ENV["JULIA_PKG_PRECOMPILE_AUTO"]=0; import Pkg; Pkg.add(;url="https://github.com/MilesCranmer/AirspeedVelocity.jl.git")'
+                julia -e 'ENV["JULIA_PKG_PRECOMPILE_AUTO"]=0; import Pkg; Pkg.build("AirspeedVelocity")'
+            - name: Add ~/.julia/bin to PATH
+              run: |
+                echo "$HOME/.julia/bin" >> $GITHUB_PATH
+            - name: Run benchmarks
+              run: |
+                echo $PATH
+                ls -l ~/.julia/bin
+                mkdir results
+                benchpkg ${{ steps.extract-package-name.outputs.package_name }} --rev="${{github.event.repository.default_branch}},${{github.event.pull_request.head.sha}}" --url=${{ github.event.repository.clone_url }} --bench-on="${{github.event.repository.default_branch}}" --output-dir=results/ --tune
+            - name: Create plots from benchmarks
+              run: |
+                mkdir -p plots
+                benchpkgplot ${{ steps.extract-package-name.outputs.package_name }} --rev="${{github.event.repository.default_branch}},${{github.event.pull_request.head.sha}}" --npart=10 --format=png --input-dir=results/ --output-dir=plots/
+            - name: Upload plot as artifact
+              uses: actions/upload-artifact@v4
+              with:
+                name: plots
+                path: plots
+            - name: Create markdown table from benchmarks
+              run: |
+                benchpkgtable ${{ steps.extract-package-name.outputs.package_name }} --rev="${{github.event.repository.default_branch}},${{github.event.pull_request.head.sha}}" --input-dir=results/ --ratio > table.md
+                echo '### Benchmark Results' > body.md
+                echo '' >> body.md
+                echo '' >> body.md
+                cat table.md >> body.md
+                echo '' >> body.md
+                echo '' >> body.md
+                echo '### Benchmark Plots' >> body.md
+                echo 'A plot of the benchmark results have been uploaded as an artifact to the workflow run for this PR.' >> body.md
+                echo 'Go to "Actions"->"Benchmark a pull request"->[the most recent run]->"Artifacts" (at the bottom).' >> body.md
+
+            - name: Find Comment
+              uses: peter-evans/find-comment@v3
+              id: fcbenchmark
+              with:
+                issue-number: ${{ github.event.pull_request.number }}
+                comment-author: 'github-actions[bot]'
+                body-includes: Benchmark Results
+
+            - name: Comment on PR
+              uses: peter-evans/create-or-update-comment@v4
+              with:
+                comment-id: ${{ steps.fcbenchmark.outputs.comment-id }}
+                issue-number: ${{ github.event.pull_request.number }}
+                body-path: body.md
+                edit-mode: replace

.gitignore

@@ -1,7 +1,34 @@
 .vscode
 *dev/
+.CondaPkg/
 
+
+# Files generated by invoking Julia with --code-coverage
+*.jl.cov
+*.jl.*.cov
+
+# Files generated by invoking Julia with --track-allocation
+*.jl.mem
+
+# System-specific files and directories generated by the BinaryProvider and BinDeps packages
+# They contain absolute paths specific to the host computer, and so should not be committed
+deps/deps.jl
+deps/build.log
+deps/downloads/
+deps/usr/
+deps/src/
+
+# Build artifacts for creating documentation generated by the Documenter package
 docs/build/
+docs/site/
 docs/Manifest.toml
 
-Manifest.toml
+# File generated by Pkg, the package manager, based on a corresponding Project.toml
+# It records a fixed state of all packages used by the project. As such, it should not be
+# committed for packages, but should be committed for applications that require a static
+# environment.
+Manifest.toml
+
+
+## python
+__pycache__/

CondaPkg.toml

@@ -0,0 +1,2 @@
+[deps]
+pot = ""

Project.toml

@@ -1,36 +1,49 @@
 name = "NetworkHistogram"
 uuid = "7806f430-7229-459c-b2e6-df35e8e4eb5d"
+version = "0.6.0"
 authors = ["Charles Dufour", "Jake Grainger"]
-version = "0.5.2"
 
 [deps]
-ArnoldiMethod = "ec485272-7323-5ecc-a04f-4719b315124d"
-Arpack = "7d9fca2a-8960-54d3-9f78-7d1dccf2cb97"
+Accessors = "7d9f7c33-5ae7-4f3b-8dc6-eff91059b697"
+ArgCheck = "dce04be8-c92d-5529-be00-80e4d2c0e197"
 BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
-CodecZstd = "6b39b394-51ab-5f42-8807-6242bab2b4c2"
-HTTP = "cd3eb016-35fb-5094-929b-558a96fad6f3"
-JLD = "4138dd39-2aa7-5051-a626-17a0bb65d9c8"
-Kronecker = "2c470bb0-bcc8-11e8-3dad-c9649493f05e"
+Clustering = "aaaa29a8-35af-508c-8bc3-b662a17a0fe5"
+Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
+Graphons = "e0c12bfd-47d7-434e-afb7-632611640ca5"
+Hungarian = "e91730f6-4275-51fb-a7a0-7064cfbd3b39"
+IntervalSets = "8197267c-284f-5f27-9208-e0e47529a953"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+LogExpFunctions = "2ab3a3ac-af41-5b50-aa03-7779005ae688"
 ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
+Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
+SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
+StatsAPI = "82ae8749-77ed-4fe6-ae5f-f523153014b0"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
-TranscodingStreams = "3bb67fe8-82b1-5028-8e26-92a6c54297fa"
-ValueHistories = "98cad3c8-aec3-5f06-8e41-884608649ab7"
+
+[weakdeps]
+Bootstrap = "e28b5b4c-05e8-5b66-bc03-6f0c0a0a06e0"
+LightMC = "b58f5c6e-c887-41d6-b553-02118416cd5d"
+Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
+PythonCall = "6099a3de-0909-46bc-b1f4-468b9a2dfc0d"
+
+[extensions]
+BootstrapExt = "Bootstrap"
+LightMCExt = "LightMC"
+MakieExt = "Makie"
+PythonOptimalTransport = "PythonCall"
 
 [compat]
-ArnoldiMethod = "0.2.0"
-Arpack = "0.5.4"
-BenchmarkTools = "1.3.2"
-CodecZstd = "0.7.2"
-HTTP = "1.7.4"
-JLD = "0.13.3"
-Kronecker = "0.5"
-ProgressMeter = "1.7.2"
-StatsBase = "0.33.21"
-TranscodingStreams = "0.9.11"
-ValueHistories = "0.5.4"
-julia = "1.8"
+Accessors = "0.1.42"
+ArgCheck = "2.5.0"
+BenchmarkTools = "1.6.3"
+Clustering = "0.15.8"
+Hungarian = "0.7.0"
+IntervalSets = "0.7.11"
+LinearAlgebra = "1.12.0"
+LogExpFunctions = "0.3.29"
+Reexport = "1.2.2"
 
 [extras]
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

README.md

@@ -12,9 +12,31 @@
 [![Stable](https://img.shields.io/badge/docs-stable-blue.svg)](https://sds-epfl.github.io/NetworkHistogram.jl/stable/)
 [![DOI](https://zenodo.org/badge/572018079.svg)](https://zenodo.org/doi/10.5281/zenodo.10212851)
 
-
-Implementation of the network histogram for graphon estimation from the paper [Network histograms and universality of blockmodel approximation](https://doi.org/10.1073/pnas.1400374111) by Sofia C. Olhede and Patrick J. Wolfe.
-
+Implementation of the network histogram for graphon estimation from the paper
+[Network histograms and universality of blockmodel approximation (2014)](https://doi.org/10.1073/pnas.1400374111)
+by Sofia C. Olhede and Patrick J. Wolfe and its extension to decorated graphs
+by Charles Dufour and Sofia C. Olhede
+[Inference for decorated graphs and application to multiplex networks (2024)](https://arxiv.org/abs/2408.12339).
+
+The network histogram is a nonparametric estimator for the generating mechanism
+of an exchangeable random graph (see graphons, decorated graphons and
+probability graphons). We assume our observed graph is
+$A \in \mathcal{K}^{n \times n}$, where $\mathcal{K}$ is a set of edge
+decorations (e.g. $\{0,1\}$ for unweighted graphs, $\mathbb{N}$ for count
+edges, $\mathbb{R}$ for real-valued edges, etc.). Using the Aldous-Hoover
+theorem, we know that $A$ is generated from a graphon
+$W: [0,1]^2 \to \mathcal{P}\left(\mathcal{K}\right)$, where
+$\mathcal{P}\left(\mathcal{K}\right)$ is the set of probability measures on
+$\mathcal{K}$ in the following way:
+
+1. Sample $U_1, \ldots, U_n \sim \text{iid } \text{Uniform}[0,1]$.
+2. For each pair of nodes $i,j$, sample the edge $A_{ij} \sim W(U_i, U_j)$
+   independently.
+
+The network histogram approximates the generating graphon
+$W: [0,1]^2 \to \mathcal{P}\left(\mathcal{K}\right)$ by a piecewise constant
+function, i.e. a stochastic block model with $k$ blocks. For details, see the
+papers mentioned above.
 
 ## Installation
 
@@ -24,7 +46,10 @@ Pkg.add("NetworkHistogram")
 
 ## Usage
 
-We fit the estimator and then extract the estimated graphon matrix and node labels.
+### Basic Usage
+
+We fit the estimator and then extract the estimated graphon matrix and node
+labels.
 
 ```julia
 using NetworkHistogram, LinearAlgebra
@@ -35,7 +60,7 @@ A[diagind(A)] .= 0
 # approximate the graphon with a network histogram
 hist = graphhist(A)
 
-# get the graphist structure
+# get the graphhist structure
 estimate = hist.graphhist
 
 # get the estimated graphon matrix
@@ -45,4 +70,105 @@ sbm_matrix = estimate.θ
 node_labels = estimate.node_labels
 ```
 
-You can control the optimization process by modifying the rules used in the optimization. Check out the docs for more information.
+### Advanced Usage with Custom Parameters
+
+You can control the optimization process by modifying the rules used in the
+optimization:
+
+```julia
+using NetworkHistogram
+
+# Binary network
+A = Symmetric(rand(0:1, 100, 100))
+A[diagind(A)] .= 0
+
+# Initial partition into k groups
+k = 3
+initial_labels = rand(1:k, 100)
+
+# Configure optimization parameters
+params = GreedyParams(
+    50_000,                      # Maximum iterations
+    RandomNodeSwap(),             # How to select nodes to swap
+    Strict(),                     # Only accept improvements
+    PreviousBestValue(5000),     # Stop after 5000 iterations without improvement
+    true                          # Show progress bar
+)
+
+# Fit the network histogram
+result = nethist(A, Bernoulli(0.5), initial_labels, params)
+
+# Extract results
+ll = loglikelihood(result)
+block_params = result.θ
+node_groups = result.node_labels
+```
+
+### Working with Different Edge Types
+
+The package supports various edge types through custom distributions:
+
+```julia
+using NetworkHistogram
+using Distributions  # For standard distributions
+
+# Example 1: Weighted networks with continuous edges
+W = Symmetric(rand(100, 100))
+W[diagind(W)] .= 0
+# You can use any distribution that implements the required interface
+
+# Example 2: Count data (e.g., number of interactions)
+C = Symmetric(rand(Poisson(2), 100, 100))
+C[diagind(C)] .= 0
+# Use appropriate count distribution
+
+# Example 3: Sparse networks with missing edges
+A_sparse = Symmetric(rand([0, 1, missing], 100, 100))
+A_sparse[diagind(A_sparse)] .= 0
+# Missing values are treated as absent edges
+```
+
+### Visualizing Results (with Makie.jl)
+
+```julia
+using NetworkHistogram
+using CairoMakie  # or GLMakie
+
+# Fit model
+A = Symmetric(rand(0:1, 100, 100))
+A[diagind(A)] .= 0
+result = nethist(A, Bernoulli(0.5), rand(1:3, 100), GreedyParams())
+
+# Create heatmap of estimated parameters
+fig = heatmap_params(result, ordering=true, colormap=:viridis)
+save("network_histogram.png", fig)
+```
+
+### Sampling from a Block Model
+
+```julia
+using NetworkHistogram
+
+# Define a 3-block model
+k = 3
+bm = BlockModel(k, Bernoulli(0.5))
+
+# Set custom edge probabilities between blocks
+bm[1, 1] = Bernoulli(0.8)  # High within-group connectivity
+bm[2, 2] = Bernoulli(0.7)
+bm[3, 3] = Bernoulli(0.6)
+bm[1, 2] = Bernoulli(0.1)  # Low between-group connectivity
+bm[1, 3] = Bernoulli(0.05)
+bm[2, 3] = Bernoulli(0.05)
+
+# Sample a network
+n_nodes = 150
+latents, A = sample(bm, n_nodes)
+
+# latents contains the true block assignments
+# A is the sampled adjacency matrix
+```
+
+Check out the
+[documentation](https://sds-epfl.github.io/NetworkHistogram.jl/dev/) for more
+examples and detailed API information.