Building a graph from a matrix

Given an upper triangular matrix B that is ordered by the atoms in each molecule such that the atoms in molecule 1 come first, then molecule 2, etc., I want to check if any atoms in any two molecules meet some criteria. Then, I want to designate an edge between those two molecules, but being careful not to count more than one edge between the same two molecules. A molecule can form edges with many other molecules, though. How do I do this efficiently?

Something like the following?

N = 200; % number of molecules
G = graph;
G = addnode(G,N);
C = ((B < 1.75) & (B > 0));
% Here I'm not sure how to efficiently determine the molecules that meet those criteria
G = addedge(G,a,b); % where a and b are all molecules a and b (I guess I need a double loop?)

Answer 1

[Matt J]

Given an upper triangular matrix B

Isn't B your adjacency matrix? If so, then why not just,

G=graph(B);

[L'O.G.]

Interesting! But I don't want to construct the graph of all of the atoms in all molecules. (B provides distances between each pair of atoms.) Instead, I want to graph just the molecules themselves as the nodes. What you suggested does, in a sense, the opposite of what I want.

[Matt J]

Maybe with conncomp, you could get the adjacency matrix of the molecules.

bins=conncomp(graph(C));

S=sparse(bins,1:numnodes(graph(B)),1);

G=graph(S*B*S'>0)

[L'O.G.]

@Matt J Thanks. That gives an error that the adjacency matrix isn't symmetric. The error is with the first line. What do you suggest?

[Matt J]

B provides distances between each pair of atoms

If B and C are distance matrices, I wonder why they aren't symmetric?

[Steven Lord]

Perhaps you should show us a concrete example with a small number of atoms. Show us a specific 6-by-6 matrix B (or if 6 isn't sufficient, some number of atoms smaller than 200) and walk us through how you determine if two of the nodes represented by that matrix should be connected.

[L'O.G.]

Thanks both. Before walking you through an example, I should say B is upper triangular to avoid double counting distances (A -> B is the same as B -> A). I'll go through a little example in a bit.

[Matt J]

I should say B is upper triangular to avoid double counting distances (A -> B is the same as B -> A)

You should just make it symmetric.

B=B+B.';

graph() is smart enough to know that A -> B is the same as B -> A. If that were not the case, you would use digraph().

[L'O.G.]

@Matt J Hm, that doesn't work: I get a graph with 1 node only. It also takes around 5 minutes to run for my B matrix

I should also clarify (I thought I did this, but maybe not very well) that B is a matrix of distances between each pair of atoms in the entire system. For the graph, however, I want to look at the nodes of just molecules (so # nodes = N rather than # nodes = N x # atoms on each molecule). The edges are defined based on whether or not a molecule is adjacent to another molecule, which is why I need the matrix B: two molecules are adjacent if any of their atoms are adjacent. Also, only one edge is possible between any two molecules, so even if more than one pair of atoms is adjacent, I still want to define only one edge. I hope that's clear now.

[Matt J]

That was always clear. We probably need to see some example B and C matrices as Steve said.

[L'O.G.]

@Matt J @Steven Lord Thanks for the feedback. Here's an example that is just slightly larger than what Steve suggested: 8 x 8 rather than 6 x 6. There are 4 molecules, each consisting of 2 sites. Atoms 1 and 2 (same molecule) are adjacent. So are atoms 2 and 3 (different molecule), 3 and 4 (same molecule), 5 and 6 (same molecule), 1 and 6 (different molecule), 2 and 6 (different molecule), and 7 and 8 (same molecule). I hope this helps. Here I would like to make a graph with 4 nodes, and there would be edges between nodes 1 and 2, 1 and 3 I believe. I hope that helps.

[Matt J]

It would be easiest if you attach the B and C matrices for the simplified problem.

[Steven Lord]

From your description it almost sounds like you should create a digraph representing the connections between your atoms and then use condensation to collapse the strongly connected components representing your molecules into individual nodes. Does the example on the documentation page look something like what you had in mind?

[L'O.G.]

@Matt J Sorry, I thought I attached the data! Ugh. Here are the matrices.

[L'O.G.]

@Steven Lord Yes, that sounds useful: representing each component as a single node. I'll play around with it. Looking at it, I don't think I can say anything about the strongly connecting components in my system, at least for this particular problem. Regardless of how many adjacent atoms there are, I define an edge as long as at least one of them is adjacent between the two molecules. I feel that condensation won't be able to capture that, but I'll dig into the help page some more. Thanks for letting me know about it.

[Matt J]

The B you have attached is neither upper triangular nor symmetric.

[Matt J]

Also, you say that atoms 2 and 3 should belong to different molecules, but this contradicts the fact that C(2,3)=1.

[L'O.G.]

@Matt J Weird, sorry, I could've sworn I sent the right files. I tripled checked these new ones: they're symmetric. Also, I took B=B+B.' as you suggested to get a symmetric matrix, where the original B was upper triangular. Indeed, C(2,3) = 1. That means atoms 2 and 3 are adjacent, which is perfectly ok because atoms can be adjacent regardless of whether or not they're on the same molecule. Here they are on different molecules. The atoms are in order for the molecules that they're in: atoms 1 and 2 = molecule 1, atoms 3 and 4 = molecule 2, and so on. I hope this is clear / makes sense.

[Matt J]

Yes, but what is the general criterion for determining whether two atoms belong to the same molecule? Does it have anything to do with the values of B and C? In your example, we don't seem to need to determine the groupings into molecules from B and C, because you've simply told us that consecutive pairs of atoms always form a separate molecule, but will this always be the case?

[L'O.G.]

Consecutive atoms will always be in the same molecule up to the number of atoms fixed for each molecule. It's not necessarily going to be pairs. Here it's pairs just due to the simple example, but you could have 100's of atoms in a molecule. Those atoms would all be in order, as I mentioned. It should always be just one kind of molecule, though, at least I think.

[Matt J]

OK, well, then in addition to C, you need an additional input which is a label vector, L, stating which molecule the atoms belong to. For the simplified example, this would be,

L=[1     1     2     2     3     3     4     4];

and adapting my approach from above leads to,

S=sparse(L,1:numel(L),1);

A=S*C*S'>0;
A(1:size(A,1)+1:end)=0;  %molecule adjacency matrix.

plot(graph(A))

In an earlier comment, you predicted that there would be additional edges 1-3 and 1-4, but I think that must have been an error. There is no adjacency, according to your provided C matrix, between the atoms of molecules 1 and 3 or 1 and 4.

[L'O.G.]

Thank you! Can you please explain the reason for A=S*C*S'>0; ? I will accept your answer in a few seconds -- thank you for your patience / walking me through this!