0001 function sbmlModel = convertCobraToSBML(model,sbmlLevel,sbmlVersion,compSymbolList,compNameList,debug_function)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
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
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
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
0077 tmp_unit_definition.id = reaction_units;
0078
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
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
0100 the_compartments = {};
0101
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
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
0119
0120
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
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
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
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
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
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
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
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
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,'&','&');
0254 str = strrep(str,'<','<');
0255 str = strrep(str,'>','>');
0256 str = strrep(str,'"','"');