helix_graph_algorithms/algorithms/
leiden.rs

1use std::cmp::Ordering;
2use std::collections::BTreeMap;
3use std::num::NonZeroUsize;
4
5use serde::{Deserialize, Serialize};
6
7use super::louvain::{
8    assignment_for_original_nodes, canonical_communities, compress_communities, modularity,
9    next_random, shuffle, validate_community_graph, Community, LevelGraph,
10};
11use crate::{Graph, GraphError, PositiveFiniteF64};
12
13/// Weighted Leiden options matching Graphify's graspologic defaults while
14/// making convergence work explicitly bounded.
15#[derive(Debug, Clone, Copy, PartialEq, Serialize, Deserialize)]
16pub struct LeidenOptions {
17    /// Positive finite modularity resolution.
18    pub resolution: PositiveFiniteF64,
19    /// Positive finite stochastic refinement temperature.
20    pub randomness: PositiveFiniteF64,
21    /// Deterministic base seed.
22    pub seed: u64,
23    /// Independent starts; the highest-quality canonical partition wins.
24    pub trials: NonZeroUsize,
25    /// Maximum local-moving sweeps at each level.
26    pub max_iterations: NonZeroUsize,
27    /// Maximum refinement and aggregation levels.
28    pub max_levels: NonZeroUsize,
29}
30
31impl Default for LeidenOptions {
32    fn default() -> Self {
33        Self {
34            resolution: PositiveFiniteF64::new(1.0).expect("one is positive"),
35            randomness: PositiveFiniteF64::new(0.001).expect("default randomness is positive"),
36            seed: 42,
37            trials: NonZeroUsize::new(1).expect("one is non-zero"),
38            max_iterations: NonZeroUsize::new(100).expect("100 is non-zero"),
39            max_levels: NonZeroUsize::new(10).expect("10 is non-zero"),
40        }
41    }
42}
43
44/// Canonical weighted Leiden partition and diagnostics.
45#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
46pub struct LeidenResult {
47    /// Canonically ordered communities.
48    pub communities: Vec<Community>,
49    /// Final weighted modularity quality.
50    pub modularity: f64,
51    /// Number of refinement and aggregation levels in the winning trial.
52    pub levels: usize,
53    /// Zero-based winning trial index.
54    pub winning_trial: usize,
55}
56
57struct TrialResult {
58    communities: Vec<Community>,
59    modularity: f64,
60    levels: usize,
61}
62
63impl Graph {
64    /// Compute weighted Leiden communities with connected refinement.
65    pub fn leiden(&self, options: LeidenOptions) -> Result<LeidenResult, GraphError> {
66        validate_community_graph(self, "Leiden")?;
67        if self.node_count() == 0 {
68            return Ok(LeidenResult {
69                communities: Vec::new(),
70                modularity: 0.0,
71                levels: 0,
72                winning_trial: 0,
73            });
74        }
75
76        let mut best: Option<(usize, TrialResult)> = None;
77        for trial in 0..options.trials.get() {
78            let seed = options.seed.wrapping_add(
79                u64::try_from(trial)
80                    .expect("usize trial fits in u64")
81                    .wrapping_mul(0x9e37_79b9_7f4a_7c15),
82            );
83            let candidate = leiden_trial(self, options, seed);
84            let replace = match &best {
85                None => true,
86                Some((_, current)) => match candidate.modularity.total_cmp(&current.modularity) {
87                    Ordering::Greater => true,
88                    Ordering::Equal => candidate.communities < current.communities,
89                    Ordering::Less => false,
90                },
91            };
92            if replace {
93                best = Some((trial, candidate));
94            }
95        }
96        let (winning_trial, best) = best.expect("non-zero trials produce a result");
97        Ok(LeidenResult {
98            communities: best.communities,
99            modularity: best.modularity,
100            levels: best.levels,
101            winning_trial,
102        })
103    }
104}
105
106fn leiden_trial(graph: &Graph, options: LeidenOptions, mut seed: u64) -> TrialResult {
107    let original = LevelGraph::from_original(graph);
108    let original_edges = original.edges.clone();
109    if original.total_weight == 0.0 {
110        return TrialResult {
111            communities: canonical_communities(graph, &(0..graph.node_count()).collect::<Vec<_>>()),
112            modularity: 0.0,
113            levels: 0,
114        };
115    }
116
117    let mut level = original;
118    let mut initial = (0..level.node_members.len()).collect::<Vec<_>>();
119    let mut final_assignment = (0..graph.node_count()).collect::<Vec<_>>();
120    let mut final_modularity = 0.0;
121    let mut levels = 0;
122    for _ in 0..options.max_levels.get() {
123        let (moved, moved_any) = level.local_move(
124            &initial,
125            options.resolution.get(),
126            &mut seed,
127            options.max_iterations.get(),
128        );
129        let refined = refine_connected(
130            &level,
131            &moved,
132            options.resolution.get(),
133            options.randomness.get(),
134            &mut seed,
135        );
136        levels += 1;
137        final_assignment =
138            assignment_for_original_nodes(graph.node_count(), &level.node_members, &refined);
139        final_modularity = modularity(
140            graph.node_count(),
141            &original_edges,
142            &final_assignment,
143            options.resolution.get(),
144        );
145
146        let refined_count = refined.iter().max().map_or(0, |maximum| maximum + 1);
147        if !moved_any || refined_count == level.node_members.len() {
148            break;
149        }
150        initial = aggregate_initial_partition(&moved, &refined);
151        level = level.aggregate(&refined);
152    }
153
154    TrialResult {
155        communities: canonical_communities(graph, &final_assignment),
156        modularity: final_modularity,
157        levels,
158    }
159}
160
161/// Refine each locally moved community from singleton, connected pieces.
162/// Nodes can only join a same-parent community reached by a positive-weight
163/// edge, so every refined community is connected by construction.
164fn refine_connected(
165    graph: &LevelGraph,
166    parent: &[usize],
167    resolution: f64,
168    randomness: f64,
169    seed: &mut u64,
170) -> Vec<usize> {
171    let node_count = graph.node_members.len();
172    if graph.total_weight == 0.0 {
173        return (0..node_count).collect();
174    }
175    let mut refined = (0..node_count).collect::<Vec<_>>();
176    let mut community_degrees = graph.degrees.clone();
177    let mut community_sizes = vec![1_usize; node_count];
178    let mut order = (0..node_count).collect::<Vec<_>>();
179    shuffle(&mut order, seed);
180
181    for node in order {
182        let old_community = refined[node];
183        if community_sizes[old_community] != 1 {
184            continue;
185        }
186        let mut edge_weight_by_community = BTreeMap::<usize, f64>::new();
187        for (neighbor, weight) in &graph.adjacency[node] {
188            if *weight > 0.0 && parent[*neighbor] == parent[node] {
189                *edge_weight_by_community
190                    .entry(refined[*neighbor])
191                    .or_default() += *weight;
192            }
193        }
194        if edge_weight_by_community.is_empty() {
195            continue;
196        }
197
198        let node_degree = graph.degrees[node];
199        let denominator = 2.0 * graph.total_weight * graph.total_weight;
200        let mut candidates = vec![(old_community, 1.0_f64)];
201        for (community, edge_weight) in edge_weight_by_community {
202            if community == old_community {
203                continue;
204            }
205            let gain = edge_weight / graph.total_weight
206                - resolution * community_degrees[community] * node_degree / denominator;
207            if gain >= 0.0 {
208                candidates.push((community, (gain / randomness).exp()));
209            }
210        }
211        let total = candidates.iter().map(|(_, weight)| weight).sum::<f64>();
212        let chosen = if total.is_finite() {
213            let random = random_unit(seed) * total;
214            let mut cumulative = 0.0;
215            candidates
216                .iter()
217                .find_map(|(community, weight)| {
218                    cumulative += weight;
219                    (random < cumulative).then_some(*community)
220                })
221                .unwrap_or(old_community)
222        } else {
223            candidates
224                .iter()
225                .max_by(|left, right| {
226                    left.1
227                        .total_cmp(&right.1)
228                        .then_with(|| right.0.cmp(&left.0))
229                })
230                .map_or(old_community, |(community, _)| *community)
231        };
232        if chosen != old_community {
233            refined[node] = chosen;
234            community_degrees[old_community] -= node_degree;
235            community_degrees[chosen] += node_degree;
236            community_sizes[old_community] -= 1;
237            community_sizes[chosen] += 1;
238        }
239    }
240    compress_communities(&refined)
241}
242
243fn aggregate_initial_partition(parent: &[usize], refined: &[usize]) -> Vec<usize> {
244    let refined_count = refined.iter().max().map_or(0, |maximum| maximum + 1);
245    let mut parents = vec![usize::MAX; refined_count];
246    for (node, refined_community) in refined.iter().enumerate() {
247        let parent_community = &mut parents[*refined_community];
248        if *parent_community == usize::MAX {
249            *parent_community = parent[node];
250        } else {
251            assert_eq!(
252                *parent_community, parent[node],
253                "refinement communities remain inside their parent"
254            );
255        }
256    }
257    compress_communities(&parents)
258}
259
260fn random_unit(seed: &mut u64) -> f64 {
261    const DENOMINATOR: f64 = (1_u64 << 53) as f64;
262    ((next_random(seed) >> 11) as f64) / DENOMINATOR
263}
264
265#[cfg(test)]
266mod tests {
267    use approx::assert_abs_diff_eq;
268
269    use super::*;
270    use crate::{Edge, GraphKind, Node};
271
272    fn bridged_triangles(kind: GraphKind) -> Graph {
273        Graph::new(
274            kind,
275            ["a", "b", "c", "d", "e", "f", "isolated"]
276                .into_iter()
277                .map(Node::new),
278            [
279                Edge::new("ab", "a", "b").with_weight(2.0),
280                Edge::new("bc", "b", "c").with_weight(2.0),
281                Edge::new("ca", "c", "a").with_weight(2.0),
282                Edge::new("de", "d", "e").with_weight(2.0),
283                Edge::new("ef", "e", "f").with_weight(2.0),
284                Edge::new("fd", "f", "d").with_weight(2.0),
285                Edge::new("cd", "c", "d").with_weight(0.1),
286            ],
287        )
288        .unwrap()
289    }
290
291    #[test]
292    fn weighted_leiden_is_seeded_and_keeps_isolates() {
293        let graph = bridged_triangles(GraphKind::Graph);
294        let first = graph.leiden(LeidenOptions::default()).unwrap();
295        let second = graph.leiden(LeidenOptions::default()).unwrap();
296        assert_eq!(first, second);
297        let options_json = serde_json::to_vec(&LeidenOptions::default()).unwrap();
298        assert_eq!(
299            serde_json::from_slice::<LeidenOptions>(&options_json).unwrap(),
300            LeidenOptions::default()
301        );
302        let result_json = serde_json::to_vec(&first).unwrap();
303        let decoded = serde_json::from_slice::<LeidenResult>(&result_json).unwrap();
304        assert_eq!(decoded.communities, first.communities);
305        assert_eq!(decoded.levels, first.levels);
306        assert_eq!(decoded.winning_trial, first.winning_trial);
307        assert_abs_diff_eq!(decoded.modularity, first.modularity, epsilon = f64::EPSILON);
308        assert_eq!(first.communities.len(), 3);
309        assert_eq!(first.communities[0].node_ids, ["a", "b", "c"]);
310        assert_eq!(first.communities[1].node_ids, ["d", "e", "f"]);
311        assert_eq!(first.communities[2].node_ids, ["isolated"]);
312        assert!(first.modularity > 0.0);
313    }
314
315    #[test]
316    fn leiden_supports_parallel_edges_self_loops_trials_and_caps() {
317        let graph = Graph::new(
318            GraphKind::MultiGraph,
319            [Node::new("a"), Node::new("b"), Node::new("c")],
320            [
321                Edge::new("ab1", "a", "b").with_weight(2.0),
322                Edge::new("ab2", "a", "b").with_weight(3.0),
323                Edge::new("bc", "b", "c").with_weight(0.1),
324                Edge::new("aa", "a", "a").with_weight(0.5),
325            ],
326        )
327        .unwrap();
328        let result = graph
329            .leiden(LeidenOptions {
330                trials: NonZeroUsize::new(3).unwrap(),
331                max_iterations: NonZeroUsize::new(1).unwrap(),
332                max_levels: NonZeroUsize::new(1).unwrap(),
333                ..LeidenOptions::default()
334            })
335            .unwrap();
336        assert!(result.winning_trial < 3);
337        assert_eq!(result.levels, 1);
338        assert!(result.communities.iter().any(|community| {
339            community.node_ids.contains(&"a".into()) && community.node_ids.contains(&"b".into())
340        }));
341    }
342
343    #[test]
344    fn refinement_never_merges_disconnected_components() {
345        let graph = Graph::new(
346            GraphKind::Graph,
347            [
348                Node::new("a"),
349                Node::new("b"),
350                Node::new("c"),
351                Node::new("d"),
352            ],
353            [Edge::new("ab", "a", "b"), Edge::new("cd", "c", "d")],
354        )
355        .unwrap();
356        let level = LevelGraph::from_original(&graph);
357        let refined = refine_connected(&level, &[0, 0, 0, 0], 1.0, 0.001, &mut 42);
358        assert_eq!(refined[0], refined[1]);
359        assert_eq!(refined[2], refined[3]);
360        assert_ne!(refined[0], refined[2]);
361    }
362
363    #[test]
364    fn leiden_rejects_directed_graphs_and_handles_empty_or_zero_weight_graphs() {
365        assert!(matches!(
366            bridged_triangles(GraphKind::DiGraph).leiden(LeidenOptions::default()),
367            Err(GraphError::InvalidOption(_))
368        ));
369        let empty = Graph::new(GraphKind::Graph, [], []).unwrap();
370        assert_eq!(empty.leiden(LeidenOptions::default()).unwrap().levels, 0);
371        let zero = Graph::new(
372            GraphKind::Graph,
373            [Node::new("a"), Node::new("b")],
374            [Edge::new("ab", "a", "b").with_weight(0.0)],
375        )
376        .unwrap()
377        .leiden(LeidenOptions::default())
378        .unwrap();
379        assert_eq!(zero.communities.len(), 2);
380        assert_abs_diff_eq!(zero.modularity, 0.0);
381    }
382
383    #[test]
384    fn randomness_resolution_and_canonical_trial_ties_are_explicit() {
385        assert!(PositiveFiniteF64::new(0.0).is_err());
386        let graph = bridged_triangles(GraphKind::Graph);
387        let high_resolution = graph
388            .leiden(LeidenOptions {
389                resolution: PositiveFiniteF64::new(2.0).unwrap(),
390                randomness: PositiveFiniteF64::new(0.01).unwrap(),
391                seed: 7,
392                trials: NonZeroUsize::new(2).unwrap(),
393                ..LeidenOptions::default()
394            })
395            .unwrap();
396        assert!(high_resolution.winning_trial < 2);
397        assert!(high_resolution.modularity.is_finite());
398    }
399}
helix_graph_algorithms/algorithms/leiden.rs

helix_graph_algorithms/algorithms/
leiden.rs
