13
$\begingroup$

Related Wolfram Community question.


I am looking for the fastest way to achieve the following:

Given a Graph, retrieve its edge list in terms of vertex indices (not actual vertices).


For example, given

Graph[{a, b, c}, {a <-> b, a <-> c}]

I am looking to get the output

{{1,2}, {1,3}}

as a packed array. Alternatively the flattened version would do just as well:

{1,2, 1,3}

Here {1,2} corresponds to UndirectedEdge[a,b] as a has vertex index 1 and b has vertex index 2.


What I have so far is the very straightforward

idxEdgeList[graph_] := 
 Developer`ToPackedArray[
  List @@@ EdgeList[graph] /. AssociationThread[VertexList[graph] -> Range@VertexCount[graph]]
 ]

g = GridGraph[{250,250}];

idxEdgeList[g]; // AbsoluteTiming

(* 0.27 seconds *)

Using undocumented features or poking inside of the Graph object is okay for as long as the method is proven to be reliable for 10.0 – 10.2 for various directed and undirected, simple and non-simple graphs. Multigraphs (multiple edges between the same vertices) must be supported, but mixed graphs (both directed and undirected edges) do not. A documented way is of course always preferred!

This is admittedly a fairly boring performance tuning problem, but this turned out to be a bottleneck in some cases, and I don't want to lose out on any possible performance improvements I may have missed.

Use case: The edge list will eventually be passed to a LibraryLink function. What hasn't occurred to me before typing up the question is that maybe I should be using sparse arrays, which are directly supported by LibraryLink.


Update:

The solution proposed by @halmir, through IndexGraph, works well for the GridGraph above. But it is not fast for all graphs. In particular:

g = GridGraph[{250, 250}];

IndexGraph[g]; // AbsoluteTiming
(* {4.*10^-6, Null} *)

g = Graph[VertexList[g], EdgeList[g]];

IndexGraph[g]; // AbsoluteTiming
(* {0.259276, Null} *)

We are now back to the same speed as the Replace method. Re-creating the graph from its vertex and edge lists somehow made IndexGraph be slow on it, and no matter what I try I cannot convert the graph back to a "fast" format.

The SparseArray-based method is much faster, and proves that it is technically possible to extract the information quickly. But it has a big problem: it does not preserve the edge order, which means that I cannot match up the edges with an EdgeWeight vector anymore. It's also difficult to handle for multigraphs, though that would be solvable if I could preserve the ordering ...

Update / 2017

@Ramble suggests using the IncidenceMatrix of the graph. The fastest way I found so far is to process the incidence matrix in C, using LibraryLink, to extract the index-based edge list.

According to the documentation, an incidence matrix uses the following values:

  • -1 represents the starting point of a directed edge
  • 1 represents the endpoint of a directed edge or an undirected egde
  • 2 represents an undirected self-loop
  • -2 represents a directed self-loop

This is not accurate. Between 10.0-11.2, both directed and undirected self-loops are represented with a positive 2. This prevents the correct representation of mixed graphs (MixedGraphQ), but I do not need that anyway. Multigraphs are easily handled by this approach.


This is now available in IGraph/M 0.3.95 as IGEdgeIndexList.

This function is actually faster than EdgeList, and can be used to implement many edge-list based operations efficiently. An index-based edge list can be used to reconstruct a graph using the undocumented syntax Graph[vertexList, indexEdgeList], e.g. Graph[{a,b}, {{1,2}}].


Here's the LTemplate code I used for this:

mma::IntTensorRef incidenceToEdgeList(mma::SparseMatrixRef<mint> im, bool directed) {
    auto edgeList = mma::makeVector<mint>(2*im.cols());
    if (directed) {
        for (auto it = im.begin(); it != im.end(); ++it) {
            switch (*it) {
            case -1:
                edgeList[2*it.col()] = it.row();
                break;
            case  1:
                edgeList[2*it.col() + 1] = it.row();
                break;
            case  2:
            case -2:
                edgeList[2*it.col()] = it.row();
                edgeList[2*it.col() + 1] = it.row();
                break;
            default:
                throw mma::LibraryError("Invalid incidence matrix.");
            }
        }
    } else {
        for (auto &el : edgeList)
            el = -1;
        for (auto it = im.begin(); it != im.end(); ++it) {
            switch (*it) {
            case  1:
                if (edgeList[2*it.col()] == -1)
                    edgeList[2*it.col()] = it.row();
                else
                    edgeList[2*it.col() + 1] = it.row();
                break;
            case  2:
                edgeList[2*it.col()] = it.row();
                edgeList[2*it.col() + 1] = it.row();
                break;
            default:
                throw mma::LibraryError("Invalid incidence matrix.");
            }
        }
    }
    return edgeList;
}
$\endgroup$
  • 2
    $\begingroup$ I think the last paragraph answers the question: pass the sparse adjacency matrix to LibraryLink and do the conversion in C, sorting out the different cases for directed/undirected and simple/multigraphs. $\endgroup$ – Szabolcs Sep 29 '15 at 9:43
11
$\begingroup$

Using IndexGraph:

 g = GridGraph[{250, 250}];

 a = Developer`ToPackedArray[
    List @@@ EdgeList[IndexGraph[g]]]; // AbsoluteTiming

{0.063857, Null}

Using AdjacencyMatrix:

b = UpperTriangularize[AdjacencyMatrix[g]][
    "NonzeroPositions"]; // AbsoluteTiming

{0.002584, Null}

c = idxEdgeList[g]; // AbsoluteTiming

{0.276563, Null}

Test results:

Developer`PackedArrayQ /@ {a, b, c}

{True, True, True}

a == b == c

True

$\endgroup$
  • $\begingroup$ Thank you, excellent! I didn't know about IndexGraph. Just taking the NonzeroPositions from the AdjacencyMatrix is not sufficient if we have multigraphs, which is why I thought of going down the C route ... but IndexGraph makes that unnecessary. $\endgroup$ – Szabolcs Sep 29 '15 at 11:17
  • $\begingroup$ Another question: I need to test if a graph is edge-weighted. Is this a fully reliable way? edgeWeightedQ = WeightedGraphQ[#] && PropertyValue[#, EdgeList] =!= Automatic & $\endgroup$ – Szabolcs Sep 29 '15 at 11:21
  • $\begingroup$ maybe WeightedGraphQ[#] && PropertyValue[#, EdgeWeight] =!= Automatic& is better? $\endgroup$ – halmir Sep 29 '15 at 13:05
  • $\begingroup$ Isn't that exactly the same that I wrote? I think you may have accidentally copy pasted the wrong line? $\endgroup$ – Szabolcs Sep 29 '15 at 14:04
  • 1
    $\begingroup$ @Szabolcs for 1, you could use WeightedAdjacencyMatrix: g = ExampleData[{"NetworkGraph", "CondensedMatterCollaborations2005"}]; {e, w} = UpperTriangularize[ WeightedAdjacencyMatrix[g]] /@ {"NonzeroPositions", "NonzeroValues"}; e is edge indices and w is weights matching with e. $\endgroup$ – halmir Oct 7 '15 at 13:10
7
+150
$\begingroup$

Or possibly:

idxEdgeList4[graph_] := 
 Module[{edge = IncidenceMatrix[graph]}, 
  Flatten[Transpose[edge]["ColumnIndices"]]] 

When comparing solution answer with your test sequence and the suggestion from Wolfram Community (Carl Woll)

g1 = GridGraph[{250, 250, 2}];
g2 = Graph@EdgeList[g1];
g3 = GraphComputation`ToGraphRepresentation[g2, "Simple"];
g4 = GraphComputation`ToGraphRepresentation[g2, "Sparse"];
wg1 = Graph[g1, EdgeWeight -> RandomReal[1, EdgeCount[g1]]];
wg2 = Graph[g2, EdgeWeight -> RandomReal[1, EdgeCount[g2]]];

h1 = RandomGraph[{10000, 400000}];
h2 = Graph@EdgeList[h1];
h3 = GraphComputation`ToGraphRepresentation[h2, "Simple"];
h4 = GraphComputation`ToGraphRepresentation[h2, "Sparse"];
wh1 = Graph[h1, EdgeWeight -> RandomReal[1, EdgeCount[h1]]];
wh2 = Graph[h2, EdgeWeight -> RandomReal[1, EdgeCount[h2]]];

e1 = ExampleData[{"NetworkGraph", "CondensedMatterCollaborations"}];
e2 = ExampleData[{"NetworkGraph", "HighEnergyTheoryCollaborations"}];

m1 = Graph[UndirectedEdge @@@ RandomInteger[{1, 2000}, {80000, 2}]];
m2 = GraphComputation`ToGraphRepresentation[m1, "Sparse"];

graphs = AssociationThread[{"g1", "g2", "g3", "g4", "wg1", "wg2", 
    "h1", "h2", "h3", "h4", "wh1", "wh2", "e1", "e2", "m1", 
    "m2"}, {g1, g2, g3, g4, wg1, wg2, h1, h2, h3, h4, wh1, wh2, e1, 
    e2, m1, m2}];  



timings = 
       Table[First@AbsoluteTiming@idxEdgeList4[#], {5}] & /@ 
        graphs; // AbsoluteTiming
    timings // Map[Min] // Normal // Map@Apply[List] // TableForm

    Table[Partition[idxEdgeList4[graphs[[idx]]], 2] == 
      elist[graphs[[idx]]], {idx, 1, Length[graphs] - 2}]

{1.65259, Null}

{
 {"g1", 0.009295},
 {"g2", 0.012675},
 {"g3", 0.008175},
 {"g4", 0.008289},
 {"wg1", 0.009168},
 {"wg2", 0.011549},
 {"h1", 0.033966},
 {"h2", 0.047896},
 {"h3", 0.033908},
 {"h4", 0.033706},
 {"wh1", 0.033781},
 {"wh2", 0.046331},
 {"e1", 0.001917},
 {"e2", 0.000581},
 {"m1", 0.00986},
 {"m2", 0.005835}
}

{True, True, True, True, True, True, True, True, True, True, True, \ True, True, True}

Comparing solutions with:

igEdgeList[graph_] :=
Developer`ToPackedArray@If[GraphComputation`GraphRepresentation[graph] === "Simple",
  Flatten[EdgeList@IndexGraph[graph, 0], 1, If[DirectedGraphQ[graph], DirectedEdge, UndirectedEdge]]
  ,
  Lookup[
    AssociationThread[VertexList[graph], Range@VertexCount[graph] - 1],
    Flatten[EdgeList[graph], 1, If[DirectedGraphQ[graph], DirectedEdge, UndirectedEdge]]
  ]
]

resulted in:

{True, False, True, False, True, False, True, False, True, False, \ True, False, False, False}

Where the differences were mainly to do with numbering (noting that my solution didn't account for Direction) and could be fixed through reversal of solutions or adding +1 to the output of igEdgeList. I haven't checked the FullForm output for these graphs but from past experience the output algorithm orders sometimes differ from those provided by functions like "EdgeList" and further vary between others. In other words unless you check there are often hidden edge cases especially with huge graphs. I have a large file of comparisons between the output of EdgeList and those generated by Matrix Functions and indexes don't always match. e.g. Graph[h2] functions index from 1-131 and skip 1-2, 1-3, etc. I wouldn't be surprised if on some cases what is being generated and shown in Mathematica is different from that generated by the functions we are creating from first principles but who has the time to actually prove that for every case?

Elaboration: In response to Szabolics "... but from past experience the output algorithm orders sometimes differ from those provided by functions like "EdgeList" ... "

Yes, my comment was in relation to edge ordering and the general problem of validating output in Mathematica. Or rather accepting at face value that the functions work as advertised. From simple inspection I can see that the ordering may be different (why does numbering differ?). Given there is this question I am not sure as to what the correct output may be and from what function it may be derived and whether \[Equal] breaks down for certain cases or large solutions. I just wanted to be upfront about the fact I wasn't sure whether things were breaking down for your test cases as I wasn't sure which output to validate against.

In the past I had huge problems when I was trying to speed up Voronoi solutions and Meshfunctions. At the time the fastest implementation was parsing the Fullform version of the graphical output, rather than the list version... It turned out I needed to validate every analysis output, as for certain large solutions the were problems with vertice order (and both outputs differed)... More recently I submitted a bug concerning clustering on datasets and the answer came back that the issue was a known bug and the function would simply work better with large datasets and I should just trust the output.

Aside from the edge ordering there appeared to be "issues" with how functions operated across datasets that contain multiple graphs or disjointedgraphs. To check the correct output, one may need to break down the set into its disparate graphs and validate against each.

But it may all be fine. For one the complexity of some graphs is so off the charts equivalency may just be near enough is good enough. Also if your just trying to export maybe the representation does't really matter so long as you can prove that reversing it results in the original graph:)

$\endgroup$
  • 1
    $\begingroup$ This is a great idea! One difficulty is self-loops, though. $\endgroup$ – Szabolcs Dec 21 '17 at 20:29
  • $\begingroup$ I guess I could transfer the incidence matrix to LibraryLink and build the edge list in pure C code. Then handling directed graphs and self-loops wouldn't be too big of a performance hit. $\endgroup$ – Szabolcs Dec 21 '17 at 20:44
  • $\begingroup$ I implemented this in C++ (the sparse array support in LTemplate made it easy). It looks to be the best approach. It can deal with multigraphs too. I'll accept after I'm confident enough that everything works as it should. $\endgroup$ – Szabolcs Dec 21 '17 at 22:05
  • 1
    $\begingroup$ Can you elaborate on the following? ... but from past experience the output algorithm orders sometimes differ from those provided by functions like "EdgeList" ... Do you mean that sometimes IncidenceMatrix messes up the edge ordering? That would be a problem. So far I haven't noticed any problems (see my update). $\endgroup$ – Szabolcs Dec 22 '17 at 9:57
  • $\begingroup$ Well, if you come across a specific example where the function doesn't return what it should, please do let me know. I now made this a critical part of IGraph/M, so I can't afford wrong results from it. I did do a lot of testing and found no issue so far—but Mathematica is Mathematica and you never know what will happen ... $\endgroup$ – Szabolcs Jan 2 '18 at 23:57

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.