Builds Directed Enzyme-Enzyme Networks Removing Currency Metabolites (based-on a Library file) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% The function reads a Metabolic Network SBML file and builds Directed Enzyme-Enzyme Networks. The Remove Currency Metabolites based-on Library (RCMLib) algorithm removes currency metabolites in the metabolic network automatically IF AND ONLY IF the currency metabolits exists in the Library file. Note: COBRA Toolbox must be installed in MATLAB before running this function [Output] = enz_cent_RCMLib_dir(fileName1,fileName2) INPUTS fileName1 The Library file includes pre-defined currency metabolites (in .txt format) Note: Library text file must include one metabolites per line (all in one column) fileName2 The metabolic Network in the SBML format OUTPUTS *_Enzyme_Cent_RCMLib_Dir_Cyt.sif Directed-Enzyme-Enzyme Network - Cytoscape Compatible Yazdan Asgari 07/16/2016 http://yazdan59.ir/scan %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0001 function [Output] = enz_cent_RCMLib_dir(fileName1,fileName2) 0002 % Builds Directed Enzyme-Enzyme Networks Removing Currency Metabolites (based-on a Library file) 0003 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0004 % The function reads a Metabolic Network SBML file and builds Directed Enzyme-Enzyme Networks. 0005 % The Remove Currency Metabolites based-on Library (RCMLib) algorithm removes currency metabolites 0006 % in the metabolic network automatically IF AND ONLY IF the currency metabolits exists in the Library file. 0007 % Note: COBRA Toolbox must be installed in MATLAB before running this function 0008 % 0009 % [Output] = enz_cent_RCMLib_dir(fileName1,fileName2) 0010 % 0011 %INPUTS 0012 % fileName1 The Library file includes pre-defined currency metabolites (in .txt format) 0013 % Note: Library text file must include one metabolites per line (all in one column) 0014 % fileName2 The metabolic Network in the SBML format 0015 % 0016 %OUTPUTS 0017 % *_Enzyme_Cent_RCMLib_Dir_Cyt.sif Directed-Enzyme-Enzyme Network - Cytoscape Compatible 0018 % 0019 % Yazdan Asgari 07/16/2016 http://yazdan59.ir/scan 0020 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0021 0022 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0023 % check validity of input files format 0024 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0025 check1=regexp(fileName1,'.txt'); 0026 assert(~isempty(check1),'Error in the first input: The fileName1 must contain .txt at its end') 0027 check2=regexp(fileName2,'.xml'); 0028 assert(~isempty(check2),'Error in the second input: The fileName2 must contain .xml at its end') 0029 0030 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0031 % start time evaluation of program 0032 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0033 tic; 0034 0035 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0036 % reading the Library text file and construct array of currency metabolites 0037 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0038 fid = fopen(fileName1); 0039 tline = fgetl(fid); 0040 i=1; 0041 Curr_met={}; 0042 while ischar(tline) 0043 Curr_met{i,1}=tline; 0044 tline = fgetl(fid); 0045 i=i+1; 0046 end 0047 fclose(fid); 0048 [h,g]=size(Curr_met); 0049 0050 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0051 % reading the SBML file using COBRA Toolbox Command, and sets size of the S matrix 0052 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0053 model=readCbModel(fileName2); 0054 [m,n]=size(model.S); 0055 0056 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0057 % calculate summation of each columns (i.e. How many metabolites each enzyme correlates) 0058 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0059 S_bin=zeros(size(model.S)); 0060 S_bin(find(model.S))=1; 0061 CB=sum(S_bin,1); 0062 A=zeros(m,n); 0063 B=zeros(m,1); 0064 N3=zeros(m,1); 0065 N_curr=zeros(m,1); 0066 0067 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0068 % reading the Metabolites array and check their availability in the library text file 0069 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0070 for q=1:m 0071 for i=1:h 0072 if strcmp(model.metNames{q},Curr_met{i,1})==1 0073 N_curr(q,1)=N_curr(q,1)+1; 0074 end 0075 end 0076 end 0077 0078 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0079 % for each binary S-matrix element, subtracts its value from the column summation and put the result in the A matrix. 0080 % A(q) means the metabolite q connects to how many other metabolites through the enzyme i. 0081 % for each row, sums the binary S-matrix over all columns. 0082 % B(q) means how many enzymes the metabolite q correlates. 0083 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0084 for q=1:m 0085 for i=1:n 0086 if S_bin(q,i)~=0 0087 A(q,i)=CB(1,i)-S_bin(q,i); 0088 end 0089 B(q,1)=B(q,1)+S_bin(q,i); 0090 end 0091 end 0092 0093 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0094 % Assumption: Generally, every metabolite is connected to the other one through a specific enzyme. 0095 % If a metabolite connects to more than one metabolite through an enzyme, this will be considered as a suspicious case. 0096 % Therefore, every A(q) value equal to 3 will be marked for further analysis. 0097 % In addition, the availability of the metabolite in the library file will be checked. 0098 % So, the metabolites which do not exist in the library file, will not select for further analysis. 0099 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0100 for q=1:m 0101 for i=1:n 0102 if A(q,i)==3 && N_curr(q,1)~=0 0103 N3(q,1)=N3(q,1)+1; 0104 end 0105 end 0106 end 0107 0108 s=0; 0109 for q=1:m 0110 if N3(q,1)~=0 0111 s=1; 0112 end 0113 end 0114 0115 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0116 % If there is any value for N3 array, the RCM algorithm will be done. 0117 % This algorithm will be deleted the most probable metabolite among all (i.e. the one with the maximum value of N3 and C) 0118 % The selected metabolite will be deleted in the binary S-Matrix, and the "WHILE LOOP" repeated. 0119 % The algorithm is ended if there is not any suspicious case. 0120 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0121 while s==1 0122 C=zeros(m,1); 0123 max1=max(N3,[],1); 0124 for q=1:m 0125 if N3(q,1)==max1 0126 C(q,1)=B(q,1); 0127 else 0128 C(q,1)=0; 0129 end 0130 end 0131 max2=max(C,[],1); 0132 for q=1:m 0133 if ( (N3(q,1)==max1) && (C(q,1)==max2) ) 0134 for i=1:n 0135 S_bin(q,i)=0; 0136 model.S(q,i)=0; 0137 end 0138 end 0139 end 0140 0141 CB=sum(S_bin,1); 0142 A=zeros(m,n); 0143 B=zeros(m,1); 0144 N3=zeros(m,1); 0145 for q=1:m 0146 for i=1:n 0147 if S_bin(q,i)~=0 0148 A(q,i)=CB(1,i)-S_bin(q,i); 0149 end 0150 B(q,1)=B(q,1)+S_bin(q,i); 0151 end 0152 end 0153 for q=1:m 0154 for i=1:n 0155 if A(q,i)==3 && N_curr(q,1)~=0 0156 N3(q,1)=N3(q,1)+1; 0157 end 0158 end 0159 end 0160 s=0; 0161 for q=1:m 0162 if N3(q,1)~=0 0163 s=1; 0164 end 0165 end 0166 end 0167 0168 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0169 % building the output file name (Cytoscape compatble file) 0170 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0171 outname=strrep(fileName2,'.xml','_Enzyme_Cent_RCMLib_Dir_Cyt.sif') 0172 fout = fopen(outname, 'w+'); 0173 0174 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0175 % construction of Directed-Enzyme-Enzyme Network based on the new binary S-matrix 0176 % finds non-zero elements of the new S-matrix (in order to make the algorithm faster), 0177 % parses through each row, and considers an edge for every unlike-signs. 0178 % It also consider Reversibility for Enzyme-Enzyme network. 0179 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0180 for j=1:m 0181 indices=find(model.S(j,:)); 0182 [a,b]=size(indices); 0183 r=0; 0184 if b~=0 0185 r=1; 0186 end 0187 while r<b 0188 i=1; 0189 while i<(b-r+1) 0190 if model.S(j,indices(1,r))>0 && model.S(j,indices(1,r+i))<0 0191 if model.rev(indices(r+i))==1 && model.rev(indices(r))==1 0192 fprintf(fout,'%s\t%s\t%s\n',model.rxns{indices(1,r)},'reaction-product',model.rxns{indices(1,r+i)}); 0193 fprintf(fout,'%s\t%s\t%s\n',model.rxns{indices(1,r+i)},'reaction-product',model.rxns{indices(1,r)}); 0194 else 0195 fprintf(fout,'%s\t%s\t%s\n',model.rxns{indices(1,r)},'reaction-product',model.rxns{indices(1,r+i)}); 0196 end 0197 elseif model.S(j,indices(1,r))<0 && model.S(j,indices(1,r+i))>0 0198 if model.rev(indices(r+i))==1 && model.rev(indices(r))==1 0199 fprintf(fout,'%s\t%s\t%s\n',model.rxns{indices(1,r+i)},'reaction-product',model.rxns{indices(1,r)}); 0200 fprintf(fout,'%s\t%s\t%s\n',model.rxns{indices(1,r)},'reaction-product',model.rxns{indices(1,r+i)}); 0201 else 0202 fprintf(fout,'%s\t%s\t%s\n',model.rxns{indices(1,r+i)},'reaction-product',model.rxns{indices(1,r)}); 0203 end 0204 end 0205 i=i+1; 0206 end 0207 r=r+1; 0208 end 0209 end 0210 fclose(fout); 0211 0212 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0213 % End of time evaluation of program 0214 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0215 toc;