addMets

 addMets
   Adds metabolites to a model

 Input:
   model       a model structure
   metsToAdd   the metabolite structure can have the following fields:
       mets            cell array with unique strings that identifies each
                       metabolite (optional, default IDs of new
                       metabolites are numbered with the prefix defined
                       below)
       metNames        cell array with the names of each metabolite
       compartments    cell array with the compartment of each
                       metabolite. Should match model.comps. If this is a
                       string rather than a cell array it is assumed that
                       all mets are in that compartment
       b               Nx1 or Nx2 matrix with equality constraints for
                       each metabolite (optional, default 0)
       unconstrained   vector describing if each metabolite is an exchange
                       metabolite (1) or not (0) (optional, default 0)
       inchis          cell array with InChI strings (optional, default '')
       metSmiles       cell array with SMILES strings (optional, default '')
       metFormulas     cell array with the formulas  (optional, default '')
       metMiriams      cell array with MIRIAM structures (optional, default [])
       metCharges      metabolite charge (optional, default NaN)
       metDeltaG       Gibbs free energy of formation at biochemical
                       standard condition in kJ/mole (optional, default NaN)
       metNotes        cell array with metabolite notes as strings
                       (optional, default '')
   copyInfo    when adding metabolites to a compartment where it
               previously did not exist, the function will copy any
               available annotation from the metabolite in another
               compartment (optional, default true)
   prefix      when metsToAdd.mets is not specified, new metabolite IDs
               are generated with the prefix specified here. If IDs with
               the prefix are already used in the model then the
               numbering will start from the highest existing integer+1
               (optional, default 'm_')

 Output:
   newModel     an updated model structure

 This function does not make extensive checks about MIRIAM formats,
 forbidden characters or such.

 If multiple metabolites are added at once, the metMiriams cell array
 should be defined as (example with ChEBI and KEGG):

 metsToAdd.metMiriams{1} = struct('name',{{'chebi';'kegg.compound'}},...
     'value',{{'CHEBI:18072';'C11821'}});
 metsToAdd.metMiriams{2} = struct('name',{{'chebi';'kegg.compound'}},...
     'value',{{'CHEBI:31132';'C12248'}});

 Usage: newModel = addMets(model, metsToAdd, copyInfo, prefix)

0001 function newModel=addMets(model,metsToAdd,copyInfo,prefix)
0002 % addMets
0003 %   Adds metabolites to a model
0004 %
0005 % Input:
0006 %   model       a model structure
0007 %   metsToAdd   the metabolite structure can have the following fields:
0008 %       mets            cell array with unique strings that identifies each
0009 %                       metabolite (optional, default IDs of new
0010 %                       metabolites are numbered with the prefix defined
0011 %                       below)
0012 %       metNames        cell array with the names of each metabolite
0013 %       compartments    cell array with the compartment of each
0014 %                       metabolite. Should match model.comps. If this is a
0015 %                       string rather than a cell array it is assumed that
0016 %                       all mets are in that compartment
0017 %       b               Nx1 or Nx2 matrix with equality constraints for
0018 %                       each metabolite (optional, default 0)
0019 %       unconstrained   vector describing if each metabolite is an exchange
0020 %                       metabolite (1) or not (0) (optional, default 0)
0021 %       inchis          cell array with InChI strings (optional, default '')
0022 %       metSmiles       cell array with SMILES strings (optional, default '')
0023 %       metFormulas     cell array with the formulas  (optional, default '')
0024 %       metMiriams      cell array with MIRIAM structures (optional, default [])
0025 %       metCharges      metabolite charge (optional, default NaN)
0026 %       metDeltaG       Gibbs free energy of formation at biochemical
0027 %                       standard condition in kJ/mole (optional, default NaN)
0028 %       metNotes        cell array with metabolite notes as strings
0029 %                       (optional, default '')
0030 %   copyInfo    when adding metabolites to a compartment where it
0031 %               previously did not exist, the function will copy any
0032 %               available annotation from the metabolite in another
0033 %               compartment (optional, default true)
0034 %   prefix      when metsToAdd.mets is not specified, new metabolite IDs
0035 %               are generated with the prefix specified here. If IDs with
0036 %               the prefix are already used in the model then the
0037 %               numbering will start from the highest existing integer+1
0038 %               (optional, default 'm_')
0039 %
0040 % Output:
0041 %   newModel     an updated model structure
0042 %
0043 % This function does not make extensive checks about MIRIAM formats,
0044 % forbidden characters or such.
0045 %
0046 % If multiple metabolites are added at once, the metMiriams cell array
0047 % should be defined as (example with ChEBI and KEGG):
0048 %
0049 % metsToAdd.metMiriams{1} = struct('name',{{'chebi';'kegg.compound'}},...
0050 %     'value',{{'CHEBI:18072';'C11821'}});
0051 % metsToAdd.metMiriams{2} = struct('name',{{'chebi';'kegg.compound'}},...
0052 %     'value',{{'CHEBI:31132';'C12248'}});
0053 %
0054 % Usage: newModel = addMets(model, metsToAdd, copyInfo, prefix)
0055 
0056 if nargin<3
0057     copyInfo=true;
0058 end
0059 if nargin<4
0060     prefix='m_';
0061 else
0062     prefix=char(prefix);
0063 end
0064 
0065 newModel=model;
0066 
0067 if isempty(metsToAdd)
0068     return;
0069 end
0070 
0071 %Check some stuff regarding the required fields
0072 if ~isfield(metsToAdd,'mets')
0073     metsToAdd.mets=generateNewIds(newModel,'mets',prefix,numel(metsToAdd.metNames));
0074 else
0075     metsToAdd.mets=convertCharArray(metsToAdd.mets);
0076 end
0077 if ~isfield(metsToAdd,'metNames')
0078     metsToAdd.metNames=metsToAdd.mets;
0079 else
0080     metsToAdd.metNames=convertCharArray(metsToAdd.metNames);
0081 end
0082 if ~isfield(metsToAdd,'compartments')
0083     EM='compartments is a required field in metsToAdd';
0084     dispEM(EM);
0085 else
0086     metsToAdd.compartments=convertCharArray(metsToAdd.compartments);
0087     %If only one compartment is given, assume it is for all metabolites
0088     if numel(metsToAdd.compartments)==1 && numel(metsToAdd.mets)>1
0089         temp=cell(numel(metsToAdd.mets),1);
0090         temp(:)=metsToAdd.compartments;
0091         metsToAdd.compartments=temp;
0092     end
0093 end
0094 
0095 %Number of metabolites
0096 nMets=numel(metsToAdd.mets);
0097 nOldMets=numel(model.mets);
0098 filler=cell(nMets,1);
0099 filler(:)={''};
0100 largeFiller=cell(nOldMets,1);
0101 largeFiller(:)={''};
0102 
0103 %Check that no metabolite ids are already present in the model
0104 I=ismember(metsToAdd.mets,model.mets);
0105 if any(I)
0106     error('One or more elements in metsToAdd.mets are already present in model.mets: %s',...
0107         strjoin(metsToAdd.mets(I),', '));
0108 else
0109     newModel.mets=[newModel.mets;metsToAdd.mets(:)];
0110 end
0111 
0112 %Check that all the compartments could be found
0113 [I, compMap]=ismember(metsToAdd.compartments,model.comps);
0114 if ~all(I)
0115     EM='metsToAdd.compartments must match model.comps';
0116     dispEM(EM);
0117 end
0118 
0119 %Check that the metabolite names aren't present in the same compartment.
0120 %Not the neatest way maybe..
0121 t1=strcat(metsToAdd.metNames(:),'***',metsToAdd.compartments(:));
0122 t2=strcat(model.metNames,'***',model.comps(model.metComps));
0123 if any(ismember(t1,t2))
0124     EM='One or more elements in metsToAdd.metNames already exist in the same compartments as the one it is being added to';
0125     dispEM(EM);
0126 end
0127 
0128 %Some more checks and if they pass then add each field to the structure
0129 if numel(metsToAdd.metNames)~=nMets
0130     EM='metsToAdd.metNames must have the same number of elements as metsToAdd.mets';
0131     dispEM(EM);
0132 else
0133     newModel.metNames=[newModel.metNames;metsToAdd.metNames(:)];
0134 end
0135 
0136 if numel(compMap)~=nMets
0137     EM='metsToAdd.compartments must have the same number of elements as metsToAdd.mets';
0138     dispEM(EM);
0139 else
0140     newModel.metComps=[newModel.metComps;compMap(:)];
0141 end
0142 
0143 if isfield(metsToAdd,'b')
0144     if size(metsToAdd.b,1)~=nMets
0145         EM='metsToAdd.b must have the same number of elements as metsToAdd.mets';
0146         dispEM(EM);
0147     else
0148         %Add empty field if it doesn't exist
0149         if ~isfield(newModel,'b')
0150             newModel.b=zeros(nOldMets,1);
0151         end
0152         
0153         %If the original is only one vector
0154         if size(metsToAdd.b,2)>size(newModel.b,2)
0155             newModel.b=[newModel.b newModel.b];
0156         end
0157         %Add the new ones
0158         newModel.b=[newModel.b;metsToAdd.b];
0159     end
0160 else
0161     if isfield(newModel,'b')
0162         %Add the default
0163         newModel.b=[newModel.b;zeros(nMets,size(newModel.b,2))];
0164     end
0165 end
0166 
0167 if isfield(metsToAdd,'unconstrained')
0168     if numel(metsToAdd.unconstrained)~=nMets
0169         EM='metsToAdd.unconstrained must have the same number of elements as metsToAdd.mets';
0170         dispEM(EM);
0171     else
0172         %Add empty field if it doesn't exist
0173         if ~isfield(newModel,'unconstrained')
0174             newModel.unconstrained=zeros(nOldMets,1);
0175         end
0176         
0177         %Add the new ones
0178         newModel.unconstrained=[newModel.unconstrained;metsToAdd.unconstrained(:)];
0179     end
0180 else
0181     if isfield(newModel,'unconstrained')
0182         %Add the default
0183         newModel.unconstrained=[newModel.unconstrained;zeros(nMets,1)];
0184     end
0185 end
0186 
0187 if isfield(metsToAdd,'inchis')
0188     metsToAdd.inchis=convertCharArray(metsToAdd.inchis);
0189     if numel(metsToAdd.inchis)~=nMets
0190         EM='metsToAdd.inchis must have the same number of elements as metsToAdd.mets';
0191         dispEM(EM);
0192     end
0193     %Add empty field if it doesn't exist
0194     if ~isfield(newModel,'inchis')
0195         newModel.inchis=largeFiller;
0196     end
0197     newModel.inchis=[newModel.inchis;metsToAdd.inchis(:)];
0198 else
0199     %Add empty strings if structure is in model
0200     if isfield(newModel,'inchis')
0201         newModel.inchis=[newModel.inchis;filler];
0202     end
0203 end
0204 
0205 
0206 if isfield(metsToAdd,'metSmiles')
0207     metsToAdd.metSmiles=convertCharArray(metsToAdd.metSmiles);
0208     if numel(metsToAdd.metSmiles)~=nMets
0209         EM='metsToAdd.metSmiles must have the same number of elements as metsToAdd.mets';
0210         dispEM(EM);
0211     end
0212     %Add empty field if it doesn't exist
0213     if ~isfield(newModel,'metSmiles')
0214         newModel.metSmiles=largeFiller;
0215     end
0216     newModel.metSmiles=[newModel.metSmiles;metsToAdd.metSmiles(:)];
0217 else
0218     %Add empty strings if structure is in model
0219     if isfield(newModel,'metSmiles')
0220         newModel.metSmiles=[newModel.metSmiles;filler];
0221     end
0222 end
0223 
0224 if isfield(metsToAdd,'metFormulas')
0225     metsToAdd.metFormulas=convertCharArray(metsToAdd.metFormulas);
0226     if numel(metsToAdd.metFormulas)~=nMets
0227         EM='metsToAdd.metFormulas must have the same number of elements as metsToAdd.mets';
0228         dispEM(EM);
0229     end
0230     %Add empty field if it doesn't exist
0231     if ~isfield(newModel,'metFormulas')
0232         newModel.metFormulas=largeFiller;
0233     end
0234     newModel.metFormulas=[newModel.metFormulas;metsToAdd.metFormulas(:)];
0235 else
0236     %Add default
0237     if isfield(newModel,'metFormulas')
0238         newModel.metFormulas=[newModel.metFormulas;filler];
0239     end
0240 end
0241 
0242 if isfield(metsToAdd,'metCharges')
0243     if numel(metsToAdd.metCharges)~=nMets
0244         EM='metsToAdd.metCharges must have the same number of elements as metsToAdd.mets';
0245         dispEM(EM);
0246     end
0247     if ~isnumeric(metsToAdd.metCharges)
0248         EM='metsToAdd.metCharges must be of type "double"';
0249         dispEM(EM);
0250     end
0251     if ~isfield(newModel,'metCharges')
0252         newModel.metCharges=NaN(numel(largeFiller),1);
0253     end
0254     newModel.metCharges=[newModel.metCharges;metsToAdd.metCharges(:)];
0255 else
0256     %Add default
0257     if isfield(newModel,'metCharges')
0258         newModel.metCharges=[newModel.metCharges;NaN(numel(filler),1)];
0259     end
0260 end
0261 
0262 if isfield(metsToAdd,'metDeltaG')
0263     if numel(metsToAdd.metDeltaG)~=nMets
0264         EM='metsToAdd.metDeltaG must have the same number of elements as metsToAdd.mets';
0265         dispEM(EM);
0266     end
0267     if ~isnumeric(metsToAdd.metDeltaG)
0268         EM='metsToAdd.metDeltaG must be of type "double"';
0269         dispEM(EM);
0270     end
0271     if ~isfield(newModel,'metDeltaG')
0272         newModel.metDeltaG=NaN(numel(largeFiller),1);
0273     end
0274     newModel.metDeltaG=[newModel.metDeltaG;metsToAdd.metDeltaG(:)];
0275 else
0276     %Add default
0277     if isfield(newModel,'metDeltaG')
0278         newModel.metDeltaG=[newModel.metDeltaG;NaN(numel(filler),1)];
0279     end
0280 end
0281 
0282 
0283 if isfield(metsToAdd,'metNotes')
0284     metsToAdd.metNotes=convertCharArray(metsToAdd.metNotes);
0285     if numel(metsToAdd.metNotes)==1 && numel(metsToAdd.mets)>1
0286         temp=cell(numel(metsToAdd.mets),1);
0287         temp(:)=metsToAdd.metNotes;
0288         metsToAdd.metNotes=temp;
0289     end
0290     if numel(metsToAdd.metNotes)~=nMets
0291         EM='metsToAdd.metNotes must have the same number of elements as metsToAdd.mets';
0292         dispEM(EM);
0293     end
0294     %Add empty field if it doesn't exist
0295     if ~isfield(newModel,'metNotes')
0296         newModel.metNotes=largeFiller;
0297     end
0298     newModel.metNotes=[newModel.metNotes;metsToAdd.metNotes(:)];
0299 else
0300     %Add empty strings if structure is in model
0301     if isfield(newModel,'metNotes')
0302         newModel.metNotes=[newModel.metNotes;filler];
0303     end
0304 end
0305 
0306 %Don't check the type of metMiriams
0307 if isfield(metsToAdd,'metMiriams')
0308     if numel(metsToAdd.metMiriams)==1 && numel(metsToAdd.mets)>1
0309         temp=cell(numel(metsToAdd.mets),1);
0310         temp(:)={metsToAdd.metMiriams};
0311         metsToAdd.metMiriams=temp;
0312     end
0313     if numel(metsToAdd.metMiriams)~=nMets
0314         EM='metsToAdd.metMiriams must have the same number of elements as metsToAdd.mets';
0315         dispEM(EM);
0316     end
0317     %Add empty field if it doesn't exist
0318     if ~isfield(newModel,'metMiriams')
0319         newModel.metMiriams=cell(nOldMets,1);
0320     end
0321     newModel.metMiriams=[newModel.metMiriams;metsToAdd.metMiriams(:)];
0322 else
0323     if isfield(newModel,'metMiriams')
0324         newModel.metMiriams=[newModel.metMiriams;cell(nMets,1)];
0325     end
0326 end
0327 
0328 if isfield(newModel,'metFrom')
0329     newModel.metFrom=[newModel.metFrom;filler];
0330 end
0331 
0332 %Expand the S matrix
0333 newModel.S=[newModel.S;sparse(nMets,size(newModel.S,2))];
0334 
0335 if copyInfo==true
0336     [I, J]=ismember(metsToAdd.metNames,model.metNames);
0337     J=J(I);
0338     %I is the indexes of the new metabolites for which a metabolite with
0339     %the same name existed
0340     I=find(I)+nOldMets;
0341     %Go through each of the added mets and copy annotation if it doesn't
0342     %exist
0343     for i=1:numel(I)
0344         if isfield(newModel,'inchis')
0345             if isempty(newModel.inchis{I(i)})
0346                 newModel.inchis(I(i))=newModel.inchis(J(i));
0347             end
0348         end
0349         if isfield(newModel,'metSmiles')
0350             if isempty(newModel.metSmiles{I(i)})
0351                 newModel.metSmiles(I(i))=newModel.metSmiles(J(i));
0352             end
0353         end        
0354         if isfield(newModel,'metFormulas')
0355             if isempty(newModel.metFormulas{I(i)})
0356                 newModel.metFormulas(I(i))=newModel.metFormulas(J(i));
0357             end
0358         end
0359         if isfield(newModel,'metMiriams')
0360             if isempty(newModel.metMiriams{I(i)})
0361                 newModel.metMiriams(I(i))=newModel.metMiriams(J(i));
0362             end
0363         end
0364         if isfield(newModel,'metCharges')
0365             newModel.metCharges(I(i))=newModel.metCharges(J(i));
0366         end
0367         if isfield(newModel,'metDeltaG')
0368             newModel.metDeltaG(I(i))=newModel.metDeltaG(J(i));
0369         end
0370     end
0371 end
0372 end