convertCobraToSBML

PURPOSE ^

convertCobraToSBML converts a cobra structure to an sbml

SYNOPSIS ^

function sbmlModel = convertCobraToSBML(model,sbmlLevel,sbmlVersion,compSymbolList,compNameList,debug_function)

DESCRIPTION ^

convertCobraToSBML converts a cobra structure to an sbml
structure using the structures provided in the SBML toolbox 3.1.0

 sbmlModel = convertCobraToSBML(model,sbmlLevel,sbmlVersion,compSymbolList,compNameList)

INPUTS
 model             COBRA model structure

OPTIONAL INPUTS
 sbmlLevel         SBML Level (default = 2)
 sbmlVersion       SBML Version (default = 1)
 compSymbolList    List of compartment symbols
 compNameList      List of copmartment names correspoding to compSymbolList

OUTPUT
 sbmlModel         SBML MATLAB structure

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function sbmlModel = convertCobraToSBML(model,sbmlLevel,sbmlVersion,compSymbolList,compNameList,debug_function)
0002 %convertCobraToSBML converts a cobra structure to an sbml
0003 %structure using the structures provided in the SBML toolbox 3.1.0
0004 %
0005 % sbmlModel = convertCobraToSBML(model,sbmlLevel,sbmlVersion,compSymbolList,compNameList)
0006 %
0007 %INPUTS
0008 % model             COBRA model structure
0009 %
0010 %OPTIONAL INPUTS
0011 % sbmlLevel         SBML Level (default = 2)
0012 % sbmlVersion       SBML Version (default = 1)
0013 % compSymbolList    List of compartment symbols
0014 % compNameList      List of copmartment names correspoding to compSymbolList
0015 %
0016 %OUTPUT
0017 % sbmlModel         SBML MATLAB structure
0018 %
0019 %
0020 
0021 %NOTE: The name mangling of reaction and metabolite ids is necessary
0022 %for compliance with the SBML sID standard.
0023 %
0024 %NOTE: Sometimes the Model_create function doesn't listen to the
0025 %sbmlVersion parameter, so it is essential that the items that
0026 %are added to the sbmlModel are defined with the sbmlModel's level
0027 %and version:  sbmlModel.SBML_level,sbmlModel.SBML_version
0028 %
0029 %NOTE:  Some of the structures are recycled to reduce to overhead for
0030 %their creation.  There's a chance this can cause bugs in the future.
0031 %
0032 %NOTE: Currently, I don't add in the boundary metabolites.
0033 %
0034 %NOTE: Speed could probably be improved by directly adding structures to
0035 %lists in a struct instead of using the SBML _addItem function, but this
0036 %could break in future versions of the SBML toolbox.
0037 %
0038 %POTENTIAL FUTURE BUG: To speed things up, sbml structs have been
0039 %recycled and are directly appended into lists instead of using _addItem
0040 if (~exist('sbmlLevel','var') || isempty(sbmlLevel))
0041     sbmlLevel = 2;
0042 end
0043 if (~exist('sbmlVersion','var') || isempty(sbmlVersion))
0044     sbmlVersion = 1;
0045 end
0046 if (~exist('debug_function','var') || isempty(debug_function))
0047   debug_function = 0;
0048 end
0049 reaction_units = 'mmol_per_gDW_per_hr';
0050 sbmlModel = Model_create(sbmlLevel, sbmlVersion);
0051 
0052 sbmlModel.namespaces = struct();
0053 sbmlModel.namespaces.prefix = '';
0054 sbmlModel.namespaces.uri = 'http://www.sbml.org/sbml/level2';
0055 if isfield(model,'description')
0056     sbmlModel.id = strrep(strrep(strrep(model.description,'.','_'), filesep, '_'), ':','_');
0057 else
0058     sbmlModel.id = '';
0059 end
0060 %POTENTIAL FUTURE BUG: Create temporary structs to speed things up.
0061 tmp_unit = Unit_create(sbmlModel.SBML_level, sbmlModel.SBML_version);
0062 tmp_species = Species_create(sbmlModel.SBML_level, sbmlModel.SBML_version);
0063 sbml_tmp_compartment = Compartment_create(sbmlModel.SBML_level, sbmlModel.SBML_version);
0064 sbml_tmp_parameter = Parameter_create(sbmlModel.SBML_level, sbmlModel.SBML_version);
0065 sbml_tmp_species_ref = SpeciesReference_create(sbmlModel.SBML_level, sbmlModel.SBML_version);
0066 sbml_tmp_reaction = Reaction_create(sbmlModel.SBML_level, sbmlModel.SBML_version);
0067 sbml_tmp_law = KineticLaw_create(sbmlModel.SBML_level, sbmlModel.SBML_version);
0068 tmp_unit_definition = UnitDefinition_create(sbmlModel.SBML_level, sbmlModel.SBML_version);
0069 
0070 %% Compartments
0071 if ~exist('compSymbolList','var') || isempty(compSymbolList)
0072     compSymbolList = {'c','m','v','x','e','t','g','r','n','p','l'};
0073     compNameList = {'Cytoplasm','Mitochondrion','Vacuole','Peroxisome','Extracellular','Pool','Golgi','Endoplasmic_reticulum','Nucleus','Periplasm','Lysosome'};
0074 end
0075 
0076 %Create and add the unit definition to the sbml model struct.
0077 tmp_unit_definition.id =  reaction_units;
0078 %The 4 following lists are in matched order for each unit.
0079 unit_kinds = {'mole','gram','second'};
0080 unit_exponents = [1 -1 -1];
0081 unit_scales = [-3 0 0];
0082 unit_multipliers = [1 1 1.0/60/60];
0083 %Add the units to the unit definition
0084 for i = 1:size(unit_kinds, 2)
0085     tmp_unit.kind = unit_kinds{ i };
0086     tmp_unit.exponent = unit_exponents(i);
0087     tmp_unit.scale = unit_scales(i);
0088     tmp_unit.multiplier = unit_multipliers(i);
0089     tmp_unit_definition = UnitDefinition_addUnit(tmp_unit_definition, tmp_unit);
0090 end
0091 if debug_function
0092   if ~isSBML_Unit(tmp_unit_definition)
0093     error('unit definition failed')
0094   end
0095 end
0096 sbmlModel = Model_addUnitDefinition(sbmlModel, tmp_unit_definition);
0097 
0098 
0099 %List to hold the compartment ids.
0100 the_compartments = {};
0101 %separate metabolite and compartment
0102 [tokens tmp_met_struct] = regexp(model.mets,'(?<met>.+)\[(?<comp>.+)\]|(?<met>.+)\((?<comp>.+)\)','tokens','names');
0103 for (i=1:size(model.mets, 1))
0104     tmp_notes='';
0105     tmp_met = tmp_met_struct{i}.met;
0106     %Change id to correspond to SBML id specifications
0107     tmp_met = strcat('M_', (tmp_met), '_', tmp_met_struct{i}.comp);
0108     model.mets{ i } = formatForSBMLID(tmp_met);
0109     tmp_species.id = formatForSBMLID(tmp_met);
0110     tmp_species.compartment = formatForSBMLID(tmp_met_struct{i}.comp);
0111     if isfield(model, 'metNames')
0112         tmp_species.name = (model.metNames{i});
0113     end
0114     if isfield(model, 'metFormulas')
0115         tmp_notes = [tmp_notes '<p>FORMULA: ' model.metFormulas{i} '</p>'];
0116     end
0117     if isfield(model, 'metCharge')
0118         %NOTE: charge is being removed in SBML level 3
0119 %         tmp_species.charge = model.metCharge(i);
0120 %         tmp_species.isSetCharge = 1;
0121         tmp_notes = [tmp_notes '<p>CHARGE: ' num2str(model.metCharge(i)) '</p>'];
0122     end
0123 
0124 
0125     if ~isempty(tmp_notes)
0126         tmp_species.notes = ['<body xmlns="http://www.w3.org/1999/xhtml">' tmp_notes '</body>'];
0127     end
0128     sbmlModel.species = [ sbmlModel.species tmp_species ];
0129     %This is where the compartment symbols are aggregated.
0130     the_compartments{ i } = tmp_species.compartment ;
0131 end
0132 if debug_function
0133   for (i = 1:size(sbmlModel.species, 2))
0134     if ~isSBML_Species(sbmlModel.species(i), sbmlLevel, sbmlVersion)
0135       error('SBML species failed to pass test')
0136     end
0137   end
0138 end
0139 %Add the unique compartments to the model struct.
0140 the_compartments = unique(the_compartments);
0141 for (i=1:size(the_compartments,2))
0142     tmp_id = the_compartments{1,i};
0143     tmp_symbol_index = find(strcmp(formatForSBMLID(compSymbolList),tmp_id));
0144     %Check that symbol is in compSymbolList
0145     if ~isempty(tmp_symbol_index)
0146         tmp_name = compNameList{tmp_symbol_index};
0147     else
0148         error(['Unknown compartment: ' tmp_id '. Be sure that ' tmp_id ' is specified in compSymbolList and compNameList.'])
0149     end
0150     tmp_id = formatForSBMLID(tmp_id);
0151     sbml_tmp_compartment.id = tmp_id;
0152     sbml_tmp_compartment.name = tmp_name;
0153     sbmlModel = Model_addCompartment(sbmlModel, sbml_tmp_compartment);
0154 end
0155 if debug_function
0156   for (i = 1:size(sbmlModel.compartment, 2))
0157     if ~isSBML_Compartment(sbmlModel.compartment(i), sbmlLevel, sbmlVersion)
0158       error('SBML compartment failed to pass test')
0159     end
0160   end
0161 end
0162 %Add the reactions to the model struct.  Use the species references.
0163 sbml_tmp_parameter.units = reaction_units;
0164 sbml_tmp_parameter.isSetValue = 1;
0165 for (i=1:size(model.rxns, 1))
0166     tmp_id =  strcat('R_', formatForSBMLID(model.rxns{i}));
0167     model.rxns{i} = tmp_id;
0168     met_idx = find(model.S(:, i));
0169     sbml_tmp_reaction.notes = '';
0170     %Reset the fields that have been filled.
0171     sbml_tmp_reaction.reactant = [];
0172     sbml_tmp_reaction.product = [];
0173     sbml_tmp_reaction.kineticLaw = [];
0174     sbml_tmp_reaction.id = tmp_id;
0175     if isfield(model, 'rxnNames')
0176         sbml_tmp_reaction.name = model.rxnNames{i};
0177     end
0178     if isfield(model, 'rev')
0179         sbml_tmp_reaction.reversible = model.rev(i);
0180     end
0181     sbml_tmp_law.parameter = [];
0182     sbml_tmp_law.formula = 'FLUX_VALUE';
0183     sbml_tmp_parameter.id = 'LOWER_BOUND';
0184     sbml_tmp_parameter.value = model.lb(i);
0185     sbml_tmp_law.parameter = [ sbml_tmp_law.parameter sbml_tmp_parameter ];
0186     sbml_tmp_parameter.id = 'UPPER_BOUND';
0187     sbml_tmp_parameter.value = model.ub(i);
0188     sbml_tmp_law.parameter = [ sbml_tmp_law.parameter sbml_tmp_parameter ];
0189     sbml_tmp_parameter.id = 'FLUX_VALUE';
0190     sbml_tmp_parameter.value = 0;
0191     sbml_tmp_law.parameter = [ sbml_tmp_law.parameter sbml_tmp_parameter ];
0192     sbml_tmp_parameter.id = 'OBJECTIVE_COEFFICIENT';
0193     sbml_tmp_parameter.value = model.c(i);
0194     sbml_tmp_law.parameter = [ sbml_tmp_law.parameter sbml_tmp_parameter ];
0195     sbml_tmp_reaction.kineticLaw = sbml_tmp_law;
0196     %Add in other notes
0197     tmp_note = '';
0198     if isfield(model, 'grRules')
0199         tmp_note = [tmp_note '<p>GENE_ASSOCIATION: ' model.grRules{i} '</p>' ];
0200     end
0201     if isfield(model, 'subSystems')
0202         tmp_note = [ tmp_note ' <p>SUBSYSTEM: ' model.subSystems{i} '</p>'];
0203     end
0204     if isfield(model, 'rxnECNumbers')
0205         tmp_note = [ tmp_note ' <p>EC Number: ' model.rxnECNumbers{i} '</p>'];
0206     end
0207     if isfield(model, 'confidenceScores')
0208         tmp_note = [ tmp_note ' <p>Confidence Level: ' model.confidenceScores{i} '</p>'];
0209     end
0210     if isfield(model, 'rxnReferences')
0211         tmp_note = [ tmp_note ' <p>AUTHORS: ' model.rxnReferences{i} '</p>'];
0212     end
0213     if isfield(model, 'rxnNotes')
0214         tmp_note = [ tmp_note ' <p>' model.rxnNotes{i} '</p>'];
0215     end
0216     if ~isempty(tmp_note)
0217         sbml_tmp_reaction.notes = ['<body xmlns="http://www.w3.org/1999/xhtml">' tmp_note '</body>'];
0218     end
0219     %Add in the reactants and products
0220     for (j_met=1:size(met_idx,1))
0221         tmp_idx = met_idx(j_met,1);
0222         sbml_tmp_species_ref.species = model.mets{tmp_idx};
0223         met_stoich = model.S(tmp_idx, i);
0224         sbml_tmp_species_ref.stoichiometry = abs(met_stoich);
0225         if (met_stoich > 0)
0226             sbml_tmp_reaction.product = [ sbml_tmp_reaction.product sbml_tmp_species_ref ];
0227         else
0228             sbml_tmp_reaction.reactant = [ sbml_tmp_reaction.reactant sbml_tmp_species_ref];
0229         end
0230     end
0231     sbmlModel.reaction = [ sbmlModel.reaction sbml_tmp_reaction ];
0232 end
0233 if debug_function
0234   for (i = 1:size(sbmlModel.reaction, 2))
0235     if ~isSBML_Reaction(sbmlModel.reaction(i), sbmlLevel, sbmlVersion)
0236       error('SBML reaction failed to pass test')
0237     end
0238   end
0239 end
0240 %% Format For SBML
0241 function str = formatForSBMLID(str)
0242 str = strrep(str,'-','_DASH_');
0243 str = strrep(str,'/','_FSLASH_');
0244 str = strrep(str,'\','_BSLASH_');
0245 str = strrep(str,'(','_LPAREN_');
0246 str = strrep(str,')','_RPAREN_');
0247 str = strrep(str,'[','_LSQBKT_');
0248 str = strrep(str,']','_RSQBKT_');
0249 str = strrep(str,',','_COMMA_');
0250 str = strrep(str,'.','_PERIOD_');
0251 str = strrep(str,'''','_APOS_');
0252 str = regexprep(str,'\(e\)$','_e');
0253 str = strrep(str,'&','&amp;');
0254 str = strrep(str,'<','&lt;');
0255 str = strrep(str,'>','&gt;');
0256 str = strrep(str,'"','&quot;');

Generated on Thu 21-Jun-2012 15:39:23 by m2html © 2003