getAllSubGraphs Get all metabolic subgraphs in a model. Two metabolites are connected if they share a reaction. model a model structure subGraphs a boolean matrix where the rows correspond to the metabolites and the columns to which subgraph they are assigned to. The columns are ordered so that larger subgraphs come first Usage: subGraphs=getAllSubGraphs(model) Rasmus Agren, 2013-11-22
0001 function subGraphs=getAllSubGraphs(model) 0002 % getAllSubGraphs 0003 % Get all metabolic subgraphs in a model. Two metabolites 0004 % are connected if they share a reaction. 0005 % 0006 % model a model structure 0007 % 0008 % subGraphs a boolean matrix where the rows correspond to the metabolites 0009 % and the columns to which subgraph they are assigned to. The 0010 % columns are ordered so that larger subgraphs come first 0011 % 0012 % Usage: subGraphs=getAllSubGraphs(model) 0013 % 0014 % Rasmus Agren, 2013-11-22 0015 % 0016 0017 %Generate the connectivity graph. Metabolites are connected through 0018 %reactions. This is not a bipartite graph with the reactions. 0019 G=model.S*model.S'; 0020 G(G~=0)=1; 0021 0022 %Keeps track of which mets have been assigned to a subgraph 0023 isAssigned=false(numel(model.mets),1); 0024 0025 %Allocate space for 100 subgraphs 0026 subGraphs=false(numel(model.mets),100); 0027 0028 %Main loop continues until all mets have been assigned to a subgraph 0029 counter=1; 0030 while ~all(isAssigned) 0031 currentSG=false(numel(model.mets),1); 0032 %Find the first non-assigned metabolite and assign it to the current SG 0033 I=find(isAssigned==false,1); 0034 currentSG(I)=true; 0035 isAssigned(I)=true; 0036 %Then iteratively add all mets that are connected to each other, until 0037 %no more such mets can be found 0038 while true 0039 J=sum(G(:,currentSG),2); 0040 0041 %Get the new mets for this SG. Also assign them to the current SG 0042 newOnes=J~=0 & currentSG==false; 0043 isAssigned(newOnes)=true; 0044 currentSG(newOnes)=true; 0045 0046 %If there are no new mets, then abort 0047 if ~any(newOnes) 0048 subGraphs(currentSG,counter)=true; 0049 counter=counter+1; 0050 break; 0051 end 0052 end 0053 end 0054 subGraphs=subGraphs(:,1:counter); 0055 0056 [crap, I]=sort(sum(subGraphs),'descend'); 0057 0058 subGraphs=subGraphs(:,I); 0059 0060 %Also remove empty subgraphs (can happen when metabolites are never used) 0061 subGraphs(:,sum(subGraphs)==0)=[]; 0062 end