spectral: new package for spectral graph analysis
Initial commit to this is just a move of the Laplacian type from network.
diff --git a/graph/network/diffusion.go b/graph/network/diffusion.go
index 6aec9c0..b9a99de 100644
--- a/graph/network/diffusion.go
+++ b/graph/network/diffusion.go
@@ -5,9 +5,7 @@
package network
import (
- "math"
-
- "gonum.org/v1/gonum/graph"
+ "gonum.org/v1/gonum/graph/spectral"
"gonum.org/v1/gonum/mat"
)
@@ -23,7 +21,7 @@
// Nodes without corresponding entries in h are given an initial heat of zero,
// and entries in h without a corresponding node in the original graph are
// not altered when written to dst.
-func Diffuse(dst, h map[int64]float64, by Laplacian, t float64) map[int64]float64 {
+func Diffuse(dst, h map[int64]float64, by spectral.Laplacian, t float64) map[int64]float64 {
heat := make([]float64, len(by.Index))
for id, i := range by.Index {
heat[i] = h[id]
@@ -58,7 +56,7 @@
// Nodes without corresponding entries in h are given an initial heat of zero,
// and entries in h without a corresponding node in the original graph are
// not altered when written to dst.
-func DiffuseToEquilibrium(dst, h map[int64]float64, by Laplacian, tol float64, iters int) (eq map[int64]float64, ok bool) {
+func DiffuseToEquilibrium(dst, h map[int64]float64, by spectral.Laplacian, tol float64, iters int) (eq map[int64]float64, ok bool) {
heat := make([]float64, len(by.Index))
for id, i := range by.Index {
heat[i] = h[id]
@@ -94,119 +92,3 @@
}
return dst, ok
}
-
-// Laplacian is a graph Laplacian matrix.
-type Laplacian struct {
- // Matrix holds the Laplacian matrix.
- mat.Matrix
-
- // Nodes holds the input graph nodes.
- Nodes []graph.Node
-
- // Index is a mapping from the graph
- // node IDs to row and column indices.
- Index map[int64]int
-}
-
-// NewLaplacian returns a Laplacian matrix for the simple undirected graph g.
-// The Laplacian is defined as D-A where D is a diagonal matrix holding the
-// degree of each node and A is the graph adjacency matrix of the input graph.
-// If g contains self edges, NewLaplacian will panic.
-func NewLaplacian(g graph.Undirected) Laplacian {
- nodes := graph.NodesOf(g.Nodes())
- indexOf := make(map[int64]int, len(nodes))
- for i, n := range nodes {
- id := n.ID()
- indexOf[id] = i
- }
-
- l := mat.NewSymDense(len(nodes), nil)
- for j, u := range nodes {
- uid := u.ID()
- to := graph.NodesOf(g.From(uid))
- l.SetSym(j, j, float64(len(to)))
- for _, v := range to {
- vid := v.ID()
- if uid == vid {
- panic("network: self edge in graph")
- }
- if uid < vid {
- l.SetSym(indexOf[vid], j, -1)
- }
- }
- }
-
- return Laplacian{Matrix: l, Nodes: nodes, Index: indexOf}
-}
-
-// NewSymNormLaplacian returns a symmetric normalized Laplacian matrix for the
-// simple undirected graph g.
-// The normalized Laplacian is defined as I-D^(-1/2)AD^(-1/2) where D is a
-// diagonal matrix holding the degree of each node and A is the graph adjacency
-// matrix of the input graph.
-// If g contains self edges, NewSymNormLaplacian will panic.
-func NewSymNormLaplacian(g graph.Undirected) Laplacian {
- nodes := graph.NodesOf(g.Nodes())
- indexOf := make(map[int64]int, len(nodes))
- for i, n := range nodes {
- id := n.ID()
- indexOf[id] = i
- }
-
- l := mat.NewSymDense(len(nodes), nil)
- for j, u := range nodes {
- uid := u.ID()
- to := graph.NodesOf(g.From(uid))
- if len(to) == 0 {
- continue
- }
- l.SetSym(j, j, 1)
- squdeg := math.Sqrt(float64(len(to)))
- for _, v := range to {
- vid := v.ID()
- if uid == vid {
- panic("network: self edge in graph")
- }
- if uid < vid {
- l.SetSym(indexOf[vid], j, -1/(squdeg*math.Sqrt(float64(g.From(vid).Len()))))
- }
- }
- }
-
- return Laplacian{Matrix: l, Nodes: nodes, Index: indexOf}
-}
-
-// NewRandomWalkLaplacian returns a damp-scaled random walk Laplacian matrix for
-// the simple graph g.
-// The random walk Laplacian is defined as I-D^(-1)A where D is a diagonal matrix
-// holding the degree of each node and A is the graph adjacency matrix of the input
-// graph.
-// If g contains self edges, NewRandomWalkLaplacian will panic.
-func NewRandomWalkLaplacian(g graph.Graph, damp float64) Laplacian {
- nodes := graph.NodesOf(g.Nodes())
- indexOf := make(map[int64]int, len(nodes))
- for i, n := range nodes {
- id := n.ID()
- indexOf[id] = i
- }
-
- l := mat.NewDense(len(nodes), len(nodes), nil)
- for j, u := range nodes {
- uid := u.ID()
- to := graph.NodesOf(g.From(uid))
- if len(to) == 0 {
- continue
- }
- l.Set(j, j, 1-damp)
- rudeg := (damp - 1) / float64(len(to))
- for _, v := range to {
- vid := v.ID()
- if uid == vid {
- panic("network: self edge in graph")
- }
- l.Set(indexOf[vid], j, rudeg)
- }
- }
-
- return Laplacian{Matrix: l, Nodes: nodes, Index: indexOf}
-}
diff --git a/graph/network/diffusion_test.go b/graph/network/diffusion_test.go
index a3f6b05..1deac06 100644
--- a/graph/network/diffusion_test.go
+++ b/graph/network/diffusion_test.go
@@ -6,15 +6,12 @@
import (
"math"
- "sort"
"testing"
"gonum.org/v1/gonum/floats"
"gonum.org/v1/gonum/graph"
- "gonum.org/v1/gonum/graph/internal/ordered"
- "gonum.org/v1/gonum/graph/iterator"
"gonum.org/v1/gonum/graph/simple"
- "gonum.org/v1/gonum/mat"
+ "gonum.org/v1/gonum/graph/spectral"
)
var diffuseTests = []struct {
@@ -160,7 +157,7 @@
g.SetEdge(simple.Edge{F: simple.Node(u), T: simple.Node(v)})
}
}
- for j, lfn := range []func(g graph.Undirected) Laplacian{NewLaplacian, NewSymNormLaplacian} {
+ for j, lfn := range []func(g graph.Undirected) spectral.Laplacian{spectral.NewLaplacian, spectral.NewSymNormLaplacian} {
normalize := j == 1
var wantTemp float64
for _, v := range test.h {
@@ -194,123 +191,6 @@
}
}
-var randomWalkLaplacianTests = []struct {
- g []set
- damp float64
-
- want *mat.Dense
-}{
- {
- g: []set{
- A: linksTo(B, C),
- B: linksTo(C),
- C: nil,
- },
-
- want: mat.NewDense(3, 3, []float64{
- 1, 0, 0,
- -0.5, 1, 0,
- -0.5, -1, 0,
- }),
- },
- {
- g: []set{
- A: linksTo(B, C),
- B: linksTo(C),
- C: nil,
- },
- damp: 0.85,
-
- want: mat.NewDense(3, 3, []float64{
- 0.15, 0, 0,
- -0.075, 0.15, 0,
- -0.075, -0.15, 0,
- }),
- },
- {
- g: []set{
- A: linksTo(B),
- B: linksTo(C),
- C: linksTo(A),
- },
- damp: 0.85,
-
- want: mat.NewDense(3, 3, []float64{
- 0.15, 0, -0.15,
- -0.15, 0.15, 0,
- 0, -0.15, 0.15,
- }),
- },
- {
- // Example graph from http://en.wikipedia.org/wiki/File:PageRanks-Example.svg 16:17, 8 July 2009
- g: []set{
- A: nil,
- B: linksTo(C),
- C: linksTo(B),
- D: linksTo(A, B),
- E: linksTo(D, B, F),
- F: linksTo(B, E),
- G: linksTo(B, E),
- H: linksTo(B, E),
- I: linksTo(B, E),
- J: linksTo(E),
- K: linksTo(E),
- },
-
- want: mat.NewDense(11, 11, []float64{
- 0, 0, 0, -0.5, 0, 0, 0, 0, 0, 0, 0,
- 0, 1, -1, -0.5, -1. / 3., -0.5, -0.5, -0.5, -0.5, 0, 0,
- 0, -1, 1, 0, 0, 0, 0, 0, 0, 0, 0,
- 0, 0, 0, 1, -1. / 3., 0, 0, 0, 0, 0, 0,
- 0, 0, 0, 0, 1, -0.5, -0.5, -0.5, -0.5, -1, -1,
- 0, 0, 0, 0, -1. / 3., 1, 0, 0, 0, 0, 0,
- 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0,
- 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0,
- 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0,
- 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0,
- 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1,
- }),
- },
-}
-
-func TestRandomWalkLaplacian(t *testing.T) {
- const tol = 1e-14
- for i, test := range randomWalkLaplacianTests {
- g := simple.NewDirectedGraph()
- for u, e := range test.g {
- // Add nodes that are not defined by an edge.
- if g.Node(int64(u)) == nil {
- g.AddNode(simple.Node(u))
- }
- for v := range e {
- g.SetEdge(simple.Edge{F: simple.Node(u), T: simple.Node(v)})
- }
- }
- l := NewRandomWalkLaplacian(g, test.damp)
- _, c := l.Dims()
- for j := 0; j < c; j++ {
- if got := mat.Sum(l.Matrix.(*mat.Dense).ColView(j)); !floats.EqualWithinAbsOrRel(got, 0, tol, tol) {
- t.Errorf("unexpected column sum for test %d, column %d: got:%v want:0", i, j, got)
- }
- }
- l = NewRandomWalkLaplacian(sortedNodeGraph{g}, test.damp)
- if !mat.EqualApprox(l, test.want, tol) {
- t.Errorf("unexpected result for test %d:\ngot:\n% .2v\nwant:\n% .2v",
- i, mat.Formatted(l), mat.Formatted(test.want))
- }
- }
-}
-
-type sortedNodeGraph struct {
- graph.Graph
-}
-
-func (g sortedNodeGraph) Nodes() graph.Nodes {
- n := graph.NodesOf(g.Graph.Nodes())
- sort.Sort(ordered.ByID(n))
- return iterator.NewOrderedNodes(n)
-}
-
var diffuseToEquilibriumTests = []struct {
g []set
builder builder
@@ -478,7 +358,7 @@
for _, v := range test.h {
wantTemp += v
}
- got, ok := DiffuseToEquilibrium(nil, test.h, NewRandomWalkLaplacian(g, test.damp), test.tol*test.tol, test.iter)
+ got, ok := DiffuseToEquilibrium(nil, test.h, spectral.NewRandomWalkLaplacian(g, test.damp), test.tol*test.tol, test.iter)
if ok != test.wantOK {
t.Errorf("unexpected success value for test %d: got:%t want:%t", i, ok, test.wantOK)
}
diff --git a/graph/spectral/doc.go b/graph/spectral/doc.go
new file mode 100644
index 0000000..b26425f
--- /dev/null
+++ b/graph/spectral/doc.go
@@ -0,0 +1,6 @@
+// Copyright ©2019 The Gonum Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+// Package spectral provides graph spectral analysis functions.
+package spectral // import "gonum.org/v1/gonum/graph/spectral"
diff --git a/graph/spectral/laplacian.go b/graph/spectral/laplacian.go
new file mode 100644
index 0000000..0c3bd41
--- /dev/null
+++ b/graph/spectral/laplacian.go
@@ -0,0 +1,128 @@
+// Copyright ©2017 The Gonum Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package spectral
+
+import (
+ "math"
+
+ "gonum.org/v1/gonum/graph"
+ "gonum.org/v1/gonum/mat"
+)
+
+// Laplacian is a graph Laplacian matrix.
+type Laplacian struct {
+ // Matrix holds the Laplacian matrix.
+ mat.Matrix
+
+ // Nodes holds the input graph nodes.
+ Nodes []graph.Node
+
+ // Index is a mapping from the graph
+ // node IDs to row and column indices.
+ Index map[int64]int
+}
+
+// NewLaplacian returns a Laplacian matrix for the simple undirected graph g.
+// The Laplacian is defined as D-A where D is a diagonal matrix holding the
+// degree of each node and A is the graph adjacency matrix of the input graph.
+// If g contains self edges, NewLaplacian will panic.
+func NewLaplacian(g graph.Undirected) Laplacian {
+ nodes := graph.NodesOf(g.Nodes())
+ indexOf := make(map[int64]int, len(nodes))
+ for i, n := range nodes {
+ id := n.ID()
+ indexOf[id] = i
+ }
+
+ l := mat.NewSymDense(len(nodes), nil)
+ for j, u := range nodes {
+ uid := u.ID()
+ to := graph.NodesOf(g.From(uid))
+ l.SetSym(j, j, float64(len(to)))
+ for _, v := range to {
+ vid := v.ID()
+ if uid == vid {
+ panic("network: self edge in graph")
+ }
+ if uid < vid {
+ l.SetSym(indexOf[vid], j, -1)
+ }
+ }
+ }
+
+ return Laplacian{Matrix: l, Nodes: nodes, Index: indexOf}
+}
+
+// NewSymNormLaplacian returns a symmetric normalized Laplacian matrix for the
+// simple undirected graph g.
+// The normalized Laplacian is defined as I-D^(-1/2)AD^(-1/2) where D is a
+// diagonal matrix holding the degree of each node and A is the graph adjacency
+// matrix of the input graph.
+// If g contains self edges, NewSymNormLaplacian will panic.
+func NewSymNormLaplacian(g graph.Undirected) Laplacian {
+ nodes := graph.NodesOf(g.Nodes())
+ indexOf := make(map[int64]int, len(nodes))
+ for i, n := range nodes {
+ id := n.ID()
+ indexOf[id] = i
+ }
+
+ l := mat.NewSymDense(len(nodes), nil)
+ for j, u := range nodes {
+ uid := u.ID()
+ to := graph.NodesOf(g.From(uid))
+ if len(to) == 0 {
+ continue
+ }
+ l.SetSym(j, j, 1)
+ squdeg := math.Sqrt(float64(len(to)))
+ for _, v := range to {
+ vid := v.ID()
+ if uid == vid {
+ panic("network: self edge in graph")
+ }
+ if uid < vid {
+ l.SetSym(indexOf[vid], j, -1/(squdeg*math.Sqrt(float64(g.From(vid).Len()))))
+ }
+ }
+ }
+
+ return Laplacian{Matrix: l, Nodes: nodes, Index: indexOf}
+}
+
+// NewRandomWalkLaplacian returns a damp-scaled random walk Laplacian matrix for
+// the simple graph g.
+// The random walk Laplacian is defined as I-D^(-1)A where D is a diagonal matrix
+// holding the degree of each node and A is the graph adjacency matrix of the input
+// graph.
+// If g contains self edges, NewRandomWalkLaplacian will panic.
+func NewRandomWalkLaplacian(g graph.Graph, damp float64) Laplacian {
+ nodes := graph.NodesOf(g.Nodes())
+ indexOf := make(map[int64]int, len(nodes))
+ for i, n := range nodes {
+ id := n.ID()
+ indexOf[id] = i
+ }
+
+ l := mat.NewDense(len(nodes), len(nodes), nil)
+ for j, u := range nodes {
+ uid := u.ID()
+ to := graph.NodesOf(g.From(uid))
+ if len(to) == 0 {
+ continue
+ }
+ l.Set(j, j, 1-damp)
+ rudeg := (damp - 1) / float64(len(to))
+ for _, v := range to {
+ vid := v.ID()
+ if uid == vid {
+ panic("network: self edge in graph")
+ }
+ l.Set(indexOf[vid], j, rudeg)
+ }
+ }
+
+ return Laplacian{Matrix: l, Nodes: nodes, Index: indexOf}
+}
diff --git a/graph/spectral/laplacian_test.go b/graph/spectral/laplacian_test.go
new file mode 100644
index 0000000..6cdaea1
--- /dev/null
+++ b/graph/spectral/laplacian_test.go
@@ -0,0 +1,177 @@
+// Copyright ©2017 The Gonum Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style
+// license that can be found in the LICENSE file.
+
+package spectral
+
+import (
+ "sort"
+ "testing"
+
+ "gonum.org/v1/gonum/floats"
+ "gonum.org/v1/gonum/graph"
+ "gonum.org/v1/gonum/graph/internal/ordered"
+ "gonum.org/v1/gonum/graph/iterator"
+ "gonum.org/v1/gonum/graph/simple"
+ "gonum.org/v1/gonum/mat"
+)
+
+var randomWalkLaplacianTests = []struct {
+ g []set
+ damp float64
+
+ want *mat.Dense
+}{
+ {
+ g: []set{
+ A: linksTo(B, C),
+ B: linksTo(C),
+ C: nil,
+ },
+
+ want: mat.NewDense(3, 3, []float64{
+ 1, 0, 0,
+ -0.5, 1, 0,
+ -0.5, -1, 0,
+ }),
+ },
+ {
+ g: []set{
+ A: linksTo(B, C),
+ B: linksTo(C),
+ C: nil,
+ },
+ damp: 0.85,
+
+ want: mat.NewDense(3, 3, []float64{
+ 0.15, 0, 0,
+ -0.075, 0.15, 0,
+ -0.075, -0.15, 0,
+ }),
+ },
+ {
+ g: []set{
+ A: linksTo(B),
+ B: linksTo(C),
+ C: linksTo(A),
+ },
+ damp: 0.85,
+
+ want: mat.NewDense(3, 3, []float64{
+ 0.15, 0, -0.15,
+ -0.15, 0.15, 0,
+ 0, -0.15, 0.15,
+ }),
+ },
+ {
+ // Example graph from http://en.wikipedia.org/wiki/File:PageRanks-Example.svg 16:17, 8 July 2009
+ g: []set{
+ A: nil,
+ B: linksTo(C),
+ C: linksTo(B),
+ D: linksTo(A, B),
+ E: linksTo(D, B, F),
+ F: linksTo(B, E),
+ G: linksTo(B, E),
+ H: linksTo(B, E),
+ I: linksTo(B, E),
+ J: linksTo(E),
+ K: linksTo(E),
+ },
+
+ want: mat.NewDense(11, 11, []float64{
+ 0, 0, 0, -0.5, 0, 0, 0, 0, 0, 0, 0,
+ 0, 1, -1, -0.5, -1. / 3., -0.5, -0.5, -0.5, -0.5, 0, 0,
+ 0, -1, 1, 0, 0, 0, 0, 0, 0, 0, 0,
+ 0, 0, 0, 1, -1. / 3., 0, 0, 0, 0, 0, 0,
+ 0, 0, 0, 0, 1, -0.5, -0.5, -0.5, -0.5, -1, -1,
+ 0, 0, 0, 0, -1. / 3., 1, 0, 0, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0,
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0,
+ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1,
+ }),
+ },
+}
+
+func TestRandomWalkLaplacian(t *testing.T) {
+ const tol = 1e-14
+ for i, test := range randomWalkLaplacianTests {
+ g := simple.NewDirectedGraph()
+ for u, e := range test.g {
+ // Add nodes that are not defined by an edge.
+ if g.Node(int64(u)) == nil {
+ g.AddNode(simple.Node(u))
+ }
+ for v := range e {
+ g.SetEdge(simple.Edge{F: simple.Node(u), T: simple.Node(v)})
+ }
+ }
+ l := NewRandomWalkLaplacian(g, test.damp)
+ _, c := l.Dims()
+ for j := 0; j < c; j++ {
+ if got := mat.Sum(l.Matrix.(*mat.Dense).ColView(j)); !floats.EqualWithinAbsOrRel(got, 0, tol, tol) {
+ t.Errorf("unexpected column sum for test %d, column %d: got:%v want:0", i, j, got)
+ }
+ }
+ l = NewRandomWalkLaplacian(sortedNodeGraph{g}, test.damp)
+ if !mat.EqualApprox(l, test.want, tol) {
+ t.Errorf("unexpected result for test %d:\ngot:\n% .2v\nwant:\n% .2v",
+ i, mat.Formatted(l), mat.Formatted(test.want))
+ }
+ }
+}
+
+type sortedNodeGraph struct {
+ graph.Graph
+}
+
+func (g sortedNodeGraph) Nodes() graph.Nodes {
+ n := graph.NodesOf(g.Graph.Nodes())
+ sort.Sort(ordered.ByID(n))
+ return iterator.NewOrderedNodes(n)
+}
+
+const (
+ A = iota
+ B
+ C
+ D
+ E
+ F
+ G
+ H
+ I
+ J
+ K
+ L
+ M
+ N
+ O
+ P
+ Q
+ R
+ S
+ T
+ U
+ V
+ W
+ X
+ Y
+ Z
+)
+
+// set is an integer set.
+type set map[int64]struct{}
+
+func linksTo(i ...int64) set {
+ if len(i) == 0 {
+ return nil
+ }
+ s := make(set)
+ for _, v := range i {
+ s[v] = struct{}{}
+ }
+ return s
+}
