Skip to content

fix: outputs and metabolite merging in replaceMets #756

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
edkerk merged 2 commits into from
Sep 30, 2026
Merged

fix: outputs and metabolite merging in replaceMets #756

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
a2ba49c
fix: outputs and metabolite merging in replaceMets
edkerk Sep 30, 2026
96723f8
Merge branch 'develop' into fix/replaceMets-outputs-develop
edkerk Sep 30, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
82 changes: 35 additions & 47 deletions core/replaceMets.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@
% Output:
% model model structure with selected metabolites replaced
% removedRxns identifiers of duplicate reactions that were removed
% idxDuplRxns index of removedRxns in original model
% idxDuplRxns indices of removedRxns in the input model
%
% Note: This function is useful when the model contains both 'oxygen' and
% 'o2' as metabolite names. If 'oxygen' and 'o2' are identifiers instead,
Expand Down Expand Up @@ -59,10 +59,11 @@
error('The to-be-replaced metabolite cannot be found in the model.');
end

rxnsWithMet = find(model.S(metIdx,:));
[~, rxnsWithMet] = find(model.S(metIdx,:));
rxnsWithMet = unique(rxnsWithMet);
if verbose==true
fprintf('\n\nThe following reactions contain the to-be-replaced metabolite as reactant:\n')
fprintf(strjoin(model.rxns(rxnsWithMet),'\n'))
fprintf('\n\nThe following reactions involve the to-be-replaced metabolite:\n')
fprintf('%s\n',model.rxns{rxnsWithMet})
fprintf('\n')
end

Expand All @@ -88,8 +89,11 @@

idxDelete=[];
if identifiers
if any(ismember(metIdx,repIdx))
error('The to-be-replaced metabolite and the replacement metabolite are the same.');
end
originalStoch = model.S(metIdx,rxnsWithMet);
model.S(repIdx,rxnsWithMet) = originalStoch;
model.S(repIdx,rxnsWithMet) = model.S(repIdx,rxnsWithMet) + originalStoch;
model.S(metIdx,rxnsWithMet) = 0;
idxDelete = metIdx;
else
Expand All @@ -99,54 +103,38 @@
% the to-be-replace metabolite ID deleted.

% Build list of metaboliteName[compartment]
metCompsN =cellstr(num2str(model.metComps));
map = containers.Map(cellstr(num2str(transpose(1:length(model.comps)))),model.comps);
metCompsN = map.values(metCompsN);
metCompsN = strcat(lower(model.metNames),'[',metCompsN,']');
comps = model.comps(:);
metCompsN = strcat(lower(model.metNames(:)),'[',comps(model.metComps(:)),']');

for i = 1:length(repIdx)
metCompsNidx=find(strcmp(metCompsN(repIdx(i)), metCompsN));
if length(metCompsNidx)>1
for j = 2:length(metCompsNidx)
model.S(metCompsNidx(1),:) = model.S(metCompsNidx(1),:) + model.S(metCompsNidx(j),:);
idxDelete=[idxDelete; metCompsNidx(j)]; % Make list of metabolite IDs to delete
end
% Within each metaboliteName[compartment] group, the first replacement
% metabolite is kept and the other members are merged into it.
groups = unique(metCompsN(repIdx),'stable');
for i = 1:numel(groups)
metCompsNidx = find(strcmp(groups{i}, metCompsN));
keepIdx = metCompsNidx(find(ismember(metCompsNidx,repIdx),1));
others = metCompsNidx(metCompsNidx ~= keepIdx);
if ~isempty(others)
model.S(keepIdx,:) = model.S(keepIdx,:) + sum(model.S(others,:),1);
idxDelete = [idxDelete; others]; % Make list of metabolite IDs to delete
end
end
end

if ~isempty(idxDelete)
model.S(idxDelete,:) =[];
model.mets(idxDelete) = [];
model.metNames(idxDelete) = [];
model.metComps(idxDelete) = [];
model.b(idxDelete) = [];
if isfield(model,'metFormulas')
model.metFormulas(idxDelete) = [];
end
if isfield(model,'unconstrained')
model.unconstrained(idxDelete) = [];
end
if isfield(model,'metMiriams')
model.metMiriams(idxDelete) = [];
end
if isfield(model,'metCharges')
model.metCharges(idxDelete) = [];
end
if isfield(model,'metDeltaG')
model.metDeltaG(idxDelete) = [];
end
if isfield(model,'inchis')
model.inchis(idxDelete) = [];
end
if isfield(model,'metSmiles')
model.metSmiles(idxDelete) = [];
end
if isfield(model,'metFrom')
model.metFrom(idxDelete) = [];
end
repMets = model.mets(repIdx);
model = removeMets(model,idxDelete);
% Indices of the replacement metabolites shift when metabolites with a
% lower index are deleted.
[~, repIdx] = ismember(repMets,model.mets);
repIdx(repIdx==0) = [];
end

% This could now have created duplicate reactions. Contract model.
model=contractModel(model,[],repIdx);
if isempty(idxDelete)
removedRxns = cell(0,1);
idxDuplRxns = zeros(0,1);
else
originalRxns = model.rxns;
[model, removedRxns] = contractModel(model,[],repIdx);
[~, idxDuplRxns] = ismember(removedRxns,originalRxns);
end
end
2 changes: 1 addition & 1 deletion doc/core/removeMets.html
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -58,7 +58,7 @@ <h2><a name="_cross"></a>CROSS-REFERENCE INFORMATION <a href="#_top"><img alt="^
<li><a href="convertCharArray.html" class="code" title="function inputConverted = convertCharArray(funcInput)">convertCharArray</a> convertCharArray</li><li><a href="getIndexes.html" class="code" title="function indexes=getIndexes(model, objects, type, returnLogical)">getIndexes</a> getIndexes</li><li><a href="removeReactions.html" class="code" title="function reducedModel=removeReactions(model,rxnsToRemove,removeUnusedMets,removeUnusedGenes,removeUnusedComps)">removeReactions</a> removeReactions</li></ul>
This function is called by:
<ul style="list-style-image:url(../matlabicon.gif)">
<li><a href="mergeCompartments.html" class="code" title="function [model, deletedRxns, duplicateRxns]=mergeCompartments(model,keepUnconstrained,deleteRxnsWithOneMet,distReverse)">mergeCompartments</a> mergeCompartments</li><li><a href="removeReactions.html" class="code" title="function reducedModel=removeReactions(model,rxnsToRemove,removeUnusedMets,removeUnusedGenes,removeUnusedComps)">removeReactions</a> removeReactions</li><li><a href="simplifyModel.html" class="code" title="function [reducedModel, deletedReactions, deletedMetabolites]=simplifyModel(model,deleteUnconstrained, deleteDuplicates, deleteZeroInterval, deleteInaccessible, deleteMinMax, groupLinear, constrainReversible, reservedRxns, suppressWarnings)">simplifyModel</a> simplifyModel</li></ul>
<li><a href="mergeCompartments.html" class="code" title="function [model, deletedRxns, duplicateRxns]=mergeCompartments(model,keepUnconstrained,deleteRxnsWithOneMet,distReverse)">mergeCompartments</a> mergeCompartments</li><li><a href="removeReactions.html" class="code" title="function reducedModel=removeReactions(model,rxnsToRemove,removeUnusedMets,removeUnusedGenes,removeUnusedComps)">removeReactions</a> removeReactions</li><li><a href="replaceMets.html" class="code" title="function [model, removedRxns, idxDuplRxns]=replaceMets(model,metabolite,replacement,verbose,identifiers)">replaceMets</a> replaceMets</li><li><a href="simplifyModel.html" class="code" title="function [reducedModel, deletedReactions, deletedMetabolites]=simplifyModel(model,deleteUnconstrained, deleteDuplicates, deleteZeroInterval, deleteInaccessible, deleteMinMax, groupLinear, constrainReversible, reservedRxns, suppressWarnings)">simplifyModel</a> simplifyModel</li></ul>
<!-- crossreference -->


Expand Down
Loading
Loading