diff --git a/curation/curateModelFromTables.m b/curation/curateModelFromTables.m index a7fd2dc0..efd628a6 100644 --- a/curation/curateModelFromTables.m +++ b/curation/curateModelFromTables.m @@ -258,8 +258,15 @@ rxnsToAdd.stoichCoeffs{i}=cell2mat(rxnsToAdd.stoichCoeffs{i}); %Check if the reaction not already exists, by stoichiometry of its %products and reactants - modelCoeffs = transpose(newModel.S(getIndexes(newModel,rxnsToAdd.mets{i},'mets'),:)); + metIdx = getIndexes(newModel,rxnsToAdd.mets{i},'mets'); + modelCoeffs = transpose(newModel.S(metIdx,:)); modelCoeffs2 = find(~all(modelCoeffs==0,2)); + % Only reactions built from exactly these metabolites can be + % duplicates. Any model reaction that also involves other + % metabolites projects onto the same row of coefficients, but is a + % different reaction and must not be overwritten. + sameSupport = full(sum(newModel.S(:,modelCoeffs2)~=0,1))' == numel(metIdx); + modelCoeffs2 = modelCoeffs2(sameSupport); modelCoeffs = modelCoeffs(modelCoeffs2,:); [~, duplicateRxn] = intersect(modelCoeffs,rxnsToAdd.stoichCoeffs{i},'rows'); if ~isempty(duplicateRxn) diff --git a/testing/function_tests/tCuration.m b/testing/function_tests/tCuration.m index d4b310fb..ae6327af 100644 --- a/testing/function_tests/tCuration.m +++ b/testing/function_tests/tCuration.m @@ -1,12 +1,95 @@ classdef tCuration < RavenTestCase % tCuration Tests for the table-driven curation functions in curation/. +% +% curateModelFromTables reads its mets/genes/rxns input from tab-delimited +% files; the tests below write those files to a temporary directory. methods (Test) - function curateModelFromTablesNeedsTableFiles(testCase) - % curateModelFromTables reads its mets/genes/rxns input from files - % (tab-delimited / Excel); without those it cannot run. - testCase.assumeFail('Requires curation input table files.'); + function curateModelFromTablesKeepsSupersetReaction(testCase) + % A model reaction that uses the candidate's metabolites *and* + % others is a different reaction, not a duplicate: adding + % A + B -> C to a model holding A + B -> C + H must leave the + % existing reaction alone and add the candidate separately. + model = tCuration.toyModel(); + [coeffsF, infoF] = tCuration.writeRxnTables(testCase, ... + {'A','B','C'}, [-1 -1 1]); + evalc(['out = curateModelFromTables(model,''none'',''rxnsCoeffs'',' ... + 'coeffsF,''rxnsInfo'',infoF);']); + + testCase.verifyNumElements(out.rxns, 2); + % The existing reaction keeps its name and its fourth metabolite + testCase.verifyEqual(out.rxnNames{1}, 'ATP hydrolysis with proton'); + testCase.verifyNumElements(find(out.S(:,1)), 4); + % The candidate is added as a new reaction + testCase.verifyEqual(out.rxnNames{2}, 'newRxn'); + testCase.verifyNumElements(find(out.S(:,2)), 3); + end + + function curateModelFromTablesOverwritesExactDuplicate(testCase) + % An exact stoichiometric match is still treated as the same + % reaction and has its annotation overwritten. + model = tCuration.toyModel(); + [coeffsF, infoF] = tCuration.writeRxnTables(testCase, ... + {'A','B','C','H'}, [-1 -1 1 1]); + evalc(['out = curateModelFromTables(model,''none'',''rxnsCoeffs'',' ... + 'coeffsF,''rxnsInfo'',infoF);']); + + testCase.verifyNumElements(out.rxns, 1); + testCase.verifyEqual(out.rxnNames{1}, 'newRxn'); + end + + end + + methods (Static, Access = private) + + function model = toyModel() + % Single reaction: A + B -> C + H, all in one compartment. + model = struct(); + model.id = 'toy'; + model.rxns = {'R1'}; + model.rxnNames = {'ATP hydrolysis with proton'}; + model.mets = {'m1';'m2';'m3';'m4'}; + model.metNames = {'A';'B';'C';'H'}; + model.comps = {'c'}; + model.compNames = {'cytosol'}; + model.metComps = [1;1;1;1]; + model.S = sparse([-1;-1;1;1]); + model.lb = 0; + model.ub = 1000; + model.rev = 0; + model.c = 0; + model.b = zeros(4,1); + model.genes = {}; + model.grRules = {''}; + model.rxnGeneMat= sparse(1,0); + model.metFormulas = {'';'';'';''}; + model.subSystems = {{''}}; + model.eccodes = {''}; + model.rxnNotes = {''}; + model.rxnReferences = {''}; + model.rxnConfidenceScores = 0; + end + + function [coeffsF, infoF] = writeRxnTables(testCase, metNames, coeffs) + % Write a one-reaction rxnsCoeffs/rxnsInfo pair for 'newRxn'. + d = tempname; mkdir(d); + testCase.addTeardown(@() rmdir(d,'s')); + coeffsF = fullfile(d,'coeffs.tsv'); + infoF = fullfile(d,'info.tsv'); + + fid = fopen(coeffsF,'w'); + fprintf(fid,'rxnIdx\trxnNames\tmetNames\tcomps\tcoefficient\n'); + for i=1:numel(metNames) + fprintf(fid,'1\tnewRxn\t%s\tc\t%d\n', metNames{i}, coeffs(i)); + end + fclose(fid); + + fid = fopen(infoF,'w'); + fprintf(fid,['rxnIdx\trxnNames\tgrRules\tlb\tub\trev\tsubSystems\t' ... + 'eccodes\trxnNotes\trxnReferences\trxnConfidenceScores\n']); + fprintf(fid,'1\tnewRxn\t\t0\t1000\t0\t\t\t\t\t\n'); + fclose(fid); end end