igraph/linalg/embedding.rs
1//! Spectral embeddings (`igraph_embedding.h`).
2
3use super::{ArpackOptions, check_finite};
4use crate::{
5 error::{Error, Result},
6 ffi::*,
7 igraph_call,
8 matrix::Matrix,
9 vector::Vector,
10};
11use std::ptr;
12
13crate::ffi_enum! {
14 /// Which eigenvalues (singular values for directed graphs) the spectral
15 /// embeddings use (the admissible subset of `igraph_eigen_which_position_t`).
16 pub enum EmbeddingWhich: igraph_eigen_which_position_t {
17 /// Largest magnitude (`IGRAPH_EIGEN_LM`).
18 LargestMagnitude = igraph_eigen_which_position_t_IGRAPH_EIGEN_LM,
19 /// Largest algebraic value (`IGRAPH_EIGEN_LA`); the same as
20 /// `LargestMagnitude` for directed graphs, whose singular values are
21 /// non-negative.
22 LargestAlgebraic = igraph_eigen_which_position_t_IGRAPH_EIGEN_LA,
23 /// Smallest algebraic value (`IGRAPH_EIGEN_SA`).
24 SmallestAlgebraic = igraph_eigen_which_position_t_IGRAPH_EIGEN_SA,
25 }
26}
27
28crate::ffi_enum! {
29 /// The Laplacian used by [`Graph::laplacian_spectral_embedding`](crate::Graph::laplacian_spectral_embedding)
30 /// (`igraph_laplacian_spectral_embedding_type_t`). `D` is the degree
31 /// (strength) matrix, `A` the adjacency matrix.
32 pub enum LaplacianEmbeddingType: igraph_laplacian_spectral_embedding_type_t {
33 /// `D - A`, the combinatorial Laplacian; undirected graphs only.
34 DA = igraph_laplacian_spectral_embedding_type_t_IGRAPH_EMBEDDING_D_A,
35 /// `I - D^-1/2 A D^-1/2`, the symmetric normalized Laplacian;
36 /// undirected graphs only.
37 IDAD = igraph_laplacian_spectral_embedding_type_t_IGRAPH_EMBEDDING_I_DAD,
38 /// `D^-1/2 A D^-1/2`; undirected graphs only.
39 DAD = igraph_laplacian_spectral_embedding_type_t_IGRAPH_EMBEDDING_DAD,
40 /// `O^-1/2 A P^-1/2` with the out- and in-degree matrices `O` and
41 /// `P`; the only choice for directed graphs.
42 OAP = igraph_laplacian_spectral_embedding_type_t_IGRAPH_EMBEDDING_OAP,
43 }
44}
45
46/// A spectral embedding of the vertices of a graph.
47#[derive(Debug, Clone, PartialEq)]
48pub struct SpectralEmbedding {
49 /// The latent positions `X` (or the singular vectors `U` when not
50 /// scaled): one row per vertex, one column per dimension.
51 pub x: Matrix,
52 /// The second half of the latent positions `Y` (or `V`), for directed
53 /// graphs; equal to [`x`](Self::x) for undirected graphs.
54 pub y: Matrix,
55 /// The eigenvalues (undirected) or singular values (directed) used; empty
56 /// for graphs without edges.
57 pub d: Vec<f64>,
58}
59
60fn check_common(graph: &igraph_t, no: usize, weights: Option<&[f64]>) -> Result<()> {
61 let n = graph.vcount();
62 if no == 0 || no > n {
63 return Err(Error::invalid(format!(
64 "cannot compute a {no}-dimensional embedding of {n} vertices"
65 )));
66 }
67 if let Some(w) = weights {
68 if w.len() != graph.ecount() {
69 return Err(Error::invalid(
70 "the weight vector length must match the number of edges",
71 ));
72 }
73 check_finite(w, "the weight vector")?;
74 }
75 Ok(())
76}
77
78/// The normalized Laplacians divide by square roots of the vertex
79/// strengths: a zero (or negative) strength would feed NaN or infinite
80/// values to ARPACK, which then aborts the process.
81fn check_strengths(
82 graph: &igraph_t,
83 weights: Option<&[f64]>,
84 kind: LaplacianEmbeddingType,
85) -> Result<()> {
86 if kind == LaplacianEmbeddingType::DA {
87 return Ok(());
88 }
89 let n = graph.vcount();
90 let (mut out, mut inn) = (vec![0.0; n], vec![0.0; n]);
91 for (e, (a, b)) in graph.edge_list().into_iter().enumerate() {
92 let w = weights.map_or(1.0, |w| w[e]);
93 out[a as usize] += w;
94 inn[b as usize] += w;
95 }
96 let ok = if graph.is_directed() {
97 out.iter().chain(&inn).all(|&s| s > 0.0)
98 } else {
99 out.iter().zip(&inn).all(|(o, i)| o + i > 0.0)
100 };
101 if ok {
102 Ok(())
103 } else {
104 Err(Error::invalid(
105 "normalized Laplacian embeddings need positive vertex strengths (no isolated vertices; for directed graphs, positive in- and out-strengths)",
106 ))
107 }
108}
109
110impl igraph_t {
111 /// Adjacency spectral embedding (`igraph_adjacency_spectral_embedding`).
112 ///
113 /// Computes a `no`-dimensional Euclidean representation of the graph from
114 /// the singular value decomposition of its adjacency matrix,
115 /// `A = U D V'`. For undirected graphs `X = U_no D^(1/2)` (and `Y = X`);
116 /// for directed graphs `X = U_no D^(1/2)` and `Y = V_no D^(1/2)`. If the
117 /// graph is a random dot product graph generated from latent position
118 /// vectors in `R^no`, the embedding estimates these latent positions.
119 ///
120 /// - `weights`: optional edge weights.
121 /// - `which`: which eigenvalues (singular values) to use.
122 /// - `scaled`: return `X`, `Y` if true, `U`, `V` otherwise.
123 /// - `cvec`: added to the diagonal of the adjacency matrix before the
124 /// decomposition; either one value per vertex, a single value for all,
125 /// or `None` for zero. A common choice is half the degrees.
126 /// - `options`: ARPACK options; only `tol` and `mxiter` are used (igraph
127 /// sets `which`, `nev` and `ncv`, and starts from a random vector of the
128 /// calling thread's RNG).
129 ///
130 /// Binds [`igraph_adjacency_spectral_embedding`](https://igraph.org/c/html/latest/igraph-Embedding.html#igraph_adjacency_spectral_embedding).
131 ///
132 /// See also [`dim_select`] to choose `no` from the
133 /// singular values, [`Graph::get_adjacency`](crate::Graph::get_adjacency)
134 /// for the matrix being decomposed, and
135 /// [`Graph::layout_mds`](crate::Graph::layout_mds) for a distance-based
136 /// spectral layout.
137 ///
138 /// # Errors
139 /// If `no` is zero or larger than the number of vertices, or the weight
140 /// or `cvec` lengths are wrong.
141 ///
142 /// # Examples
143 ///
144 /// Two disjoint 4-cliques: the 2-dimensional embedding places the two
145 /// groups on orthogonal axes.
146 ///
147 /// ```
148 /// use igraph::{linalg::*, prelude::*};
149 /// let k4 = Graph::full(4, false, false).unwrap();
150 /// let g = k4.disjoint_union(&k4).unwrap(); // vertices 0..4 and 4..8
151 /// let e = g
152 /// .adjacency_spectral_embedding(2, None, EmbeddingWhich::LargestAlgebraic, true, None, &ArpackOptions::default())
153 /// .unwrap();
154 /// assert!((e.d[0] - 3.0).abs() < 1e-8 && (e.d[1] - 3.0).abs() < 1e-8);
155 /// let dot = |a: usize, b: usize| e.x.row(a).iter().zip(e.x.row(b)).map(|(p, q)| p * q).sum::<f64>();
156 /// assert!(dot(0, 1) > 0.5); // same clique: similar positions
157 /// assert!(dot(0, 5).abs() < 1e-8); // different cliques: orthogonal
158 /// ```
159 #[allow(clippy::too_many_arguments)]
160 pub fn adjacency_spectral_embedding(
161 &self,
162 no: usize,
163 weights: Option<&[f64]>,
164 which: EmbeddingWhich,
165 scaled: bool,
166 cvec: Option<&[f64]>,
167 options: &ArpackOptions,
168 ) -> Result<SpectralEmbedding> {
169 check_common(self, no, weights)?;
170 let n = self.vcount();
171 // igraph accepts a single-element `cvec` but then reads it as if it
172 // had one element per vertex: always expand it.
173 let cvec: Vec<f64> = match cvec {
174 None => vec![0.0; n],
175 Some([c]) => vec![*c; n],
176 Some(c) if c.len() == n => c.to_vec(),
177 Some(_) => {
178 return Err(Error::invalid(
179 "cvec must have one element per vertex, or a single element",
180 ));
181 }
182 };
183 check_finite(&cvec, "cvec")?;
184 let w = weights.map(Vector::view);
185 let wp = w.as_ref().map_or(ptr::null(), |v| v.as_ptr());
186 let cv = Vector::view(&cvec);
187 let mut raw_opts = options.to_raw_unchecked_which(n)?;
188 let (mut x, mut y, mut d) = (Matrix::new(), Matrix::new(), Vector::new());
189 let _arpack = super::arpack::ArpackGuard::enter()?;
190 igraph_call!(igraph_adjacency_spectral_embedding(
191 self,
192 no as igraph_int_t,
193 wp,
194 which.into(),
195 scaled,
196 &mut x,
197 &mut y,
198 &mut d,
199 cv.as_ptr(),
200 &mut raw_opts
201 ))?;
202 Ok(SpectralEmbedding { x, y, d: d.into() })
203 }
204
205 /// Laplacian spectral embedding (`igraph_laplacian_spectral_embedding`).
206 ///
207 /// Like [`adjacency_spectral_embedding`](Self::adjacency_spectral_embedding),
208 /// but decomposes a Laplacian of the graph, chosen with `kind` (see
209 /// [`LaplacianEmbeddingType`]; directed graphs need
210 /// [`LaplacianEmbeddingType::OAP`]). Use
211 /// [`EmbeddingWhich::SmallestAlgebraic`] with `D - A` for the classic
212 /// spectral clustering / Fiedler vector embedding.
213 ///
214 /// Binds [`igraph_laplacian_spectral_embedding`](https://igraph.org/c/html/latest/igraph-Embedding.html#igraph_laplacian_spectral_embedding).
215 ///
216 /// The eigenvectors come from ARPACK, started from a random vector of
217 /// the calling thread's igraph RNG (seed it with
218 /// [`rng::seed`](crate::rng::seed) for reproducible runs). On tiny graphs
219 /// with repeated eigenvalues ARPACK may miss an eigenvalue or fail to
220 /// converge for some start vectors; there, diagonalize the dense
221 /// Laplacian with [`lapack_dsyevr`](super::lapack_dsyevr) instead.
222 ///
223 /// See also [`Graph::get_laplacian`](crate::Graph::get_laplacian) and
224 /// [`Graph::get_laplacian_sparse`](crate::Graph::get_laplacian_sparse)
225 /// (structural module) for the Laplacian matrices themselves; `D - A` is
226 /// [`LaplacianNormalization::Unnormalized`](crate::structural::LaplacianNormalization::Unnormalized)
227 /// and `I - D^-1/2 A D^-1/2` is
228 /// [`LaplacianNormalization::Symmetric`](crate::structural::LaplacianNormalization::Symmetric).
229 ///
230 /// # Errors
231 /// As for the adjacency embedding, plus an invalid `kind` for the
232 /// directedness of the graph, and — for the normalized Laplacians — a
233 /// vertex with zero strength (isolated vertex; for directed graphs, zero
234 /// in- or out-strength).
235 ///
236 /// # Examples
237 ///
238 /// Spectral bisection of Zachary's karate club: the Fiedler vector
239 /// (eigenvalue `~0.4685` of `D - A`, the algebraic connectivity) puts the
240 /// instructor (vertex 0) and the administrator (vertex 33) on opposite
241 /// sides.
242 ///
243 /// ```
244 /// use igraph::{linalg::*, prelude::*};
245 /// let g = Graph::famous("Zachary").unwrap();
246 /// let e = g
247 /// .laplacian_spectral_embedding(
248 /// 2,
249 /// None,
250 /// EmbeddingWhich::SmallestAlgebraic,
251 /// LaplacianEmbeddingType::DA,
252 /// false,
253 /// &ArpackOptions::default(),
254 /// )
255 /// .unwrap();
256 /// // One eigenvalue is 0 (connected graph), the other is the algebraic connectivity.
257 /// let fiedler = if e.d[0] > e.d[1] { 0 } else { 1 };
258 /// assert!(e.d[1 - fiedler].abs() < 1e-8);
259 /// assert!((e.d[fiedler] - 0.4685).abs() < 1e-4);
260 /// let f = e.x.column(fiedler);
261 /// assert!(f[0] * f[33] < 0.0);
262 /// ```
263 #[allow(clippy::too_many_arguments)]
264 pub fn laplacian_spectral_embedding(
265 &self,
266 no: usize,
267 weights: Option<&[f64]>,
268 which: EmbeddingWhich,
269 kind: LaplacianEmbeddingType,
270 scaled: bool,
271 options: &ArpackOptions,
272 ) -> Result<SpectralEmbedding> {
273 check_common(self, no, weights)?;
274 check_strengths(self, weights, kind)?;
275 let w = weights.map(Vector::view);
276 let wp = w.as_ref().map_or(ptr::null(), |v| v.as_ptr());
277 let mut raw_opts = options.to_raw_unchecked_which(self.vcount())?;
278 let (mut x, mut y, mut d) = (Matrix::new(), Matrix::new(), Vector::new());
279 let _arpack = super::arpack::ArpackGuard::enter()?;
280 igraph_call!(igraph_laplacian_spectral_embedding(
281 self,
282 no as igraph_int_t,
283 wp,
284 which.into(),
285 kind.into(),
286 scaled,
287 &mut x,
288 &mut y,
289 &mut d,
290 &mut raw_opts
291 ))?;
292 Ok(SpectralEmbedding { x, y, d: d.into() })
293 }
294}
295
296/// Dimensionality selection by profile likelihood (`igraph_dim_select`).
297///
298/// Given the (decreasingly) ordered "importance" values of the dimensions —
299/// typically the singular values of an adjacency spectral embedding — it
300/// models them as a mixture of two Gaussians with different means and equal
301/// variance, and returns the number `d` of leading values that maximizes the
302/// likelihood when the first `d` values are assigned to one component and
303/// the rest to the other. Also usable for any "where is the gap" problem.
304/// Time complexity: O(n).
305///
306/// Binds [`igraph_dim_select`](https://igraph.org/c/html/latest/igraph-Embedding.html#igraph_dim_select).
307///
308/// When no profile likelihood is a finite number — e.g. when all the values
309/// are equal, so that every one of them is NaN — igraph 1.0.0 and 1.0.1
310/// leave the result unset; this wrapper then returns `sv.len()` (all values
311/// in one group).
312///
313/// # Errors
314/// If `sv` is empty or contains NaN or infinite values.
315///
316/// # Examples
317///
318/// ```
319/// use igraph::linalg::dim_select;
320/// assert_eq!(dim_select(&[10.0, 9.8, 9.9, 1.0, 1.1, 0.9, 1.0]).unwrap(), 3);
321/// assert_eq!(dim_select(&[5.0]).unwrap(), 1);
322/// // No gap at all: everything is one group.
323/// assert_eq!(dim_select(&[2.0, 2.0, 2.0]).unwrap(), 3);
324/// ```
325pub fn dim_select(sv: &[f64]) -> Result<usize> {
326 if sv.is_empty() {
327 return Err(Error::invalid(
328 "dimensionality selection needs at least one value",
329 ));
330 }
331 check_finite(sv, "the values")?;
332 let v = Vector::view(sv);
333 let mut dim: igraph_int_t = 0;
334 igraph_call!(igraph_dim_select(v.as_ptr(), &mut dim))?;
335 Ok(if dim >= 1 { dim as usize } else { sv.len() })
336}