From 9f7a12f0bd7fc3bb79de31e9caaa2a3eae015dce Mon Sep 17 00:00:00 2001 From: Eduard Kerkhoven Date: Fri, 16 Oct 2026 12:43:44 +0200 Subject: [PATCH 1/2] Keep gathered kcats as kcat source files, one per source Gathered kcats are written to data/kcats/.tsv in the project, keyed by conventional reaction and enzyme set rather than model.ec.rxns, so a file stays valid when the ecModel is rebuilt and applies to full and light models alike. geckopy reads and writes the same files. - kcatsFromBrenda, kcatsFromDLKcat and kcatsFromOKP gather one source each, with that source's own options, and write its file. OKP keeps a job and a result per method, and only uses a stored result of the current job. - saveKcatSource writes any kcat list; readKcatSource reads a file; loadKcatSource maps a file to a model and reports rows with an invalid kcat, rows that fit no reaction, and rows that match only part of an enzyme set. - getEcRxnKeys gives the conventional reaction and enzyme set of every model.ec reaction. - fetchOpenKineticsPredictor takes the file to store the result in. Original-commit: https://github.com/SysBioChalmers/GECKO/commit/f3a696cade6ea3129f036723954cff859212dbc1 --- .gitattributes | 3 + src/gather_kcats/fetchOpenKineticsPredictor.m | 11 +- src/gather_kcats/getEcRxnKeys.m | 51 ++++ src/gather_kcats/kcatsFromBrenda.m | 97 ++++++ src/gather_kcats/kcatsFromDLKcat.m | 114 +++++++ src/gather_kcats/kcatsFromOKP.m | 214 ++++++++++++++ src/gather_kcats/loadKcatSource.m | 218 ++++++++++++++ src/gather_kcats/private/checkNameValue.m | 21 ++ src/gather_kcats/private/countRxns.m | 8 + src/gather_kcats/private/formatKcatNumber.m | 24 ++ src/gather_kcats/private/kcatSourceFile.m | 58 ++++ src/gather_kcats/private/kcatSourceJson.m | 61 ++++ src/gather_kcats/private/parseHeaderLine.m | 14 + src/gather_kcats/private/parseKcatNumber.m | 11 + src/gather_kcats/private/unquoteField.m | 11 + src/gather_kcats/private/warnSubset.m | 7 + src/gather_kcats/readKcatSource.m | 137 +++++++++ src/gather_kcats/saveKcatSource.m | 226 ++++++++++++++ test/VerificationMatrix.tsv | 7 + .../ecTestGEM/data/kcatSourceFixture.tsv | 16 + test/unit_tests/kcatSourceTests.m | 279 ++++++++++++++++++ 21 files changed, 1586 insertions(+), 2 deletions(-) create mode 100644 .gitattributes create mode 100644 src/gather_kcats/getEcRxnKeys.m create mode 100644 src/gather_kcats/kcatsFromBrenda.m create mode 100644 src/gather_kcats/kcatsFromDLKcat.m create mode 100644 src/gather_kcats/kcatsFromOKP.m create mode 100644 src/gather_kcats/loadKcatSource.m create mode 100644 src/gather_kcats/private/checkNameValue.m create mode 100644 src/gather_kcats/private/countRxns.m create mode 100644 src/gather_kcats/private/formatKcatNumber.m create mode 100644 src/gather_kcats/private/kcatSourceFile.m create mode 100644 src/gather_kcats/private/kcatSourceJson.m create mode 100644 src/gather_kcats/private/parseHeaderLine.m create mode 100644 src/gather_kcats/private/parseKcatNumber.m create mode 100644 src/gather_kcats/private/unquoteField.m create mode 100644 src/gather_kcats/private/warnSubset.m create mode 100644 src/gather_kcats/readKcatSource.m create mode 100644 src/gather_kcats/saveKcatSource.m create mode 100644 test/unit_tests/ecTestGEM/data/kcatSourceFixture.tsv create mode 100644 test/unit_tests/kcatSourceTests.m diff --git a/.gitattributes b/.gitattributes new file mode 100644 index 00000000..814b1d13 --- /dev/null +++ b/.gitattributes @@ -0,0 +1,3 @@ + +# A fixture that must keep its byte-order mark and CRLF line endings. +test/unit_tests/ecTestGEM/data/kcatSourceFixture.tsv -text diff --git a/src/gather_kcats/fetchOpenKineticsPredictor.m b/src/gather_kcats/fetchOpenKineticsPredictor.m index a992f856..1666e498 100644 --- a/src/gather_kcats/fetchOpenKineticsPredictor.m +++ b/src/gather_kcats/fetchOpenKineticsPredictor.m @@ -29,6 +29,9 @@ % return (default false). % pollInterval : double % seconds between polls when wait=true (default 30). +% outFile : char +% where the result is stored, and read from with useStored=true +% (default: data/OKP_output.csv). % modelAdapter : ModelAdapter % a loaded model adapter (default: the current default model adapter). % @@ -56,7 +59,8 @@ 'jobId', []; ... 'wait', []; ... 'pollInterval', []; ... - 'modelAdapter', []}); + 'modelAdapter', []; ... + 'outFile', []}); useStored = p.useStored; jobId = p.jobId; wait = p.wait; @@ -74,7 +78,10 @@ if isempty(wait); wait = false; end if isempty(pollInterval); pollInterval = 30; end -outFile = fullfile(params.path,'data','OKP_output.csv'); +outFile = p.outFile; +if isempty(outFile) + outFile = fullfile(params.path,'data','OKP_output.csv'); +end %% Stored-result mode: parse the saved file, no API call if useStored diff --git a/src/gather_kcats/getEcRxnKeys.m b/src/gather_kcats/getEcRxnKeys.m new file mode 100644 index 00000000..652750ea --- /dev/null +++ b/src/gather_kcats/getEcRxnKeys.m @@ -0,0 +1,51 @@ +function [baseRxns, enzymeSets] = getEcRxnKeys(model) +% getEcRxnKeys Build-independent keys of the model.ec reactions. +% +% A kcat source file (see saveKcatSource) identifies a reaction by its +% conventional reaction and its enzyme set, not by model.ec.rxns, whose +% identifiers depend on how the ecModel was built. This function gives +% both keys for every model.ec row. +% +% Parameters +% ---------- +% model : struct +% an ecModel in GECKO 3 format (with ecModel.ec structure). +% +% Returns +% ------- +% baseRxns : cell +% per model.ec.rxns entry, the conventional reaction: without a +% trailing _EXP_ (full model) or a leading _ (light model). A +% reverse direction keeps its _REV. +% enzymeSets : cell +% per model.ec.rxns entry, the UniProt accessions of its enzymes +% (model.ec.rxnEnzMat), sorted and joined by ' + '; empty without +% enzymes. +% +% Examples +% -------- +% [baseRxns, enzymeSets] = getEcRxnKeys(ecModel); +% +% See also +% -------- +% saveKcatSource, loadKcatSource + +if model.ec.geckoLight + baseRxns = regexprep(model.ec.rxns, '^\d+_', ''); +else + baseRxns = regexprep(model.ec.rxns, '_EXP_\d+$', ''); +end +baseRxns = baseRxns(:); + +n = numel(model.ec.rxns); +enzymeSets = repmat({''}, n, 1); +[rows, cols] = find(model.ec.rxnEnzMat); +if ~isempty(rows) + % Enzymes per reaction in one pass, rather than a search per reaction. + perRxn = accumarray(rows(:), cols(:), [n 1], @(c) {c}); + for i = find(~cellfun(@isempty, perRxn))' + enzymes = unique(model.ec.enzymes(perRxn{i})); + enzymeSets{i} = strjoin(enzymes(:)', ' + '); + end +end +end diff --git a/src/gather_kcats/kcatsFromBrenda.m b/src/gather_kcats/kcatsFromBrenda.m new file mode 100644 index 00000000..39721861 --- /dev/null +++ b/src/gather_kcats/kcatsFromBrenda.m @@ -0,0 +1,97 @@ +function kcatList = kcatsFromBrenda(model, varargin) +% kcatsFromBrenda Gather kcats from BRENDA and keep them as a source file. +% +% Matches the model's EC numbers and substrates to BRENDA by fuzzy +% matching (see fuzzyKcatMatching), and writes the result to +% data/kcats/brenda.tsv in the project folder (see saveKcatSource). The +% matching runs locally on the BRENDA files of GECKO. +% +% Parameters +% ---------- +% model : struct +% an ecModel in GECKO 3 format (with ecModel.ec structure), with EC +% numbers (getECfromDatabase or getECfromGEM). +% +% Name-Value Arguments +% -------------------- +% aggregation : char +% which BRENDA value per EC number, substrate and organism is used: +% 'max' or 'median' (default: params.kcatAggregation of the model +% adapter, or 'max'). +% forceWClvl : double +% force a minimum wildcard level in the EC numbers (default 0). +% ecRxns : logical +% for which model.ec.rxns kcats are gathered (default: all). +% save : logical +% write the source file (default true). +% saveAs : char +% name of the source file, without .tsv (default 'brenda'). With +% ecRxns, the file then only holds those reactions. +% modelAdapter : ModelAdapter +% a loaded model adapter (default: the current default model adapter). +% +% Returns +% ------- +% kcatList : struct +% the kcat values, as returned by fuzzyKcatMatching. +% +% Examples +% -------- +% kcatsFromBrenda(ecModel); +% kcatList = kcatsFromBrenda(ecModel, 'aggregation', 'median'); +% +% See also +% -------- +% fuzzyKcatMatching, kcatsFromDLKcat, kcatsFromOKP, saveKcatSource, loadKcatSource + +checkNameValue(varargin, {'aggregation','forceWClvl','ecRxns','save','saveAs','modelAdapter'}, 'kcatsFromBrenda'); +p = parseGECKOargs(varargin, { ... + 'aggregation', []; ... + 'forceWClvl', 0; ... + 'ecRxns', []; ... + 'save', true; ... + 'saveAs', 'brenda'; ... + 'modelAdapter', []}); +modelAdapter = p.modelAdapter; +if isempty(modelAdapter) + modelAdapter = ModelAdapterManager.getDefault(); + if isempty(modelAdapter) + error('Either send in a modelAdapter or set the default model adapter in the ModelAdapterManager.') + end +end +aggregation = resolveKcatAggregation(p.aggregation, modelAdapter); +if p.save + file = kcatSourceFile(p.saveAs, [], modelAdapter); % fails before the matching, not after + warnSubset(p.ecRxns, numel(model.ec.rxns), file); +end + +kcatList = fuzzyKcatMatching(model, 'ecRxns', p.ecRxns, 'modelAdapter', modelAdapter, ... + 'forceWClvl', p.forceWClvl, 'kcatAggregation', aggregation); + +if p.save + parameters = struct('aggregation', aggregation, 'force_wildcard_level', p.forceWClvl); + if ~isempty(p.ecRxns) + parameters.ec_rxns = countRxns(p.ecRxns); + end + saveKcatSource(model, kcatList, p.saveAs, 'modelAdapter', modelAdapter, ... + 'generator', 'kcatsFromBrenda', 'parameters', parameters, ... + 'data', brendaRelease(modelAdapter)); + fprintf('BRENDA kcats for %d reactions written to %s\n', sum(kcatList.kcats > 0), file); +end +end + +function release = brendaRelease(modelAdapter) +% The release line at the top of kcat.tsv, e.g. 'BRENDA release 2026.1'. +release = ''; +file = fullfile(modelAdapter.getBrendaDBFolder(), 'kcat.tsv'); +fid = fopen(file, 'r'); +if fid < 0 + return +end +line = fgetl(fid); +fclose(fid); +tok = regexp(line, '^#\s*(BRENDA release [^ ]+)', 'tokens', 'once'); +if ~isempty(tok) + release = tok{1}; +end +end diff --git a/src/gather_kcats/kcatsFromDLKcat.m b/src/gather_kcats/kcatsFromDLKcat.m new file mode 100644 index 00000000..fb50a39e --- /dev/null +++ b/src/gather_kcats/kcatsFromDLKcat.m @@ -0,0 +1,114 @@ +function kcatList = kcatsFromDLKcat(model, varargin) +% kcatsFromDLKcat Gather kcats from DLKcat and keep them as a source file. +% +% Predicts kcats with DLKcat and writes them to data/kcats/dlkcat.tsv in +% the project folder (see saveKcatSource). With 'run', true it writes the +% DLKcat input (writeDLKcatInput) and runs DLKcat through its Docker image +% (runDLKcat), which needs Docker; otherwise it reads the output of an +% earlier run (readDLKcatOutput). A run records how many reactions it +% covered in data/DLKcat_run.txt, so that the source file says so also when +% the output is read later. +% +% Parameters +% ---------- +% model : struct +% an ecModel in GECKO 3 format (with ecModel.ec structure), with +% metabolite SMILES (findMetSmiles) and protein sequences. +% +% Name-Value Arguments +% -------------------- +% run : logical +% write the input and run DLKcat; false reads the existing output in +% data/DLKcat.tsv (default false). +% ecRxns : logical +% for which model.ec.rxns DLKcat is run (default: all). Only used with +% 'run', true. +% save : logical +% write the source file (default true). +% saveAs : char +% name of the source file, without .tsv (default 'dlkcat'). +% modelAdapter : ModelAdapter +% a loaded model adapter (default: the current default model adapter). +% +% Returns +% ------- +% kcatList : struct +% the kcat values, as returned by readDLKcatOutput. +% +% Examples +% -------- +% kcatsFromDLKcat(ecModel, 'run', true); +% kcatList = kcatsFromDLKcat(ecModel); +% +% See also +% -------- +% writeDLKcatInput, runDLKcat, readDLKcatOutput, kcatsFromBrenda, saveKcatSource + +checkNameValue(varargin, {'run','ecRxns','save','saveAs','modelAdapter'}, 'kcatsFromDLKcat'); +p = parseGECKOargs(varargin, { ... + 'run', false; ... + 'ecRxns', []; ... + 'save', true; ... + 'saveAs', 'dlkcat'; ... + 'modelAdapter', []}); +modelAdapter = p.modelAdapter; +if isempty(modelAdapter) + modelAdapter = ModelAdapterManager.getDefault(); + if isempty(modelAdapter) + error('Either send in a modelAdapter or set the default model adapter in the ModelAdapterManager.') + end +end +if ~p.run && ~isempty(p.ecRxns) + error('kcatsFromDLKcat:ecRxns', ['ecRxns selects the reactions DLKcat is run for, so ' ... + 'it needs ''run'', true.']) +end +if p.save + file = kcatSourceFile(p.saveAs, [], modelAdapter); + if p.run + warnSubset(p.ecRxns, numel(model.ec.rxns), file); + end +end + +runFile = fullfile(modelAdapter.params.path, 'data', 'DLKcat_run.txt'); +if p.run + % As runDLKcat: on macOS, MATLAB started from the launcher may not + % have Docker on its path. + if ismac + setenv('PATH', strcat('/usr/local/bin', ':', getenv("PATH"))); + end + [status, ~] = system('docker --version'); + if status ~= 0 + error('kcatsFromDLKcat:noDocker', ['DLKcat runs in a Docker image, but Docker was ' ... + 'not found. Install it from https://www.docker.com, or call kcatsFromDLKcat ' ... + 'without ''run'' to read an earlier DLKcat output.']) + end + writeDLKcatInput(model, 'ecRxns', p.ecRxns, 'modelAdapter', modelAdapter, 'overwrite', true); + runDLKcat('modelAdapter', modelAdapter); + fid = fopen(runFile, 'w'); + if ~isempty(p.ecRxns) + fprintf(fid, 'ecRxns: %d\n', countRxns(p.ecRxns)); + end + fclose(fid); +else + outFile = fullfile(modelAdapter.params.path, 'data', 'DLKcat.tsv'); + if ~isfile(outFile) + error('kcatsFromDLKcat:noOutput', ['No DLKcat output at %s. Run ' ... + 'kcatsFromDLKcat(model, ''run'', true) first.'], outFile) + end +end +kcatList = readDLKcatOutput(model, 'modelAdapter', modelAdapter); + +if p.save + parameters = struct('run', logical(p.run)); + if isfile(runFile) + tok = regexp(fileread(runFile), 'ecRxns:\s*(\d+)', 'tokens', 'once'); + if ~isempty(tok) + parameters.ec_rxns = str2double(tok{1}); + end + end + saveKcatSource(model, kcatList, p.saveAs, 'modelAdapter', modelAdapter, ... + 'generator', 'kcatsFromDLKcat', 'parameters', parameters); + fprintf('DLKcat kcats for %d reactions written to %s\n', ... + numel(unique(kcatList.rxns(kcatList.kcats > 0))), file); +end +end diff --git a/src/gather_kcats/kcatsFromOKP.m b/src/gather_kcats/kcatsFromOKP.m new file mode 100644 index 00000000..4e0ed5d1 --- /dev/null +++ b/src/gather_kcats/kcatsFromOKP.m @@ -0,0 +1,214 @@ +function [kcatList, done] = kcatsFromOKP(model, varargin) +% kcatsFromOKP Gather kcats from OpenKineticsPredictor and keep them as a source file. +% +% Predicts kcats with one OpenKineticsPredictor (OKP) method and writes +% them to data/kcats/okp_.tsv in the project folder (see +% saveKcatSource). OKP is a web service with a job queue, so this function +% can be called repeatedly: +% +% - without a job for this method, it submits one +% (submitOpenKineticsPredictor); +% - with a job, it fetches the result (fetchOpenKineticsPredictor). When +% the job is done the file is written; otherwise it reports the job's +% status and writes nothing. +% +% Each method keeps its own job in data/OKP_job_.txt and its own +% result in data/OKP_output_.csv, so jobs of several methods can +% run at the same time. The job file records which job the stored result +% came from; useStored only uses a result of the current job. An API key +% is needed (see resolveOkpApiKey). +% +% Parameters +% ---------- +% model : struct +% an ecModel in GECKO 3 format (with ecModel.ec structure), with +% metabolite SMILES (findMetSmiles) and protein sequences. +% +% Name-Value Arguments +% -------------------- +% method : char +% OKP method, e.g. 'CataPro', 'CatPred', 'CatRange', 'DLKcat', +% 'EITLEM', 'KinForm-H', 'KinForm-L', 'OmniESI', 'TurNuP', 'UniKP' +% (default: params.okp.method of the model adapter, or 'CataPro'). +% wait : logical +% wait until the job is done (default false). +% resubmit : logical +% submit a new job even if one exists for this method; its earlier +% result is removed (default false). +% useStored : logical +% read the stored result of this method's job, without contacting OKP +% (default false). +% ecRxns : logical +% for which model.ec.rxns kcats are predicted (default: all). Only +% used when a job is submitted. +% save : logical +% write the source file (default true). +% saveAs : char +% name of the source file, without .tsv (default 'okp_', in +% lower case, with other characters than letters and digits as '_'). +% modelAdapter : ModelAdapter +% a loaded model adapter (default: the current default model adapter). +% +% Returns +% ------- +% kcatList : struct +% the kcat values, as returned by readOpenKineticsPredictorOutput; +% empty while the job is not done. +% done : logical +% whether the result was available. +% +% Examples +% -------- +% kcatsFromOKP(ecModel, 'method', 'CataPro'); % submits +% kcatsFromOKP(ecModel, 'method', 'CataPro'); % later: fetches and writes +% kcatsFromOKP(ecModel, 'method', 'UniKP', 'wait', true); +% +% See also +% -------- +% submitOpenKineticsPredictor, fetchOpenKineticsPredictor, kcatsFromBrenda, saveKcatSource + +checkNameValue(varargin, {'method','wait','resubmit','useStored','ecRxns','save','saveAs','modelAdapter'}, 'kcatsFromOKP'); +p = parseGECKOargs(varargin, { ... + 'method', []; ... + 'wait', false; ... + 'resubmit', false; ... + 'useStored', false; ... + 'ecRxns', []; ... + 'save', true; ... + 'saveAs', []; ... + 'modelAdapter', []}); +modelAdapter = p.modelAdapter; +if isempty(modelAdapter) + modelAdapter = ModelAdapterManager.getDefault(); + if isempty(modelAdapter) + error('Either send in a modelAdapter or set the default model adapter in the ModelAdapterManager.') + end +end +params = modelAdapter.params; +method = char(p.method); +if isempty(method) + if isfield(params, 'okp') && isfield(params.okp, 'method') && ~isempty(params.okp.method) + method = char(params.okp.method); + else + method = 'CataPro'; + end +end +tag = lower(regexprep(method, '[^A-Za-z0-9]+', '_')); +dataDir = fullfile(params.path, 'data'); +jobFile = fullfile(dataDir, ['OKP_job_' tag '.txt']); +outFile = fullfile(dataDir, ['OKP_output_' tag '.csv']); +name = char(p.saveAs); +if isempty(name) + name = ['okp_' tag]; +end +job = readJobFile(jobFile); +submitting = ~p.useStored && (p.resubmit || ~isfield(job, 'jobId')); +if ~isempty(p.ecRxns) && ~submitting + error('kcatsFromOKP:ecRxns', ['ecRxns selects the reactions of a new job, but this call ' ... + 'submits none (a job for %s exists; use ''resubmit'', true for a new one).'], method) +end +if p.save + file = kcatSourceFile(name, [], modelAdapter); % fails before submitting, not after + if submitting + warnSubset(p.ecRxns, numel(model.ec.rxns), file); + end +end + +if p.useStored + if ~isfield(job, 'jobId') + adoptSharedResult(model, modelAdapter, dataDir, method, jobFile, outFile); + job = readJobFile(jobFile); + end + if ~isfile(outFile) || ~isfield(job, 'resultOf') || ~strcmp(job.resultOf, job.jobId) + error('kcatsFromOKP:noStored', ['No stored result of the current %s job (%s). Call ' ... + 'kcatsFromOKP without ''useStored'' to fetch it.'], method, jobFile) + end + kcatList = readOpenKineticsPredictorOutput(model, 'outFile', outFile, 'modelAdapter', modelAdapter); + done = true; +else + if submitting + if isfile(outFile) + delete(outFile); + end + jobId = submitOpenKineticsPredictor(model, 'method', method, 'ecRxns', p.ecRxns, ... + 'overwrite', true, 'modelAdapter', modelAdapter); + job = struct('jobId', jobId, 'method', method, 'submittedAt', ... + char(datetime('now', 'TimeZone', 'UTC', 'Format', 'yyyy-MM-dd''T''HH:mm:ss''Z'''))); + if ~isempty(p.ecRxns) + job.ecRxns = sprintf('%d', countRxns(p.ecRxns)); + end + writeJobFile(jobFile, job); + end + [done, kcatList] = fetchOpenKineticsPredictor(model, 'jobId', job.jobId, 'wait', p.wait, ... + 'outFile', outFile, 'modelAdapter', modelAdapter); + if done + job.resultOf = job.jobId; + writeJobFile(jobFile, job); + end +end + +if ~done + fprintf('Call kcatsFromOKP(model, ''method'', ''%s'') again later to collect the result.\n', method); + return +end +if p.save + parameters = struct('method', method); + if isfield(job, 'ecRxns') + parameters.ec_rxns = str2double(job.ecRxns); + end + saveKcatSource(model, kcatList, name, 'modelAdapter', modelAdapter, ... + 'generator', 'kcatsFromOKP', 'parameters', parameters, 'data', ['OKP job ' job.jobId]); + fprintf('OKP (%s) kcats for %d reactions written to %s\n', method, ... + numel(unique(kcatList.rxns(kcatList.kcats > 0))), file); +end +end + +function adoptSharedResult(model, modelAdapter, dataDir, method, jobFile, outFile) +% A result that fetchOpenKineticsPredictor stored by itself, in +% data/OKP_output.csv with data/OKP_job.txt, is taken over for this method +% only when its job is for this method and its predictions are labelled +% with this method. +shared = fullfile(dataDir, 'OKP_output.csv'); +sharedJob = readJobFile(fullfile(dataDir, 'OKP_job.txt')); +ok = isfile(shared) && isfield(sharedJob, 'jobId') && isfield(sharedJob, 'method') && ... + strcmpi(sharedJob.method, method); +if ok + list = readOpenKineticsPredictorOutput(model, 'outFile', shared, 'modelAdapter', modelAdapter); + predicted = list.kcatSource(startsWith(list.kcatSource, 'OKP-') & ... + ~ismember(lower(list.kcatSource), {'okp-brenda','okp-sabio-rk','okp-uniprot'})); + ok = ~isempty(predicted) && all(strcmpi(predicted, ['OKP-' method])); +end +if ~ok + error('kcatsFromOKP:noStored', ['No stored OKP result for %s (%s). Call kcatsFromOKP ' ... + 'without ''useStored'' to submit or fetch a job.'], method, outFile) +end +copyfile(shared, outFile); +sharedJob.method = method; +sharedJob.resultOf = sharedJob.jobId; +writeJobFile(jobFile, sharedJob); +end + +function job = readJobFile(file) +% Fields of an OKP job file ('key: value' per line), as a struct. +job = struct(); +if ~isfile(file) + return +end +lines = strsplit(fileread(file), newline); +for i = 1:numel(lines) + tok = regexp(lines{i}, '^([A-Za-z]\w*):\s*(.*?)\s*$', 'tokens', 'once'); + if ~isempty(tok) + job.(tok{1}) = tok{2}; + end +end +end + +function writeJobFile(file, job) +fid = fopen(file, 'w'); +for f = {'jobId','method','submittedAt','ecRxns','resultOf'} + if isfield(job, f{1}) + fprintf(fid, '%s: %s\n', f{1}, job.(f{1})); + end +end +fclose(fid); +end diff --git a/src/gather_kcats/loadKcatSource.m b/src/gather_kcats/loadKcatSource.m new file mode 100644 index 00000000..87d99867 --- /dev/null +++ b/src/gather_kcats/loadKcatSource.m @@ -0,0 +1,218 @@ +function [kcatList, report] = loadKcatSource(model, name, varargin) +% loadKcatSource Read a kcat source file and map its rows to an ecModel. +% +% Reads data/kcats/.tsv (see saveKcatSource) and maps every row to +% the model.ec reactions it applies to, by conventional reaction and +% enzyme set (see getEcRxnKeys): +% +% 1. A row with rxn and enzymes applies to the model.ec rows with that +% conventional reaction and that enzyme set. +% 2. A row with only enzymes applies to every model.ec row with that +% enzyme set; a row with only rxn to every model.ec row of that +% conventional reaction (curated rows; saveKcatSource always writes +% both). +% 3. A row whose kcat is missing, not a number, infinite, or not above 0 +% is reported as invalid and not applied. +% 4. A row that matches nothing is reported, never dropped silently. If a +% model.ec row has at least half of the row's enzymes, that partial +% match is reported as well, but not applied. +% +% The enzymes of a row may be listed in any order. +% +% Parameters +% ---------- +% model : struct +% an ecModel in GECKO 3 format (with ecModel.ec structure). +% name : char +% name of the source (data/kcats/.tsv in the project folder), or +% the path to a kcat source file. +% +% Name-Value Arguments +% -------------------- +% modelAdapter : ModelAdapter +% a loaded model adapter (default: the current default model adapter). +% +% Returns +% ------- +% kcatList : struct +% kcat values per model.ec reaction, as used by mergeKcats and +% assignKcatValues: source (the source name), rxns, kcats, +% substrates, eccodes, wildcardLvl and origin (from the detail +% column), kcatSource (the per-row source), log10sd, and line (the +% line of the file each entry comes from). +% report : struct +% - invalid: table of the file rows with a missing or invalid kcat; +% - unmatched: table of the file rows that apply to no model.ec row; +% - partial: table of rows that only partly match an enzyme set, with +% the model.ec reaction and its enzymes. +% +% Examples +% -------- +% kcatList = loadKcatSource(ecModel, 'brenda'); +% [kcatList, report] = loadKcatSource(ecModel, 'custom'); +% +% See also +% -------- +% saveKcatSource, readKcatSource, kcatsFromBrenda, getEcRxnKeys + +checkNameValue(varargin, {'modelAdapter'}, 'loadKcatSource'); +p = parseGECKOargs(varargin, {'modelAdapter', []}); + +if isfile(name) + filename = name; +else + modelAdapter = p.modelAdapter; + if isempty(modelAdapter) + modelAdapter = ModelAdapterManager.getDefault(); + if isempty(modelAdapter) + error('Either send in a modelAdapter or set the default model adapter in the ModelAdapterManager.') + end + end + filename = fullfile(modelAdapter.params.path, 'data', 'kcats', [name '.tsv']); + if ~isfile(filename) + error('loadKcatSource:missing', ['No kcat source file %s. Make it with the ' ... + 'kcatsFrom function of this source (kcatsFromBrenda, kcatsFromDLKcat, ' ... + 'kcatsFromOKP), or with saveKcatSource.'], filename) + end +end +[records, meta] = readKcatSource(filename); +sourceName = ''; +if isfield(meta, 'source') + sourceName = meta.source; +end + +[baseRxns, enzymeSets] = getEcRxnKeys(model); +nEc = numel(baseRxns); +nRec = numel(records.kcat); + +% Lookup tables: key -> model.ec rows. +byPair = groupIndex(strcat(baseRxns, '|', enzymeSets)); +byBase = groupIndex(baseRxns); +byEnz = groupIndex(enzymeSets); +% accession -> model.ec rows containing it, for partial matches. +accRows = cell(nEc, 1); +for i = 1:nEc + if ~isempty(enzymeSets{i}) + accRows{i} = strsplit(enzymeSets{i}, ' + '); + end +end +accIdx = repelem((1:nEc)', cellfun(@numel, accRows)); +byAcc = groupIndex([accRows{:}]', accIdx); + +invalid = ~(records.kcat > 0) | ~isfinite(records.kcat); +matchRec = cell(nRec, 1); +matchEc = cell(nRec, 1); +unmatched = false(nRec, 1); +partialRec = zeros(0, 1); partialEc = zeros(0, 1); +for r = find(~invalid)' + rxn = records.rxn{r}; + enz = normalizeEnzymes(records.enzymes{r}); + if ~isempty(rxn) && ~isempty(enz) + hit = lookup(byPair, [rxn '|' enz]); + elseif ~isempty(enz) + hit = lookup(byEnz, enz); + elseif ~isempty(rxn) + hit = lookup(byBase, rxn); + else + hit = []; + end + if isempty(hit) + unmatched(r) = true; + if ~isempty(enz) + wanted = strsplit(enz, ' + '); + candidates = unique(cell2mat(cellfun(@(a) lookup(byAcc, a), wanted(:), 'UniformOutput', false))); + if ~isempty(rxn) + candidates = intersect(candidates, lookup(byBase, rxn)); + end + for c = candidates(:)' + if mean(ismember(wanted, accRows{c})) >= 0.5 + partialRec(end+1, 1) = r; %#ok + partialEc(end+1, 1) = c; %#ok + end + end + end + else + matchRec{r} = repmat(r, numel(hit), 1); + matchEc{r} = hit(:); + end +end +rec = vertcat(matchRec{:}, zeros(0, 1)); +ec = vertcat(matchEc{:}, zeros(0, 1)); + +kcatList.source = sourceName; +kcatList.rxns = reshape(model.ec.rxns(ec), [], 1); +kcatList.kcats = records.kcat(rec); +kcatList.substrates = records.substrate(rec); +kcatList.eccodes = records.ec_code(rec); +kcatList.kcatSource = records.source(rec); +kcatList.log10sd = records.log10_sd(rec); +[kcatList.wildcardLvl, kcatList.origin] = parseDetail(records.detail(rec)); +kcatList.line = records.line(rec); + +report.invalid = table(records.line(invalid), records.rxn(invalid), records.enzymes(invalid), ... + records.kcat(invalid), 'VariableNames', {'line','rxn','enzymes','kcat'}); +report.unmatched = table(records.line(unmatched), records.rxn(unmatched), records.enzymes(unmatched), ... + records.kcat(unmatched), 'VariableNames', {'line','rxn','enzymes','kcat'}); +report.partial = table(records.line(partialRec), records.rxn(partialRec), records.enzymes(partialRec), ... + reshape(model.ec.rxns(partialEc), [], 1), enzymeSets(partialEc), ... + 'VariableNames', {'line','rxn','enzymes','ecRxn','ecEnzymes'}); +if any(invalid) + printOrange(sprintf('%d of %d rows of %s have no valid kcat; see report.invalid.\n', ... + sum(invalid), nRec, filename)); +end +if any(unmatched) + printOrange(sprintf('%d of %d rows of %s apply to no reaction of this model; see report.unmatched.\n', ... + sum(unmatched), nRec, filename)); +end +end + +function map = groupIndex(keys, rows) +% containers.Map from each key to the rows (default 1:n) that have it. +if nargin < 2 + rows = (1:numel(keys))'; +end +map = containers.Map('KeyType', 'char', 'ValueType', 'any'); +if isempty(keys) + return +end +[u, ~, g] = unique(keys(:)); +groups = accumarray(g, rows(:), [], @(x) {sort(x)}); +for i = 1:numel(u) + if ~isempty(u{i}) + map(u{i}) = groups{i}; + end +end +end + +function rows = lookup(map, key) +if isKey(map, key) + rows = map(key); +else + rows = zeros(0, 1); +end +end + +function enz = normalizeEnzymes(enz) +% Sorted and joined as getEcRxnKeys does, so hand-written rows match too. +if isempty(enz) + return +end +parts = strtrim(strsplit(enz, '+')); +parts = unique(parts(~cellfun(@isempty, parts))); +enz = strjoin(sort(parts(:)'), ' + '); +end + +function [wc, origin] = parseDetail(detail) +wc = nan(numel(detail), 1); +origin = nan(numel(detail), 1); +for i = 1:numel(detail) + t = regexp(detail{i}, 'wc=(\d+)', 'tokens', 'once'); + if ~isempty(t) + wc(i) = str2double(t{1}); + end + t = regexp(detail{i}, 'origin=(\d+)', 'tokens', 'once'); + if ~isempty(t) + origin(i) = str2double(t{1}); + end +end +end diff --git a/src/gather_kcats/private/checkNameValue.m b/src/gather_kcats/private/checkNameValue.m new file mode 100644 index 00000000..ed2e7bb0 --- /dev/null +++ b/src/gather_kcats/private/checkNameValue.m @@ -0,0 +1,21 @@ +function checkNameValue(args, names, functionName) +% checkNameValue Accept only name-value pairs with known names. +% +% The kcat source functions have no positional form, so a misspelled +% option name must be an error rather than being read as a positional +% argument by parseGECKOargs. +id = [functionName ':arguments']; +if mod(numel(args), 2) ~= 0 + error(id, 'Options of %s are given as name-value pairs.', functionName) +end +for k = 1:2:numel(args) + name = args{k}; + if ~(ischar(name) || (isstring(name) && isscalar(name))) + error(id, 'Option names of %s are text, not %s. Options: %s.', functionName, ... + class(name), strjoin(names, ', ')) + elseif ~any(strcmp(char(name), names)) + error(id, 'Unknown option of %s: %s. Options: %s.', functionName, char(name), ... + strjoin(names, ', ')) + end +end +end diff --git a/src/gather_kcats/private/countRxns.m b/src/gather_kcats/private/countRxns.m new file mode 100644 index 00000000..fd812fab --- /dev/null +++ b/src/gather_kcats/private/countRxns.m @@ -0,0 +1,8 @@ +function n = countRxns(ecRxns) +% countRxns Number of reactions an ecRxns argument selects (logical or indices). +if islogical(ecRxns) + n = sum(ecRxns); +else + n = numel(unique(ecRxns)); +end +end diff --git a/src/gather_kcats/private/formatKcatNumber.m b/src/gather_kcats/private/formatKcatNumber.m new file mode 100644 index 00000000..f7251a3e --- /dev/null +++ b/src/gather_kcats/private/formatKcatNumber.m @@ -0,0 +1,24 @@ +function t = formatKcatNumber(x) +% formatKcatNumber Number text of kcat source files. +% +% Whole numbers as integers, others in the shortest decimal form that reads +% back as the same double; 'Inf' or '-Inf' for infinity, empty for NaN. +% geckopy's format_number writes the same text. +if isnan(x) + t = ''; +elseif isinf(x) + t = 'Inf'; + if x < 0 + t = '-Inf'; + end +elseif x == round(x) && abs(x) < 1e15 + t = sprintf('%d', x); +else + for digits = 1:17 + t = sprintf('%.*g', digits, x); + if str2double(t) == x + break + end + end +end +end diff --git a/src/gather_kcats/private/kcatSourceFile.m b/src/gather_kcats/private/kcatSourceFile.m new file mode 100644 index 00000000..3d112cf3 --- /dev/null +++ b/src/gather_kcats/private/kcatSourceFile.m @@ -0,0 +1,58 @@ +function filename = kcatSourceFile(name, folder, modelAdapter) +% kcatSourceFile Path of a kcat source file that may be written. +% +% Used by saveKcatSource and the kcatsFrom functions before any +% work is done: checks the name, resolves the folder (data/kcats in the +% project folder by default), and refuses a file that is curated by hand: +% named custom.tsv, or whose header says 'kind: curated' or 'source: +% custom' (in any case, also when a spreadsheet quoted the line). + +name = char(name); +if isempty(regexp(name, '^[A-Za-z0-9_.-]+$', 'once')) + error('saveKcatSource:invalidName', ... + 'The source name "%s" may only contain letters, digits, ''_'', ''-'' and ''.''.', name) +end +if isempty(folder) + if isempty(modelAdapter) + modelAdapter = ModelAdapterManager.getDefault(); + if isempty(modelAdapter) + error('Either send in a modelAdapter or set the default model adapter in the ModelAdapterManager.') + end + end + folder = fullfile(modelAdapter.params.path, 'data', 'kcats'); +end +filename = fullfile(char(folder), [name '.tsv']); +if strcmpi(name, 'custom') + error('saveKcatSource:curated', ['custom.tsv holds curated kcats and is not ' ... + 'written by GECKO; edit it by hand.']) +end +if isfile(filename) + meta = readHeaderLines(filename); + if (isfield(meta, 'kind') && strcmpi(meta.kind, 'curated')) || ... + (isfield(meta, 'source') && strcmpi(meta.source, 'custom')) + error('saveKcatSource:curated', '%s is curated by hand and is not overwritten.', filename) + end +end +end + +function meta = readHeaderLines(filename) +% Only the header lines at the top, so that a damaged table does not stop +% a gathered file from being written again. +meta = struct(); +fid = fopen(filename, 'r', 'n', 'UTF-8'); +if fid < 0 + return +end +cleanup = onCleanup(@() fclose(fid)); +line = fgetl(fid); +if ischar(line) + line = regexprep(line, ['^' char(65279)], ''); +end +while ischar(line) && (isempty(strtrim(line)) || startsWith(unquoteField(line), '#')) + [key, value] = parseHeaderLine(line); + if ~isempty(key) + meta.(key) = value; + end + line = fgetl(fid); +end +end diff --git a/src/gather_kcats/private/kcatSourceJson.m b/src/gather_kcats/private/kcatSourceJson.m new file mode 100644 index 00000000..86808bd1 --- /dev/null +++ b/src/gather_kcats/private/kcatSourceJson.m @@ -0,0 +1,61 @@ +function j = kcatSourceJson(v) +% kcatSourceJson The parameters line of kcat source files as JSON. +% +% Compact JSON with keys sorted at every level. Numbers are written with +% formatKcatNumber, and NaN and infinity as null; logical values as true +% or false; char as a string (non-ASCII characters as they are); a cell +% array as a list; a numeric or logical vector with more than one element +% as a list. geckopy's save_kcat_source writes the same text. +if isstruct(v) && isscalar(v) + f = sort(fieldnames(v)); + parts = cell(1, numel(f)); + for k = 1:numel(f) + parts{k} = [jsonString(f{k}) ':' kcatSourceJson(v.(f{k}))]; + end + j = ['{' strjoin(parts, ',') '}']; +elseif ischar(v) || (isstring(v) && isscalar(v)) + j = jsonString(char(v)); +elseif iscell(v) + j = ['[' strjoin(cellfun(@kcatSourceJson, v(:)', 'UniformOutput', false), ',') ']']; +elseif (isnumeric(v) || islogical(v)) && numel(v) ~= 1 + j = ['[' strjoin(arrayfun(@kcatSourceJson, v(:)', 'UniformOutput', false), ',') ']']; +elseif islogical(v) + j = 'false'; + if v + j = 'true'; + end +elseif isnumeric(v) + if ~isfinite(v) + j = 'null'; + else + j = formatKcatNumber(double(v)); + end +elseif isempty(v) + j = 'null'; +else + error('kcatSourceJson:type', 'Parameters of class %s cannot be written.', class(v)) +end +end + +function s = jsonString(s) +% JSON escapes, as Python's json.dumps(ensure_ascii=False). +out = ''; +for c = s + switch c + case '"', out = [out '\"']; %#ok + case '\', out = [out '\\']; %#ok + case char(8), out = [out '\b']; %#ok + case char(9), out = [out '\t']; %#ok + case char(10), out = [out '\n']; %#ok + case char(12), out = [out '\f']; %#ok + case char(13), out = [out '\r']; %#ok + otherwise + if double(c) < 32 + out = [out sprintf('\\u%04x', double(c))]; %#ok + else + out = [out c]; %#ok + end + end +end +s = ['"' out '"']; +end diff --git a/src/gather_kcats/private/parseHeaderLine.m b/src/gather_kcats/private/parseHeaderLine.m new file mode 100644 index 00000000..ee17bbc4 --- /dev/null +++ b/src/gather_kcats/private/parseHeaderLine.m @@ -0,0 +1,14 @@ +function [key, value] = parseHeaderLine(line) +% parseHeaderLine Key and value of a '# key: value' header line. +% +% The line may be wrapped in quotes by a spreadsheet. The key is returned +% in lower case; it must start with a letter and be at most 63 letters, +% digits or underscores long. Other lines give an empty key. geckopy reads +% header lines the same way. +key = ''; value = ''; +tok = regexp(unquoteField(line), '^\s*#\s*([A-Za-z][A-Za-z0-9_]{0,62})\s*:\s?(.*)$', 'tokens', 'once'); +if ~isempty(tok) + key = lower(tok{1}); + value = strtrim(tok{2}); +end +end diff --git a/src/gather_kcats/private/parseKcatNumber.m b/src/gather_kcats/private/parseKcatNumber.m new file mode 100644 index 00000000..f17b49a0 --- /dev/null +++ b/src/gather_kcats/private/parseKcatNumber.m @@ -0,0 +1,11 @@ +function x = parseKcatNumber(texts) +% parseKcatNumber Read numbers written in kcat source files. +% +% Accepts plain decimal numbers with an optional exponent (12, 0.5, 1e-05, +% 2.5E+3) and Inf; anything else, such as '5,5' from a spreadsheet with a +% decimal comma, becomes NaN instead of being misread. geckopy's +% read_kcat_source accepts the same text. +x = nan(numel(texts), 1); +ok = ~cellfun(@isempty, regexp(texts(:), '^[-+]?(\d+\.?\d*|\.\d+)([eE][-+]?\d+)?$|^[-+]?Inf$', 'once')); +x(ok) = str2double(texts(ok)); +end diff --git a/src/gather_kcats/private/unquoteField.m b/src/gather_kcats/private/unquoteField.m new file mode 100644 index 00000000..36c95aa1 --- /dev/null +++ b/src/gather_kcats/private/unquoteField.m @@ -0,0 +1,11 @@ +function t = unquoteField(t) +% unquoteField A field without surrounding spaces and double quotes. +% +% A field wrapped in double quotes, as a spreadsheet may save it, is +% unwrapped, with "" for a quote inside; spaces inside the quotes are +% trimmed too. +t = strtrim(t); +if numel(t) >= 2 && t(1) == '"' && t(end) == '"' + t = strtrim(strrep(t(2:end-1), '""', '"')); +end +end diff --git a/src/gather_kcats/private/warnSubset.m b/src/gather_kcats/private/warnSubset.m new file mode 100644 index 00000000..00064a5b --- /dev/null +++ b/src/gather_kcats/private/warnSubset.m @@ -0,0 +1,7 @@ +function warnSubset(ecRxns, nRxns, file) +% warnSubset Say that gathering for some reactions replaces the whole file. +if ~isempty(ecRxns) && countRxns(ecRxns) < nRxns && isfile(file) + printOrange(sprintf(['%s is replaced by kcats for %d reactions only. Give saveAs ' ... + 'to keep it.\n'], file, countRxns(ecRxns))); +end +end diff --git a/src/gather_kcats/readKcatSource.m b/src/gather_kcats/readKcatSource.m new file mode 100644 index 00000000..063ab8b0 --- /dev/null +++ b/src/gather_kcats/readKcatSource.m @@ -0,0 +1,137 @@ +function [records, meta] = readKcatSource(filename) +% readKcatSource Read a kcat source file as it is, without a model. +% +% Reads a file written by saveKcatSource or a kcatsFrom function, +% or curated by hand in the same format. Use loadKcatSource to map the +% rows to an ecModel. +% +% The file is read as UTF-8 (a leading byte-order mark is ignored). The +% '#' lines before the column header are the header (keys in lower case); +% '#' lines and blank lines further down are skipped. Fields wrapped in double quotes, as a +% spreadsheet may save them, are unwrapped. kcat and log10_sd must be plain +% numbers (12, 0.5, 1e-05); other text, such as '5,5' with a decimal comma, +% is read as NaN, which loadKcatSource reports. A row with more filled +% fields than there are columns is an error. +% +% Parameters +% ---------- +% filename : char +% path to the kcat source file. +% +% Returns +% ------- +% records : struct +% one field per column, each with one entry per data row: kcat and +% log10_sd as double (NaN when empty or not a number), every other +% column as a cell of char. Columns other than the standard ones are +% kept. records.line is the line number of each row in the file. +% meta : struct +% the '#' header lines, as field: value (for instance meta.source, +% meta.kind). +% +% Examples +% -------- +% [records, meta] = readKcatSource(fullfile(projectPath,'data','kcats','brenda.tsv')); +% +% See also +% -------- +% loadKcatSource, saveKcatSource + +text = readUtf8(filename); +lines = regexp(text, '\r?\n', 'split'); + +meta = struct(); +k = 1; +while k <= numel(lines) && (isempty(strtrim(lines{k})) || startsWith(unquoteField(lines{k}), '#')) + [key, value] = parseHeaderLine(lines{k}); + if ~isempty(key) + meta.(key) = value; + end + k = k + 1; +end +if k > numel(lines) + error('readKcatSource:noHeader', '%s has no column header.', filename) +end +columns = cellfun(@unquoteField, strsplit(lines{k}, sprintf('\t'), ... + 'CollapseDelimiters', false), 'UniformOutput', false); +if ~all(ismember({'kcat','source'}, columns)) || ~any(ismember({'rxn','enzymes'}, columns)) + error('readKcatSource:columns', ['%s needs the columns kcat and source, and rxn ' ... + 'and/or enzymes.'], filename) +end +named = columns(~cellfun(@isempty, columns)); +if numel(unique(named)) < numel(named) || any(strcmp(columns, 'line')) + error('readKcatSource:columns', '%s has a column name twice, or a column named line.', filename) +end +nCol = numel(columns); + +trimmed = strtrim(lines); +skip = cellfun(@isempty, trimmed) | startsWith(trimmed, '#') | startsWith(trimmed, '"#'); +skip(1:k) = true; +kept = find(~skip); +L = string(lines(kept)); +L = L(:); +F = strings(numel(kept), nCol); +nTabs = count(L, sprintf('\t')); +regular = nTabs == nCol - 1; +if any(regular) + % Rows with exactly one field per column, split in one go. + F(regular, :) = reshape(split(L(regular), sprintf('\t'), 2), [], nCol); +end +for r = find(~regular)' + fields = strsplit(char(L(r)), sprintf('\t'), 'CollapseDelimiters', false); + if numel(fields) > nCol + if any(~cellfun(@isempty, strtrim(fields(nCol+1:end)))) + error('readKcatSource:fields', '%s, line %d has more fields than the %d columns.', ... + filename, kept(r), nCol) + end + fields = fields(1:nCol); + end + fields(end+1:nCol) = {''}; + F(r, :) = string(fields); +end +F = strtrim(F); +quoted = find(startsWith(F, '"')); +for q = quoted(:)' + F(q) = unquoteField(char(F(q))); +end +filled = any(F ~= "", 2); % a row of tabs only is skipped +values = cellstr(F(filled, :)); +if isempty(values) + values = cell(0, nCol); +end +lineNo = kept(filled); +lineNo = lineNo(:); +n = numel(lineNo); + +records = struct(); +for c = 1:nCol + if isempty(columns{c}) + continue + end + name = matlab.lang.makeValidName(columns{c}); + v = values(:, c); + if any(strcmp(columns{c}, {'kcat','log10_sd'})) + v = parseKcatNumber(v); + end + records.(name) = v; +end +for c = {'rxn','enzymes','substrate','source','detail','ec_code','note'} + if ~isfield(records, c{1}) + records.(c{1}) = repmat({''}, n, 1); + end +end +if ~isfield(records, 'log10_sd') + records.log10_sd = nan(n, 1); +end +records.line = lineNo; +end + +function text = readUtf8(filename) +fid = fopen(filename, 'r', 'n', 'UTF-8'); +if fid < 0 + error('readKcatSource:read', 'Cannot read %s', filename) +end +cleanup = onCleanup(@() fclose(fid)); +text = fread(fid, [1 Inf], '*char'); +text = regexprep(text, ['^' char(65279)], ''); +end diff --git a/src/gather_kcats/saveKcatSource.m b/src/gather_kcats/saveKcatSource.m new file mode 100644 index 00000000..37101553 --- /dev/null +++ b/src/gather_kcats/saveKcatSource.m @@ -0,0 +1,226 @@ +function filename = saveKcatSource(model, kcatList, name, varargin) +% saveKcatSource Write a kcatList to a kcat source file of the project. +% +% Writes the kcat values of kcatList to data/kcats/.tsv in the +% project folder, the file that loadKcatSource reads back. Rows are keyed +% by their conventional reaction and enzyme set (see getEcRxnKeys), not +% by model.ec.rxns, so the file stays valid when the ecModel is rebuilt. +% kcatsFromBrenda, kcatsFromDLKcat and kcatsFromOKP write these files for +% their sources; use this function for a kcatList made in another way. +% +% The file is UTF-8 text. It starts with '#' header lines (source, +% generator, parameters, data, kind, created), then a tab-separated table +% with the columns rxn, enzymes, substrate, kcat, log10_sd, source, +% detail, ec_code and note, sorted by all columns from rxn onwards. +% Left out are: +% - rows whose kcat is not a positive finite number (0 means no value +% was found); +% - rows of model.ec reactions without enzymes, which have no enzyme +% constraint and could not be told apart from their isozymes. +% +% Parameters +% ---------- +% model : struct +% the ecModel the kcatList refers to. +% kcatList : struct +% kcat values as returned by kcatsFromBrenda, fuzzyKcatMatching, +% readDLKcatOutput, readOpenKineticsPredictorOutput or mergeKcats: +% fields rxns and kcats, optionally source, substrates, eccodes, +% wildcardLvl, origin, kcatSource and log10sd. +% name : char +% name of the source; the file is .tsv. +% +% Name-Value Arguments +% -------------------- +% modelAdapter : ModelAdapter +% a loaded model adapter (default: the current default model adapter). +% folder : char +% folder to write to (default: data/kcats in the project folder). +% generator : char +% what produced the values (default: 'saveKcatSource'). +% parameters : struct +% settings used to produce the values, written as compact JSON with +% sorted keys (default: none). +% data : char +% version of the data the values come from, e.g. a BRENDA release +% (default: empty). +% kind : char +% 'gathered' (default) or 'curated'. +% +% A file whose header says 'kind: curated' or 'source: custom', or named +% custom.tsv, is never written. +% +% Returns +% ------- +% filename : char +% the file written. +% +% Examples +% -------- +% kcatList = fuzzyKcatMatching(ecModel); +% saveKcatSource(ecModel, kcatList, 'mybrenda'); +% +% See also +% -------- +% kcatsFromBrenda, loadKcatSource, readKcatSource, getEcRxnKeys + +names = {'modelAdapter','folder','generator','parameters','data','kind'}; +checkNameValue(varargin, names, 'saveKcatSource'); +p = parseGECKOargs(varargin, { ... + 'modelAdapter', []; ... + 'folder', []; ... + 'generator', 'saveKcatSource'; ... + 'parameters', struct(); ... + 'data', ''; ... + 'kind', 'gathered'}); +kind = validatestring(p.kind, {'gathered','curated'}); +name = char(name); +filename = kcatSourceFile(name, p.folder, p.modelAdapter); + +[baseRxns, enzymeSets] = getEcRxnKeys(model); +[found, idx] = ismember(kcatList.rxns(:), model.ec.rxns); +if ~all(found) + error('saveKcatSource:unknownRxns', '%s', ravenList( ... + 'These kcatList reactions are not in model.ec.rxns:', kcatList.rxns(~found), false)); +end + +n = numel(kcatList.rxns); +kcats = double(kcatList.kcats(:)); +if numel(kcats) ~= n + error('saveKcatSource:fields', 'kcatList.kcats and kcatList.rxns differ in length.') +end +substrate = textColumn(kcatList, 'substrates', n); +ecCode = textColumn(kcatList, 'eccodes', n); +if isfield(kcatList, 'kcatSource') + source = textColumn(kcatList, 'kcatSource', n); +elseif isfield(kcatList, 'source') + source = repmat({char(kcatList.source)}, n, 1); +else + source = repmat({name}, n, 1); +end +wc = numColumn(kcatList, 'wildcardLvl', n); +origin = numColumn(kcatList, 'origin', n); +sd = numColumn(kcatList, 'log10sd', n); + +keep = kcats > 0 & isfinite(kcats); +noEnzymes = keep & cellfun(@isempty, enzymeSets(idx)); +if any(noEnzymes) + fprintf('%d kcats of reactions without enzymes are not written to %s.\n', sum(noEnzymes), filename); +end +keep = find(keep & ~noEnzymes); + +clean = @(c) regexprep(string(c(:)), '[\t\r\n]', ' '); +rows = strings(numel(keep), 9); +if ~isempty(keep) + wcText = strings(numel(keep), 1); + hasWc = ~isnan(wc(keep)); + wcText(hasWc) = "wc=" + string(fix(wc(keep(hasWc)))); + originText = strings(numel(keep), 1); + hasOrigin = ~isnan(origin(keep)); + originText(hasOrigin) = "origin=" + string(fix(origin(keep(hasOrigin)))); + detail = wcText + "; " + originText; + detail(~hasWc) = originText(~hasWc); + detail(~hasOrigin) = wcText(~hasOrigin); + rows(:, 1) = clean(baseRxns(idx(keep))); + rows(:, 2) = clean(enzymeSets(idx(keep))); + rows(:, 3) = clean(substrate(keep)); + rows(:, 4) = formatNumbers(kcats(keep)); + rows(:, 5) = formatNumbers(sd(keep)); + rows(:, 6) = clean(source(keep)); + rows(:, 7) = detail; + rows(:, 8) = clean(ecCode(keep)); + % Sorted by the columns in this order, as geckopy sorts tuples: one key + % per row, the fields joined by a character below any in the text. + [~, order] = sort(join(rows(:, [1 2 3 6 4 5 7 8 9]), char(1), 2)); + rows = rows(order, :); +end + +folder = fileparts(filename); +if ~isfolder(folder) + mkdir(folder); +end +fid = fopen(filename, 'w', 'n', 'UTF-8'); +if fid < 0 + error('saveKcatSource:write', 'Cannot write %s', filename) +end +cleanup = onCleanup(@() fclose(fid)); +fprintf(fid, '# source: %s\n', name); +fprintf(fid, '# generator: %s (GECKO %s)\n', cleanText(p.generator), geckoVersion()); +fprintf(fid, '# parameters: %s\n', kcatSourceJson(p.parameters)); +fprintf(fid, '# data: %s\n', cleanText(p.data)); +fprintf(fid, '# kind: %s\n', kind); +fprintf(fid, '# created: %s\n', char(datetime('now', 'TimeZone', 'UTC', 'Format', 'yyyy-MM-dd''T''HH:mm:ss''Z'''))); +fprintf(fid, 'rxn\tenzymes\tsubstrate\tkcat\tlog10_sd\tsource\tdetail\tec_code\tnote\n'); +if ~isempty(rows) + fprintf(fid, '%s\n', join(rows, sprintf('\t'), 2)); +end +end + +function col = textColumn(s, field, n) +% A kcatList field as an n x 1 cell of char; a cell of several entries +% (e.g. all substrates of a reaction) is joined by '; '. +col = repmat({''}, n, 1); +if ~isfield(s, field) || isempty(s.(field)) + return +end +v = s.(field); +if ischar(v) || isstring(v) + v = cellstr(v); +end +if numel(v) ~= n + error('saveKcatSource:fields', 'kcatList.%s and kcatList.rxns differ in length.', field) +end +for i = 1:n + e = v{i}; + if iscell(e) + e = strjoin(cellfun(@char, e(:)', 'UniformOutput', false), '; '); + elseif isstring(e) + e = char(e); + elseif ~ischar(e) + e = ''; + end + col{i} = e; +end +end + +function col = numColumn(s, field, n) +col = nan(n, 1); +if isfield(s, field) && ~isempty(s.(field)) + if numel(s.(field)) ~= n + error('saveKcatSource:fields', 'kcatList.%s and kcatList.rxns differ in length.', field) + end + col = double(s.(field)(:)); +end +end + +function t = formatNumbers(x) +% formatKcatNumber for a vector, in a few vectorized passes. +x = x(:); +t = strings(numel(x), 1); +whole = isfinite(x) & x == round(x) & abs(x) < 1e15; +t(whole) = compose('%d', x(whole)); +t(isinf(x) & x > 0) = "Inf"; +t(isinf(x) & x < 0) = "-Inf"; +open = find(isfinite(x) & ~whole); +for digits = 1:17 + if isempty(open) + break + end + candidate = compose(sprintf('%%.%dg', digits), x(open)); + hit = str2double(candidate) == x(open); + t(open(hit)) = candidate(hit); + open = open(~hit); +end +end + +function t = cleanText(t) +t = regexprep(char(t), '[\t\r\n]', ' '); +end + +function v = geckoVersion() +v = ''; +file = fullfile(findGECKOroot, 'version.txt'); +if isfile(file) + v = strtrim(fileread(file)); +end +end diff --git a/test/VerificationMatrix.tsv b/test/VerificationMatrix.tsv index fad92a5a..56eca3ff 100644 --- a/test/VerificationMatrix.tsv +++ b/test/VerificationMatrix.tsv @@ -15,6 +15,13 @@ gather_kcats/readDLKcatOutput.m tcm0001, tcm0002 - The test doesn't really check gather_kcats/runDLKcat.m tcm0001, tcm0002 - The test doesn't really check the function, only that we get something gather_kcats/assignKcatValues.m tc0011 - we did not test all possible parameterizations here gather_kcats/writeDLKcatInput.m partly by tc0011 +gather_kcats/getEcRxnKeys.m tc0070, tc0076 +gather_kcats/saveKcatSource.m tc0071, tc0073, tc0076, tc0078 +gather_kcats/readKcatSource.m tc0071, tc0077 +gather_kcats/loadKcatSource.m tc0071, tc0072, tc0073, tc0076, tc0077, tc0078 +gather_kcats/kcatsFromBrenda.m tc0074, tc0078 +gather_kcats/kcatsFromDLKcat.m tc0074, tc0078 - reads an existing DLKcat output; running DLKcat needs Docker +gather_kcats/kcatsFromOKP.m tc0075 - reads a stored result; submitting needs the OKP service gather_kcats/reviewKcatAssignment.m Has only been tested by running it manually, not tested by test cases gather_kcats/ecStatsFromBrenda.m Not verified against real BRENDA data in this environment -- see its own TODO comment get_enzyme_data/calculateMW.m Not explicitly tested - The code has existed for a long time diff --git a/test/unit_tests/ecTestGEM/data/kcatSourceFixture.tsv b/test/unit_tests/ecTestGEM/data/kcatSourceFixture.tsv new file mode 100644 index 00000000..4cfdc9f1 --- /dev/null +++ b/test/unit_tests/ecTestGEM/data/kcatSourceFixture.tsv @@ -0,0 +1,16 @@ +# source: fixture +# generator: hand-written +"# parameters: {""aggregation"":""median""}" +# kind: curated + +rxn enzymes substrate kcat log10_sd source detail ec_code note +R2 P2 + P1 m1 12.5 0.3 brenda wc=0; origin=1 1.1.1.1 enzymes in another order +# a comment between rows +R3 P4 "2,3-bisphospho-D-glycerate" 1e-05 custom quoted substrate + +R5 P5 m2 5,5 custom decimal comma: invalid +R2_REV m2 7 custom reaction only + P3 3 custom enzymes only +R9 X1 -1 custom negative: invalid +R9 X1 4 custom unmatched + diff --git a/test/unit_tests/kcatSourceTests.m b/test/unit_tests/kcatSourceTests.m new file mode 100644 index 00000000..4228cd74 --- /dev/null +++ b/test/unit_tests/kcatSourceTests.m @@ -0,0 +1,279 @@ +%run this test case with the command +%results = runtests('kcatSourceTests.m') + +function tests = kcatSourceTests + tests = functiontests(localfunctions); +end + +function testGetEcRxnKeys_tc0070(testCase) + % Full and light models give the same conventional reactions and + % enzyme sets. + [full, light] = buildKcatModels(); + [baseF, enzF] = getEcRxnKeys(full); + [baseL, enzL] = getEcRxnKeys(light); + verifyEqual(testCase, baseF(strcmp(full.ec.rxns, 'R2_EXP_1')), {'R2'}) + verifyEqual(testCase, baseF(strcmp(full.ec.rxns, 'R2_REV_EXP_2')), {'R2_REV'}) + verifyEqual(testCase, sort(unique(strcat(baseF, '|', enzF))), sort(unique(strcat(baseL, '|', enzL)))) + for i = 1:numel(enzF) + parts = strsplit(enzF{i}, ' + '); + verifyEqual(testCase, parts, sort(parts)) + end +end + +function testSaveLoadRoundTrip_tc0071(testCase) + % What saveKcatSource writes, loadKcatSource maps back to the same + % model.ec rows and values. + [full, ~, adapter] = buildKcatModels(); + list = fuzzyKcatMatching(full, [], adapter); + folder = tempname; + file = saveKcatSource(full, list, 'brenda', 'folder', folder, ... + 'parameters', struct('aggregation','max'), 'data', 'BRENDA release test'); + [records, meta] = readKcatSource(file); + verifyEqual(testCase, meta.source, 'brenda') + verifyEqual(testCase, meta.kind, 'gathered') + verifyEqual(testCase, meta.data, 'BRENDA release test') + verifyEqual(testCase, meta.parameters, '{"aggregation":"max"}') + keep = list.kcats > 0; + verifyEqual(testCase, numel(records.kcat), sum(keep)) + loaded = loadKcatSource(full, file); + verifyEqual(testCase, sort(loaded.rxns), sort(list.rxns(keep))) + [~, a] = sort(list.rxns(keep)); [~, b] = sort(loaded.rxns); + kept = find(keep); + verifyEqual(testCase, loaded.kcats(b), list.kcats(kept(a))) + verifyEqual(testCase, loaded.wildcardLvl(b), list.wildcardLvl(kept(a))) + verifyEqual(testCase, loaded.origin(b), list.origin(kept(a))) + verifyEqual(testCase, loaded.source, 'brenda') + subs = list.substrates(kept(a)); + verifyEqual(testCase, loaded.substrates(b), cellfun(@(c) strjoin(c, '; '), subs, 'UniformOutput', false)) + verifyEqual(testCase, loaded.eccodes(b), list.eccodes(kept(a))) + % Rows are sorted, so a rebuilt model writes the same file. + sorted = sortrows([records.rxn records.enzymes records.substrate records.source]); + verifyEqual(testCase, [records.rxn records.enzymes records.substrate records.source], sorted) +end + +function testFullFileOnLightModel_tc0072(testCase) + % A file written from the full model applies to the light model. + [full, light, adapter] = buildKcatModels(); + list = fuzzyKcatMatching(full, [], adapter); + file = saveKcatSource(full, list, 'brenda', 'folder', tempname); + [loaded, report] = loadKcatSource(light, file); + verifyEmpty(testCase, report.unmatched) + [baseL, enzL] = getEcRxnKeys(light); + [records] = readKcatSource(file); + for r = 1:numel(records.kcat) + hit = strcmp(baseL, records.rxn{r}) & strcmp(enzL, records.enzymes{r}); + verifyEqual(testCase, sort(loaded.rxns(loaded.line == records.line(r))), sort(light.ec.rxns(hit))) + end +end + +function testMatchingRules_tc0073(testCase) + % Rows with only enzymes or only a reaction, unmatched and partial rows. + [full, ~, adapter] = buildKcatModels(); + [base, enz] = getEcRxnKeys(full); + multi = find(contains(enz, ' + '), 1); + parts = strsplit(enz{multi}, ' + '); + folder = tempname; mkdir(folder); + file = fullfile(folder, 'custom.tsv'); + fid = fopen(file, 'w'); + fprintf(fid, '# source: custom\n# kind: curated\n'); + fprintf(fid, 'rxn\tenzymes\tsubstrate\tkcat\tlog10_sd\tsource\tdetail\tec_code\tnote\n'); + fprintf(fid, '\t%s\t\t10\t\tcustom\t\t\tenzymes only, reversed order\n', strjoin(fliplr(parts), ' + ')); + fprintf(fid, '%s\t\t\t20\t\tcustom\t\t\treaction only\n', base{multi}); + fprintf(fid, 'R999\tX1\t\t30\t\tcustom\t\t\tunmatched\n'); + fprintf(fid, '%s\t%s + XNEW\t\t40\t\tcustom\t\t\tpartial\n', base{multi}, enz{multi}); + fclose(fid); + [loaded, report] = loadKcatSource(full, file); + verifyEqual(testCase, sort(loaded.rxns(loaded.line == 4)), sort(full.ec.rxns(strcmp(enz, enz{multi})))) + verifyEqual(testCase, sort(loaded.rxns(loaded.line == 5)), sort(full.ec.rxns(strcmp(base, base{multi})))) + verifyEqual(testCase, report.unmatched.line, [6; 7]) + verifyEqual(testCase, unique(report.partial.line), 7) + % Curated files are never written: by name, by 'kind: curated' and by + % 'source: custom', in any case, and the file is left as it was. + list = fuzzyKcatMatching(full, [], adapter); + before = fileread(file); + verifyError(testCase, @() saveKcatSource(full, list, 'custom', 'folder', folder), 'saveKcatSource:curated') + verifyEqual(testCase, fileread(file), before) + writeLines(fullfile(folder, 'mine.tsv'), {'# source: Custom', 'rxn kcat source'}); + verifyError(testCase, @() saveKcatSource(full, list, 'mine', 'folder', folder), 'saveKcatSource:curated') + writeLines(fullfile(folder, 'mine2.tsv'), {'# kind: Curated', 'rxn kcat source'}); + verifyError(testCase, @() saveKcatSource(full, list, 'mine2', 'folder', folder), 'saveKcatSource:curated') +end + +function testKcatsFromBrendaAndDLKcat_tc0074(testCase) + % The source functions write their files into data/kcats. + [full, ~, adapter] = buildKcatModels(); + project = tempname; + copyfile(fullfile(adapter.params.path, 'data'), fullfile(project, 'data')); + adapter.params.path = project; + list = kcatsFromBrenda(full, 'aggregation', 'median', 'modelAdapter', adapter); + [records, meta] = readKcatSource(fullfile(project, 'data', 'kcats', 'brenda.tsv')); + verifyEqual(testCase, meta.generator(1:15), 'kcatsFromBrenda') + verifyEqual(testCase, meta.parameters, '{"aggregation":"median","force_wildcard_level":0}') + verifyEqual(testCase, meta.data, 'BRENDA release test-fixture') + verifyEqual(testCase, numel(records.kcat), sum(list.kcats > 0)) + + % DLKcat output for this model, read without running DLKcat. + rows = {'R2_EXP_1','G1','m1','C','MRAL','12.5'; 'R3','G4','m1','C','MDFM','3.25'}; + fid = fopen(fullfile(project, 'data', 'DLKcat.tsv'), 'w'); + fprintf(fid, 'rxns\tgenes\tsubstrates\tsmiles\tsequence\tkcat\n'); + for i = 1:size(rows, 1) + fprintf(fid, '%s\t%s\t%s\t%s\t%s\t%s\n', rows{i,:}); + end + fclose(fid); + kcatsFromDLKcat(full, 'modelAdapter', adapter); + [records, meta] = readKcatSource(fullfile(project, 'data', 'kcats', 'dlkcat.tsv')); + verifyEqual(testCase, meta.source, 'dlkcat') + verifyEqual(testCase, sort(records.kcat), [3.25; 12.5]) + verifyEqual(testCase, unique(records.source), {'DLKcat'}) +end + +function testKcatsFromOKPStored_tc0075(testCase) + % Stored OKP results: a method's own result of its current job, or a + % result stored by fetchOpenKineticsPredictor alone when its predictions + % are of the requested method. + [full, ~, adapter] = buildKcatModels(); + full.metSmiles = repmat({''}, numel(full.mets), 1); + full.metSmiles(strcmp(full.mets, 'm1c')) = {'C(C1C)O'}; + project = tempname; + copyfile(fullfile(adapter.params.path, 'data'), fullfile(project, 'data')); + adapter.params.path = project; + data = fullfile(project, 'data'); + okpCsv = {'kcat (1/s),Source kcat,Extra Info kcat,Protein Sequence,Substrate', ... + '12.5,Prediction from CataPro,,MRAL,C(C1C)O'}; + writeLines(fullfile(data, 'OKP_output.csv'), okpCsv); + writeLines(fullfile(data, 'OKP_job.txt'), {'jobId: abc123', 'method: CataPro'}); + [list, done] = kcatsFromOKP(full, 'method', "CataPro", 'useStored', true, 'modelAdapter', adapter); + verifyTrue(testCase, done) + [records, meta] = readKcatSource(fullfile(data, 'kcats', 'okp_catapro.tsv')); + verifyEqual(testCase, meta.data, 'OKP job abc123') + verifyEqual(testCase, numel(records.kcat), sum(list.kcats > 0)) + verifyEqual(testCase, unique(records.source), {'OKP-CataPro'}) + verifyTrue(testCase, isfile(fullfile(data, 'OKP_output_catapro.csv'))) + + % The shared job file names UniKP, but the shared result holds CataPro + % predictions (fetched after UniKP was submitted): refused. + writeLines(fullfile(data, 'OKP_job.txt'), {'jobId: uni1', 'method: UniKP'}); + verifyError(testCase, @() kcatsFromOKP(full, 'method', 'UniKP', 'useStored', true, ... + 'modelAdapter', adapter), 'kcatsFromOKP:noStored') + verifyFalse(testCase, isfile(fullfile(data, 'kcats', 'okp_unikp.tsv'))) + % A method's stored result of an earlier job is refused too. + writeLines(fullfile(data, 'OKP_job_kinform_h.txt'), {'jobId: new2', 'method: KinForm-H', 'resultOf: old1'}); + writeLines(fullfile(data, 'OKP_output_kinform_h.csv'), okpCsv); + verifyError(testCase, @() kcatsFromOKP(full, 'method', 'KinForm-H', 'useStored', true, ... + 'modelAdapter', adapter), 'kcatsFromOKP:noStored') + % ecRxns when no job is submitted, and a misspelled option, are errors. + verifyError(testCase, @() kcatsFromOKP(full, 'method', 'CataPro', 'useStored', true, ... + 'ecRxns', true(numel(full.ec.rxns), 1), 'modelAdapter', adapter), 'kcatsFromOKP:ecRxns') + verifyError(testCase, @() kcatsFromOKP(full, 'method', 'UniKP', 'wiat', true, ... + 'modelAdapter', adapter), 'kcatsFromOKP:arguments') +end + +function testGuardsNumbersAndJson_tc0078(testCase) + [full, ~, adapter] = buildKcatModels(); + project = tempname; + copyfile(fullfile(adapter.params.path, 'data'), fullfile(project, 'data')); + adapter.params.path = project; + folder = fullfile(project, 'data', 'kcats'); + mkdir(folder); + % The curated guard reads quoted header lines and ignores the case of keys. + writeLines(fullfile(folder, 'a.tsv'), {'"# kind: curated"', 'rxn kcat source'}); + writeLines(fullfile(folder, 'b.tsv'), {'# Kind: Curated', 'rxn kcat source'}); + writeLines(fullfile(folder, 'c.tsv'), {'# SOURCE: custom', 'rxn kcat source'}); + list = fuzzyKcatMatching(full, [], adapter); + for f = {'a','b','c'} + verifyError(testCase, @() saveKcatSource(full, list, f{1}, 'modelAdapter', adapter), 'saveKcatSource:curated') + end + % The source functions check the file before they gather. + verifyError(testCase, @() kcatsFromBrenda(full, 'saveAs', "custom", 'modelAdapter', adapter), 'saveKcatSource:curated') + verifyError(testCase, @() kcatsFromDLKcat(full, 'ecRxns', true(numel(full.ec.rxns), 1), ... + 'modelAdapter', adapter), 'kcatsFromDLKcat:ecRxns') + + % Infinite, NaN and non-positive kcats are not written, and not applied. + list.kcats(1) = Inf; list.kcats(2) = NaN; list.kcats(3) = -2; + file = saveKcatSource(full, list, "numbers", 'modelAdapter', adapter, 'parameters', ... + struct('z', 1, 'a', struct('y', ['caf' char(233)], 'b', true), 'n', NaN, ... + 'f', 1/3, 'g', 1e-7, 'h', 1e15, 'k', 1e20, 'l', {{5}}, 's', 'tab quote"')); + [records, meta] = readKcatSource(file); + verifyTrue(testCase, all(isfinite(records.kcat) & records.kcat > 0)) + % The same text as geckopy's test_parameters_json_matches_gecko. + verifyEqual(testCase, meta.parameters, ['{"a":{"b":true,"y":"caf' char(233) '"},"f":0.3333333333333333,' ... + '"g":1e-07,"h":1e+15,"k":1e+20,"l":[5],"n":null,"s":"tab\tquote\"","z":1}']) + writeLines(fullfile(folder, 'inf.tsv'), {'rxn enzymes kcat source', 'R3 P4 Inf custom', 'R3 P4 1e999 custom'}); + [loaded, report] = loadKcatSource(full, fullfile(folder, 'inf.tsv')); + verifyEmpty(testCase, loaded.rxns) + verifyEqual(testCase, report.invalid.line, [2; 3]) + + % Rows are sorted by every column, in the same order as geckopy. + writeLines(fullfile(folder, 'unsorted.tsv'), {'rxn enzymes substrate kcat source detail', ... + 'R3 P4 m1 2 b ', 'R3 P4 m1 2 a wc=1', 'R3 P4 5 a ', 'R2_REV P3 m2 1 a '}); + l2 = loadKcatSource(full, fullfile(folder, 'unsorted.tsv')); + file = saveKcatSource(full, l2, 'sorted', 'modelAdapter', adapter); + records = readKcatSource(file); + verifyEqual(testCase, records.rxn', {'R2_REV','R3','R3','R3'}) + verifyEqual(testCase, records.substrate', {'m2','','m1','m1'}) + verifyEqual(testCase, records.source', {'a','a','a','b'}) +end + +function testEnzymelessRowsAndLightToFull_tc0076(testCase) + % Rows of reactions without enzymes are not written, so they cannot + % spread to their isozymes; a light-model file applies to the full model. + [full, light, adapter] = buildKcatModels(); + list = fuzzyKcatMatching(full, [], adapter); + i = find(strcmp(full.ec.rxns, 'R2_EXP_2')); + full.ec.rxnEnzMat(i, :) = 0; + file = saveKcatSource(full, list, 'brenda', 'folder', tempname); + loaded = loadKcatSource(full, file); + verifyFalse(testCase, any(strcmp(loaded.rxns, 'R2_EXP_2'))) + verifyEqual(testCase, sum(strcmp(loaded.rxns, 'R2_EXP_1')), 1) + + listL = fuzzyKcatMatching(light, [], adapter); + fileL = saveKcatSource(light, listL, 'brenda', 'folder', tempname); + fullAgain = buildKcatModels(); + [~, report] = loadKcatSource(fullAgain, fileL); + verifyEmpty(testCase, report.unmatched) +end + +function testSharedFixture_tc0077(testCase) + % The fixture file, read the same way by geckopy: byte-order mark, CRLF, + % a quoted header line and field, comments and blank lines between + % rows, a row of tabs, a decimal comma, a negative kcat. + file = fullfile(findGECKOroot, 'test', 'unit_tests', 'ecTestGEM', 'data', 'kcatSourceFixture.tsv'); + [records, meta] = readKcatSource(file); + verifyEqual(testCase, meta.source, 'fixture') + verifyEqual(testCase, meta.parameters, '{"aggregation":"median"}') + verifyEqual(testCase, meta.kind, 'curated') + verifyEqual(testCase, records.line', [7 9 11 12 13 14 15]) + verifyEqual(testCase, records.rxn', {'R2','R3','R5','R2_REV','','R9','R9'}) + verifyEqual(testCase, records.substrate{2}, '2,3-bisphospho-D-glycerate') + verifyEqual(testCase, records.kcat', [12.5 1e-05 NaN 7 3 -1 4]) + verifyEqual(testCase, records.log10_sd(1), 0.3) + + [full] = buildKcatModels(); + [loaded, report] = loadKcatSource(full, file); + verifyEqual(testCase, loaded.rxns(loaded.line == 7), {'R2_EXP_1'}) + verifyEqual(testCase, sort(loaded.rxns(loaded.line == 12)), {'R2_REV_EXP_1';'R2_REV_EXP_2'}) + verifyEqual(testCase, sort(loaded.rxns(loaded.line == 13)), {'R2_EXP_2';'R2_REV_EXP_2'}) + verifyEqual(testCase, report.invalid.line, [11; 14]) + verifyEqual(testCase, report.unmatched.line, 15) + + % A row with more filled fields than columns is an error. + bad = [tempname '.tsv']; + writeLines(bad, {'rxn kcat source', 'R2 1 brenda extra'}); + verifyError(testCase, @() readKcatSource(bad), 'readKcatSource:fields') +end + +function writeLines(file, lines) +fid = fopen(file, 'w'); +fprintf(fid, '%s\n', lines{:}); +fclose(fid); +end + +function [full, light, adapter] = buildKcatModels() + geckoPath = findGECKOroot; + adapter = ModelAdapterManager.getAdapter(fullfile(geckoPath,'test','unit_tests','ecTestGEM','TestGEMAdapter.m')); + model = getGeckoTestModel(); + full = makeEcModel(model, false, adapter); + full = getECfromGEM(full); + light = makeEcModel(model, true, adapter); + light = getECfromGEM(light); +end From 4f06b59bb4fbdd47a0a38e5d15c25ac5b73c6566 Mon Sep 17 00:00:00 2001 From: Eduard Kerkhoven Date: Sat, 17 Oct 2026 22:10:47 +0200 Subject: [PATCH 2/2] Read kcat source files in loadKcatSource; no private helper folder loadKcatSource reads the file itself and returns the header lines and the rows as read in report.header and report.records; readKcatSource is gone. The helpers of src/gather_kcats/private are local functions of the functions that use them. kcatsFrom no longer check the file name before gathering; saveKcatSource does when it writes. --- src/gather_kcats/kcatsFromBrenda.m | 35 ++++- src/gather_kcats/kcatsFromDLKcat.m | 35 ++++- src/gather_kcats/kcatsFromOKP.m | 35 ++++- src/gather_kcats/loadKcatSource.m | 162 +++++++++++++++++++- src/gather_kcats/private/checkNameValue.m | 21 --- src/gather_kcats/private/countRxns.m | 8 - src/gather_kcats/private/formatKcatNumber.m | 24 --- src/gather_kcats/private/kcatSourceFile.m | 58 ------- src/gather_kcats/private/kcatSourceJson.m | 61 -------- src/gather_kcats/private/parseHeaderLine.m | 14 -- src/gather_kcats/private/parseKcatNumber.m | 11 -- src/gather_kcats/private/unquoteField.m | 11 -- src/gather_kcats/private/warnSubset.m | 7 - src/gather_kcats/readKcatSource.m | 137 ----------------- src/gather_kcats/saveKcatSource.m | 161 ++++++++++++++++++- test/VerificationMatrix.tsv | 1 - test/unit_tests/kcatSourceTests.m | 27 ++-- 17 files changed, 433 insertions(+), 375 deletions(-) delete mode 100644 src/gather_kcats/private/checkNameValue.m delete mode 100644 src/gather_kcats/private/countRxns.m delete mode 100644 src/gather_kcats/private/formatKcatNumber.m delete mode 100644 src/gather_kcats/private/kcatSourceFile.m delete mode 100644 src/gather_kcats/private/kcatSourceJson.m delete mode 100644 src/gather_kcats/private/parseHeaderLine.m delete mode 100644 src/gather_kcats/private/parseKcatNumber.m delete mode 100644 src/gather_kcats/private/unquoteField.m delete mode 100644 src/gather_kcats/private/warnSubset.m delete mode 100644 src/gather_kcats/readKcatSource.m diff --git a/src/gather_kcats/kcatsFromBrenda.m b/src/gather_kcats/kcatsFromBrenda.m index 39721861..64d29428 100644 --- a/src/gather_kcats/kcatsFromBrenda.m +++ b/src/gather_kcats/kcatsFromBrenda.m @@ -44,7 +44,7 @@ % -------- % fuzzyKcatMatching, kcatsFromDLKcat, kcatsFromOKP, saveKcatSource, loadKcatSource -checkNameValue(varargin, {'aggregation','forceWClvl','ecRxns','save','saveAs','modelAdapter'}, 'kcatsFromBrenda'); +checkOptions(varargin, {'aggregation','forceWClvl','ecRxns','save','saveAs','modelAdapter'}, 'kcatsFromBrenda'); p = parseGECKOargs(varargin, { ... 'aggregation', []; ... 'forceWClvl', 0; ... @@ -61,7 +61,7 @@ end aggregation = resolveKcatAggregation(p.aggregation, modelAdapter); if p.save - file = kcatSourceFile(p.saveAs, [], modelAdapter); % fails before the matching, not after + file = fullfile(modelAdapter.params.path, 'data', 'kcats', [char(p.saveAs) '.tsv']); warnSubset(p.ecRxns, numel(model.ec.rxns), file); end @@ -95,3 +95,34 @@ release = tok{1}; end end + +function n = countRxns(ecRxns) +if islogical(ecRxns) + n = sum(ecRxns); +else + n = numel(unique(ecRxns)); +end +end + +function warnSubset(ecRxns, nRxns, file) +% A subset written under the default name replaces the full file. +if ~isempty(ecRxns) && countRxns(ecRxns) < nRxns && isfile(file) + printOrange(sprintf(['%s is replaced by kcats for %d reactions only. Give saveAs ' ... + 'to keep it.\n'], file, countRxns(ecRxns))); +end +end + +function checkOptions(args, names, functionName) +% Options are name-value pairs with known names, so a misspelled option is +% an error instead of being taken as a positional argument. +id = [functionName ':arguments']; +given = args(1:2:end); +if mod(numel(args), 2) ~= 0 || ~all(cellfun(@(x) ischar(x) || (isstring(x) && isscalar(x)), given)) + error(id, 'Options of %s are name-value pairs. Options: %s.', functionName, strjoin(names, ', ')) +end +unknown = setdiff(cellfun(@char, given, 'UniformOutput', false), names); +if ~isempty(unknown) + error(id, 'Unknown option of %s: %s. Options: %s.', functionName, strjoin(unknown, ', '), ... + strjoin(names, ', ')) +end +end diff --git a/src/gather_kcats/kcatsFromDLKcat.m b/src/gather_kcats/kcatsFromDLKcat.m index fb50a39e..4608e3d2 100644 --- a/src/gather_kcats/kcatsFromDLKcat.m +++ b/src/gather_kcats/kcatsFromDLKcat.m @@ -44,7 +44,7 @@ % -------- % writeDLKcatInput, runDLKcat, readDLKcatOutput, kcatsFromBrenda, saveKcatSource -checkNameValue(varargin, {'run','ecRxns','save','saveAs','modelAdapter'}, 'kcatsFromDLKcat'); +checkOptions(varargin, {'run','ecRxns','save','saveAs','modelAdapter'}, 'kcatsFromDLKcat'); p = parseGECKOargs(varargin, { ... 'run', false; ... 'ecRxns', []; ... @@ -63,7 +63,7 @@ 'it needs ''run'', true.']) end if p.save - file = kcatSourceFile(p.saveAs, [], modelAdapter); + file = fullfile(modelAdapter.params.path, 'data', 'kcats', [char(p.saveAs) '.tsv']); if p.run warnSubset(p.ecRxns, numel(model.ec.rxns), file); end @@ -112,3 +112,34 @@ numel(unique(kcatList.rxns(kcatList.kcats > 0))), file); end end + +function n = countRxns(ecRxns) +if islogical(ecRxns) + n = sum(ecRxns); +else + n = numel(unique(ecRxns)); +end +end + +function warnSubset(ecRxns, nRxns, file) +% A subset written under the default name replaces the full file. +if ~isempty(ecRxns) && countRxns(ecRxns) < nRxns && isfile(file) + printOrange(sprintf(['%s is replaced by kcats for %d reactions only. Give saveAs ' ... + 'to keep it.\n'], file, countRxns(ecRxns))); +end +end + +function checkOptions(args, names, functionName) +% Options are name-value pairs with known names, so a misspelled option is +% an error instead of being taken as a positional argument. +id = [functionName ':arguments']; +given = args(1:2:end); +if mod(numel(args), 2) ~= 0 || ~all(cellfun(@(x) ischar(x) || (isstring(x) && isscalar(x)), given)) + error(id, 'Options of %s are name-value pairs. Options: %s.', functionName, strjoin(names, ', ')) +end +unknown = setdiff(cellfun(@char, given, 'UniformOutput', false), names); +if ~isempty(unknown) + error(id, 'Unknown option of %s: %s. Options: %s.', functionName, strjoin(unknown, ', '), ... + strjoin(names, ', ')) +end +end diff --git a/src/gather_kcats/kcatsFromOKP.m b/src/gather_kcats/kcatsFromOKP.m index 4e0ed5d1..5fe891d4 100644 --- a/src/gather_kcats/kcatsFromOKP.m +++ b/src/gather_kcats/kcatsFromOKP.m @@ -67,7 +67,7 @@ % -------- % submitOpenKineticsPredictor, fetchOpenKineticsPredictor, kcatsFromBrenda, saveKcatSource -checkNameValue(varargin, {'method','wait','resubmit','useStored','ecRxns','save','saveAs','modelAdapter'}, 'kcatsFromOKP'); +checkOptions(varargin, {'method','wait','resubmit','useStored','ecRxns','save','saveAs','modelAdapter'}, 'kcatsFromOKP'); p = parseGECKOargs(varargin, { ... 'method', []; ... 'wait', false; ... @@ -108,7 +108,7 @@ 'submits none (a job for %s exists; use ''resubmit'', true for a new one).'], method) end if p.save - file = kcatSourceFile(name, [], modelAdapter); % fails before submitting, not after + file = fullfile(modelAdapter.params.path, 'data', 'kcats', [char(name) '.tsv']); if submitting warnSubset(p.ecRxns, numel(model.ec.rxns), file); end @@ -212,3 +212,34 @@ function writeJobFile(file, job) end fclose(fid); end + +function n = countRxns(ecRxns) +if islogical(ecRxns) + n = sum(ecRxns); +else + n = numel(unique(ecRxns)); +end +end + +function warnSubset(ecRxns, nRxns, file) +% A subset written under the default name replaces the full file. +if ~isempty(ecRxns) && countRxns(ecRxns) < nRxns && isfile(file) + printOrange(sprintf(['%s is replaced by kcats for %d reactions only. Give saveAs ' ... + 'to keep it.\n'], file, countRxns(ecRxns))); +end +end + +function checkOptions(args, names, functionName) +% Options are name-value pairs with known names, so a misspelled option is +% an error instead of being taken as a positional argument. +id = [functionName ':arguments']; +given = args(1:2:end); +if mod(numel(args), 2) ~= 0 || ~all(cellfun(@(x) ischar(x) || (isstring(x) && isscalar(x)), given)) + error(id, 'Options of %s are name-value pairs. Options: %s.', functionName, strjoin(names, ', ')) +end +unknown = setdiff(cellfun(@char, given, 'UniformOutput', false), names); +if ~isempty(unknown) + error(id, 'Unknown option of %s: %s. Options: %s.', functionName, strjoin(unknown, ', '), ... + strjoin(names, ', ')) +end +end diff --git a/src/gather_kcats/loadKcatSource.m b/src/gather_kcats/loadKcatSource.m index 87d99867..56e2d950 100644 --- a/src/gather_kcats/loadKcatSource.m +++ b/src/gather_kcats/loadKcatSource.m @@ -19,6 +19,14 @@ % % The enzymes of a row may be listed in any order. % +% The file is read as UTF-8 (a leading byte-order mark is ignored). The +% '#' lines before the column header are the header; '#' lines and blank +% lines further down are skipped. Fields wrapped in double quotes, as a +% spreadsheet may save them, are unwrapped. kcat and log10_sd must be plain +% numbers (12, 0.5, 1e-05); other text, such as '5,5' with a decimal comma, +% is reported as invalid. A row with more filled fields than there are +% columns is an error. +% % Parameters % ---------- % model : struct @@ -44,7 +52,12 @@ % - invalid: table of the file rows with a missing or invalid kcat; % - unmatched: table of the file rows that apply to no model.ec row; % - partial: table of rows that only partly match an enzyme set, with -% the model.ec reaction and its enzymes. +% the model.ec reaction and its enzymes; +% - header: the '#' header lines of the file, as field: value (for +% instance header.source, header.parameters); +% - records: the rows of the file as read, one field per column (kcat +% and log10_sd as double, other columns as cell of char) and line, +% the line number of each row. % % Examples % -------- @@ -53,9 +66,9 @@ % % See also % -------- -% saveKcatSource, readKcatSource, kcatsFromBrenda, getEcRxnKeys +% saveKcatSource, kcatsFromBrenda, getEcRxnKeys -checkNameValue(varargin, {'modelAdapter'}, 'loadKcatSource'); +checkOptions(varargin, {'modelAdapter'}, 'loadKcatSource'); p = parseGECKOargs(varargin, {'modelAdapter', []}); if isfile(name) @@ -75,7 +88,7 @@ 'kcatsFromOKP), or with saveKcatSource.'], filename) end end -[records, meta] = readKcatSource(filename); +[records, meta] = readSourceFile(filename); sourceName = ''; if isfield(meta, 'source') sourceName = meta.source; @@ -164,6 +177,8 @@ printOrange(sprintf('%d of %d rows of %s apply to no reaction of this model; see report.unmatched.\n', ... sum(unmatched), nRec, filename)); end +report.header = meta; +report.records = records; end function map = groupIndex(keys, rows) @@ -216,3 +231,142 @@ end end end + +function [records, meta] = readSourceFile(filename) +% The rows and '#' header lines of a kcat source file, as they are. +text = readUtf8(filename); +lines = regexp(text, '\r?\n', 'split'); + +meta = struct(); +k = 1; +while k <= numel(lines) && (isempty(strtrim(lines{k})) || startsWith(unquoteField(lines{k}), '#')) + [key, value] = parseHeaderLine(lines{k}); + if ~isempty(key) + meta.(key) = value; + end + k = k + 1; +end +if k > numel(lines) + error('loadKcatSource:noHeader', '%s has no column header.', filename) +end +columns = cellfun(@unquoteField, strsplit(lines{k}, sprintf('\t'), ... + 'CollapseDelimiters', false), 'UniformOutput', false); +if ~all(ismember({'kcat','source'}, columns)) || ~any(ismember({'rxn','enzymes'}, columns)) + error('loadKcatSource:columns', ['%s needs the columns kcat and source, and rxn ' ... + 'and/or enzymes.'], filename) +end +named = columns(~cellfun(@isempty, columns)); +if numel(unique(named)) < numel(named) || any(strcmp(columns, 'line')) + error('loadKcatSource:columns', '%s has a column name twice, or a column named line.', filename) +end +nCol = numel(columns); + +trimmed = strtrim(lines); +skip = cellfun(@isempty, trimmed) | startsWith(trimmed, '#') | startsWith(trimmed, '"#'); +skip(1:k) = true; +kept = find(~skip); +L = string(lines(kept)); +L = L(:); +F = strings(numel(kept), nCol); +nTabs = count(L, sprintf('\t')); +regular = nTabs == nCol - 1; +if any(regular) + % Rows with exactly one field per column, split in one go. + F(regular, :) = reshape(split(L(regular), sprintf('\t'), 2), [], nCol); +end +for r = find(~regular)' + fields = strsplit(char(L(r)), sprintf('\t'), 'CollapseDelimiters', false); + if numel(fields) > nCol + if any(~cellfun(@isempty, strtrim(fields(nCol+1:end)))) + error('loadKcatSource:fields', '%s, line %d has more fields than the %d columns.', ... + filename, kept(r), nCol) + end + fields = fields(1:nCol); + end + fields(end+1:nCol) = {''}; + F(r, :) = string(fields); +end +F = strtrim(F); +quoted = find(startsWith(F, '"')); +for q = quoted(:)' + F(q) = unquoteField(char(F(q))); +end +filled = any(F ~= "", 2); % a row of tabs only is skipped +values = cellstr(F(filled, :)); +if isempty(values) + values = cell(0, nCol); +end +lineNo = kept(filled); +lineNo = lineNo(:); +n = numel(lineNo); + +records = struct(); +for c = 1:nCol + if isempty(columns{c}) + continue + end + name = matlab.lang.makeValidName(columns{c}); + v = values(:, c); + if any(strcmp(columns{c}, {'kcat','log10_sd'})) + v = parseKcatNumber(v); + end + records.(name) = v; +end +for c = {'rxn','enzymes','substrate','source','detail','ec_code','note'} + if ~isfield(records, c{1}) + records.(c{1}) = repmat({''}, n, 1); + end +end +if ~isfield(records, 'log10_sd') + records.log10_sd = nan(n, 1); +end +records.line = lineNo; +end + +function text = readUtf8(filename) +fid = fopen(filename, 'r', 'n', 'UTF-8'); +if fid < 0 + error('loadKcatSource:read', 'Cannot read %s', filename) +end +cleanup = onCleanup(@() fclose(fid)); +text = fread(fid, [1 Inf], '*char'); +text = regexprep(text, ['^' char(65279)], ''); +end + +function x = parseKcatNumber(texts) +% Plain numbers only: '5,5' or '12 /s' read as NaN. +x = nan(numel(texts), 1); +ok = ~cellfun(@isempty, regexp(texts(:), '^[-+]?(\d+\.?\d*|\.\d+)([eE][-+]?\d+)?$|^[-+]?Inf$', 'once')); +x(ok) = str2double(texts(ok)); +end + +function [key, value] = parseHeaderLine(line) +key = ''; value = ''; +tok = regexp(unquoteField(line), '^\s*#\s*([A-Za-z][A-Za-z0-9_]{0,62})\s*:\s?(.*)$', 'tokens', 'once'); +if ~isempty(tok) + key = lower(tok{1}); + value = strtrim(tok{2}); +end +end + +function t = unquoteField(t) +t = strtrim(t); +if numel(t) >= 2 && t(1) == '"' && t(end) == '"' + t = strtrim(strrep(t(2:end-1), '""', '"')); +end +end + +function checkOptions(args, names, functionName) +% Options are name-value pairs with known names, so a misspelled option is +% an error instead of being taken as a positional argument. +id = [functionName ':arguments']; +given = args(1:2:end); +if mod(numel(args), 2) ~= 0 || ~all(cellfun(@(x) ischar(x) || (isstring(x) && isscalar(x)), given)) + error(id, 'Options of %s are name-value pairs. Options: %s.', functionName, strjoin(names, ', ')) +end +unknown = setdiff(cellfun(@char, given, 'UniformOutput', false), names); +if ~isempty(unknown) + error(id, 'Unknown option of %s: %s. Options: %s.', functionName, strjoin(unknown, ', '), ... + strjoin(names, ', ')) +end +end diff --git a/src/gather_kcats/private/checkNameValue.m b/src/gather_kcats/private/checkNameValue.m deleted file mode 100644 index ed2e7bb0..00000000 --- a/src/gather_kcats/private/checkNameValue.m +++ /dev/null @@ -1,21 +0,0 @@ -function checkNameValue(args, names, functionName) -% checkNameValue Accept only name-value pairs with known names. -% -% The kcat source functions have no positional form, so a misspelled -% option name must be an error rather than being read as a positional -% argument by parseGECKOargs. -id = [functionName ':arguments']; -if mod(numel(args), 2) ~= 0 - error(id, 'Options of %s are given as name-value pairs.', functionName) -end -for k = 1:2:numel(args) - name = args{k}; - if ~(ischar(name) || (isstring(name) && isscalar(name))) - error(id, 'Option names of %s are text, not %s. Options: %s.', functionName, ... - class(name), strjoin(names, ', ')) - elseif ~any(strcmp(char(name), names)) - error(id, 'Unknown option of %s: %s. Options: %s.', functionName, char(name), ... - strjoin(names, ', ')) - end -end -end diff --git a/src/gather_kcats/private/countRxns.m b/src/gather_kcats/private/countRxns.m deleted file mode 100644 index fd812fab..00000000 --- a/src/gather_kcats/private/countRxns.m +++ /dev/null @@ -1,8 +0,0 @@ -function n = countRxns(ecRxns) -% countRxns Number of reactions an ecRxns argument selects (logical or indices). -if islogical(ecRxns) - n = sum(ecRxns); -else - n = numel(unique(ecRxns)); -end -end diff --git a/src/gather_kcats/private/formatKcatNumber.m b/src/gather_kcats/private/formatKcatNumber.m deleted file mode 100644 index f7251a3e..00000000 --- a/src/gather_kcats/private/formatKcatNumber.m +++ /dev/null @@ -1,24 +0,0 @@ -function t = formatKcatNumber(x) -% formatKcatNumber Number text of kcat source files. -% -% Whole numbers as integers, others in the shortest decimal form that reads -% back as the same double; 'Inf' or '-Inf' for infinity, empty for NaN. -% geckopy's format_number writes the same text. -if isnan(x) - t = ''; -elseif isinf(x) - t = 'Inf'; - if x < 0 - t = '-Inf'; - end -elseif x == round(x) && abs(x) < 1e15 - t = sprintf('%d', x); -else - for digits = 1:17 - t = sprintf('%.*g', digits, x); - if str2double(t) == x - break - end - end -end -end diff --git a/src/gather_kcats/private/kcatSourceFile.m b/src/gather_kcats/private/kcatSourceFile.m deleted file mode 100644 index 3d112cf3..00000000 --- a/src/gather_kcats/private/kcatSourceFile.m +++ /dev/null @@ -1,58 +0,0 @@ -function filename = kcatSourceFile(name, folder, modelAdapter) -% kcatSourceFile Path of a kcat source file that may be written. -% -% Used by saveKcatSource and the kcatsFrom functions before any -% work is done: checks the name, resolves the folder (data/kcats in the -% project folder by default), and refuses a file that is curated by hand: -% named custom.tsv, or whose header says 'kind: curated' or 'source: -% custom' (in any case, also when a spreadsheet quoted the line). - -name = char(name); -if isempty(regexp(name, '^[A-Za-z0-9_.-]+$', 'once')) - error('saveKcatSource:invalidName', ... - 'The source name "%s" may only contain letters, digits, ''_'', ''-'' and ''.''.', name) -end -if isempty(folder) - if isempty(modelAdapter) - modelAdapter = ModelAdapterManager.getDefault(); - if isempty(modelAdapter) - error('Either send in a modelAdapter or set the default model adapter in the ModelAdapterManager.') - end - end - folder = fullfile(modelAdapter.params.path, 'data', 'kcats'); -end -filename = fullfile(char(folder), [name '.tsv']); -if strcmpi(name, 'custom') - error('saveKcatSource:curated', ['custom.tsv holds curated kcats and is not ' ... - 'written by GECKO; edit it by hand.']) -end -if isfile(filename) - meta = readHeaderLines(filename); - if (isfield(meta, 'kind') && strcmpi(meta.kind, 'curated')) || ... - (isfield(meta, 'source') && strcmpi(meta.source, 'custom')) - error('saveKcatSource:curated', '%s is curated by hand and is not overwritten.', filename) - end -end -end - -function meta = readHeaderLines(filename) -% Only the header lines at the top, so that a damaged table does not stop -% a gathered file from being written again. -meta = struct(); -fid = fopen(filename, 'r', 'n', 'UTF-8'); -if fid < 0 - return -end -cleanup = onCleanup(@() fclose(fid)); -line = fgetl(fid); -if ischar(line) - line = regexprep(line, ['^' char(65279)], ''); -end -while ischar(line) && (isempty(strtrim(line)) || startsWith(unquoteField(line), '#')) - [key, value] = parseHeaderLine(line); - if ~isempty(key) - meta.(key) = value; - end - line = fgetl(fid); -end -end diff --git a/src/gather_kcats/private/kcatSourceJson.m b/src/gather_kcats/private/kcatSourceJson.m deleted file mode 100644 index 86808bd1..00000000 --- a/src/gather_kcats/private/kcatSourceJson.m +++ /dev/null @@ -1,61 +0,0 @@ -function j = kcatSourceJson(v) -% kcatSourceJson The parameters line of kcat source files as JSON. -% -% Compact JSON with keys sorted at every level. Numbers are written with -% formatKcatNumber, and NaN and infinity as null; logical values as true -% or false; char as a string (non-ASCII characters as they are); a cell -% array as a list; a numeric or logical vector with more than one element -% as a list. geckopy's save_kcat_source writes the same text. -if isstruct(v) && isscalar(v) - f = sort(fieldnames(v)); - parts = cell(1, numel(f)); - for k = 1:numel(f) - parts{k} = [jsonString(f{k}) ':' kcatSourceJson(v.(f{k}))]; - end - j = ['{' strjoin(parts, ',') '}']; -elseif ischar(v) || (isstring(v) && isscalar(v)) - j = jsonString(char(v)); -elseif iscell(v) - j = ['[' strjoin(cellfun(@kcatSourceJson, v(:)', 'UniformOutput', false), ',') ']']; -elseif (isnumeric(v) || islogical(v)) && numel(v) ~= 1 - j = ['[' strjoin(arrayfun(@kcatSourceJson, v(:)', 'UniformOutput', false), ',') ']']; -elseif islogical(v) - j = 'false'; - if v - j = 'true'; - end -elseif isnumeric(v) - if ~isfinite(v) - j = 'null'; - else - j = formatKcatNumber(double(v)); - end -elseif isempty(v) - j = 'null'; -else - error('kcatSourceJson:type', 'Parameters of class %s cannot be written.', class(v)) -end -end - -function s = jsonString(s) -% JSON escapes, as Python's json.dumps(ensure_ascii=False). -out = ''; -for c = s - switch c - case '"', out = [out '\"']; %#ok - case '\', out = [out '\\']; %#ok - case char(8), out = [out '\b']; %#ok - case char(9), out = [out '\t']; %#ok - case char(10), out = [out '\n']; %#ok - case char(12), out = [out '\f']; %#ok - case char(13), out = [out '\r']; %#ok - otherwise - if double(c) < 32 - out = [out sprintf('\\u%04x', double(c))]; %#ok - else - out = [out c]; %#ok - end - end -end -s = ['"' out '"']; -end diff --git a/src/gather_kcats/private/parseHeaderLine.m b/src/gather_kcats/private/parseHeaderLine.m deleted file mode 100644 index ee17bbc4..00000000 --- a/src/gather_kcats/private/parseHeaderLine.m +++ /dev/null @@ -1,14 +0,0 @@ -function [key, value] = parseHeaderLine(line) -% parseHeaderLine Key and value of a '# key: value' header line. -% -% The line may be wrapped in quotes by a spreadsheet. The key is returned -% in lower case; it must start with a letter and be at most 63 letters, -% digits or underscores long. Other lines give an empty key. geckopy reads -% header lines the same way. -key = ''; value = ''; -tok = regexp(unquoteField(line), '^\s*#\s*([A-Za-z][A-Za-z0-9_]{0,62})\s*:\s?(.*)$', 'tokens', 'once'); -if ~isempty(tok) - key = lower(tok{1}); - value = strtrim(tok{2}); -end -end diff --git a/src/gather_kcats/private/parseKcatNumber.m b/src/gather_kcats/private/parseKcatNumber.m deleted file mode 100644 index f17b49a0..00000000 --- a/src/gather_kcats/private/parseKcatNumber.m +++ /dev/null @@ -1,11 +0,0 @@ -function x = parseKcatNumber(texts) -% parseKcatNumber Read numbers written in kcat source files. -% -% Accepts plain decimal numbers with an optional exponent (12, 0.5, 1e-05, -% 2.5E+3) and Inf; anything else, such as '5,5' from a spreadsheet with a -% decimal comma, becomes NaN instead of being misread. geckopy's -% read_kcat_source accepts the same text. -x = nan(numel(texts), 1); -ok = ~cellfun(@isempty, regexp(texts(:), '^[-+]?(\d+\.?\d*|\.\d+)([eE][-+]?\d+)?$|^[-+]?Inf$', 'once')); -x(ok) = str2double(texts(ok)); -end diff --git a/src/gather_kcats/private/unquoteField.m b/src/gather_kcats/private/unquoteField.m deleted file mode 100644 index 36c95aa1..00000000 --- a/src/gather_kcats/private/unquoteField.m +++ /dev/null @@ -1,11 +0,0 @@ -function t = unquoteField(t) -% unquoteField A field without surrounding spaces and double quotes. -% -% A field wrapped in double quotes, as a spreadsheet may save it, is -% unwrapped, with "" for a quote inside; spaces inside the quotes are -% trimmed too. -t = strtrim(t); -if numel(t) >= 2 && t(1) == '"' && t(end) == '"' - t = strtrim(strrep(t(2:end-1), '""', '"')); -end -end diff --git a/src/gather_kcats/private/warnSubset.m b/src/gather_kcats/private/warnSubset.m deleted file mode 100644 index 00064a5b..00000000 --- a/src/gather_kcats/private/warnSubset.m +++ /dev/null @@ -1,7 +0,0 @@ -function warnSubset(ecRxns, nRxns, file) -% warnSubset Say that gathering for some reactions replaces the whole file. -if ~isempty(ecRxns) && countRxns(ecRxns) < nRxns && isfile(file) - printOrange(sprintf(['%s is replaced by kcats for %d reactions only. Give saveAs ' ... - 'to keep it.\n'], file, countRxns(ecRxns))); -end -end diff --git a/src/gather_kcats/readKcatSource.m b/src/gather_kcats/readKcatSource.m deleted file mode 100644 index 063ab8b0..00000000 --- a/src/gather_kcats/readKcatSource.m +++ /dev/null @@ -1,137 +0,0 @@ -function [records, meta] = readKcatSource(filename) -% readKcatSource Read a kcat source file as it is, without a model. -% -% Reads a file written by saveKcatSource or a kcatsFrom function, -% or curated by hand in the same format. Use loadKcatSource to map the -% rows to an ecModel. -% -% The file is read as UTF-8 (a leading byte-order mark is ignored). The -% '#' lines before the column header are the header (keys in lower case); -% '#' lines and blank lines further down are skipped. Fields wrapped in double quotes, as a -% spreadsheet may save them, are unwrapped. kcat and log10_sd must be plain -% numbers (12, 0.5, 1e-05); other text, such as '5,5' with a decimal comma, -% is read as NaN, which loadKcatSource reports. A row with more filled -% fields than there are columns is an error. -% -% Parameters -% ---------- -% filename : char -% path to the kcat source file. -% -% Returns -% ------- -% records : struct -% one field per column, each with one entry per data row: kcat and -% log10_sd as double (NaN when empty or not a number), every other -% column as a cell of char. Columns other than the standard ones are -% kept. records.line is the line number of each row in the file. -% meta : struct -% the '#' header lines, as field: value (for instance meta.source, -% meta.kind). -% -% Examples -% -------- -% [records, meta] = readKcatSource(fullfile(projectPath,'data','kcats','brenda.tsv')); -% -% See also -% -------- -% loadKcatSource, saveKcatSource - -text = readUtf8(filename); -lines = regexp(text, '\r?\n', 'split'); - -meta = struct(); -k = 1; -while k <= numel(lines) && (isempty(strtrim(lines{k})) || startsWith(unquoteField(lines{k}), '#')) - [key, value] = parseHeaderLine(lines{k}); - if ~isempty(key) - meta.(key) = value; - end - k = k + 1; -end -if k > numel(lines) - error('readKcatSource:noHeader', '%s has no column header.', filename) -end -columns = cellfun(@unquoteField, strsplit(lines{k}, sprintf('\t'), ... - 'CollapseDelimiters', false), 'UniformOutput', false); -if ~all(ismember({'kcat','source'}, columns)) || ~any(ismember({'rxn','enzymes'}, columns)) - error('readKcatSource:columns', ['%s needs the columns kcat and source, and rxn ' ... - 'and/or enzymes.'], filename) -end -named = columns(~cellfun(@isempty, columns)); -if numel(unique(named)) < numel(named) || any(strcmp(columns, 'line')) - error('readKcatSource:columns', '%s has a column name twice, or a column named line.', filename) -end -nCol = numel(columns); - -trimmed = strtrim(lines); -skip = cellfun(@isempty, trimmed) | startsWith(trimmed, '#') | startsWith(trimmed, '"#'); -skip(1:k) = true; -kept = find(~skip); -L = string(lines(kept)); -L = L(:); -F = strings(numel(kept), nCol); -nTabs = count(L, sprintf('\t')); -regular = nTabs == nCol - 1; -if any(regular) - % Rows with exactly one field per column, split in one go. - F(regular, :) = reshape(split(L(regular), sprintf('\t'), 2), [], nCol); -end -for r = find(~regular)' - fields = strsplit(char(L(r)), sprintf('\t'), 'CollapseDelimiters', false); - if numel(fields) > nCol - if any(~cellfun(@isempty, strtrim(fields(nCol+1:end)))) - error('readKcatSource:fields', '%s, line %d has more fields than the %d columns.', ... - filename, kept(r), nCol) - end - fields = fields(1:nCol); - end - fields(end+1:nCol) = {''}; - F(r, :) = string(fields); -end -F = strtrim(F); -quoted = find(startsWith(F, '"')); -for q = quoted(:)' - F(q) = unquoteField(char(F(q))); -end -filled = any(F ~= "", 2); % a row of tabs only is skipped -values = cellstr(F(filled, :)); -if isempty(values) - values = cell(0, nCol); -end -lineNo = kept(filled); -lineNo = lineNo(:); -n = numel(lineNo); - -records = struct(); -for c = 1:nCol - if isempty(columns{c}) - continue - end - name = matlab.lang.makeValidName(columns{c}); - v = values(:, c); - if any(strcmp(columns{c}, {'kcat','log10_sd'})) - v = parseKcatNumber(v); - end - records.(name) = v; -end -for c = {'rxn','enzymes','substrate','source','detail','ec_code','note'} - if ~isfield(records, c{1}) - records.(c{1}) = repmat({''}, n, 1); - end -end -if ~isfield(records, 'log10_sd') - records.log10_sd = nan(n, 1); -end -records.line = lineNo; -end - -function text = readUtf8(filename) -fid = fopen(filename, 'r', 'n', 'UTF-8'); -if fid < 0 - error('readKcatSource:read', 'Cannot read %s', filename) -end -cleanup = onCleanup(@() fclose(fid)); -text = fread(fid, [1 Inf], '*char'); -text = regexprep(text, ['^' char(65279)], ''); -end diff --git a/src/gather_kcats/saveKcatSource.m b/src/gather_kcats/saveKcatSource.m index 37101553..330025c0 100644 --- a/src/gather_kcats/saveKcatSource.m +++ b/src/gather_kcats/saveKcatSource.m @@ -62,10 +62,10 @@ % % See also % -------- -% kcatsFromBrenda, loadKcatSource, readKcatSource, getEcRxnKeys +% kcatsFromBrenda, loadKcatSource, getEcRxnKeys names = {'modelAdapter','folder','generator','parameters','data','kind'}; -checkNameValue(varargin, names, 'saveKcatSource'); +checkOptions(varargin, names, 'saveKcatSource'); p = parseGECKOargs(varargin, { ... 'modelAdapter', []; ... 'folder', []; ... @@ -224,3 +224,160 @@ v = strtrim(fileread(file)); end end + +function filename = kcatSourceFile(name, folder, modelAdapter) +name = char(name); +if isempty(regexp(name, '^[A-Za-z0-9_.-]+$', 'once')) + error('saveKcatSource:invalidName', ... + 'The source name "%s" may only contain letters, digits, ''_'', ''-'' and ''.''.', name) +end +if isempty(folder) + if isempty(modelAdapter) + modelAdapter = ModelAdapterManager.getDefault(); + if isempty(modelAdapter) + error('Either send in a modelAdapter or set the default model adapter in the ModelAdapterManager.') + end + end + folder = fullfile(modelAdapter.params.path, 'data', 'kcats'); +end +filename = fullfile(char(folder), [name '.tsv']); +if strcmpi(name, 'custom') + error('saveKcatSource:curated', ['custom.tsv holds curated kcats and is not ' ... + 'written by GECKO; edit it by hand.']) +end +if isfile(filename) + meta = readHeaderLines(filename); + if (isfield(meta, 'kind') && strcmpi(meta.kind, 'curated')) || ... + (isfield(meta, 'source') && strcmpi(meta.source, 'custom')) + error('saveKcatSource:curated', '%s is curated by hand and is not overwritten.', filename) + end +end +end + +function meta = readHeaderLines(filename) +% Only the header lines at the top, so that a damaged table does not stop +% a gathered file from being written again. +meta = struct(); +fid = fopen(filename, 'r', 'n', 'UTF-8'); +if fid < 0 + return +end +cleanup = onCleanup(@() fclose(fid)); +line = fgetl(fid); +if ischar(line) + line = regexprep(line, ['^' char(65279)], ''); +end +while ischar(line) && (isempty(strtrim(line)) || startsWith(unquoteField(line), '#')) + [key, value] = parseHeaderLine(line); + if ~isempty(key) + meta.(key) = value; + end + line = fgetl(fid); +end +end + +function j = kcatSourceJson(v) +if isstruct(v) && isscalar(v) + f = sort(fieldnames(v)); + parts = cell(1, numel(f)); + for k = 1:numel(f) + parts{k} = [jsonString(f{k}) ':' kcatSourceJson(v.(f{k}))]; + end + j = ['{' strjoin(parts, ',') '}']; +elseif ischar(v) || (isstring(v) && isscalar(v)) + j = jsonString(char(v)); +elseif iscell(v) + j = ['[' strjoin(cellfun(@kcatSourceJson, v(:)', 'UniformOutput', false), ',') ']']; +elseif (isnumeric(v) || islogical(v)) && numel(v) ~= 1 + j = ['[' strjoin(arrayfun(@kcatSourceJson, v(:)', 'UniformOutput', false), ',') ']']; +elseif islogical(v) + j = 'false'; + if v + j = 'true'; + end +elseif isnumeric(v) + if ~isfinite(v) + j = 'null'; + else + j = formatKcatNumber(double(v)); + end +elseif isempty(v) + j = 'null'; +else + error('kcatSourceJson:type', 'Parameters of class %s cannot be written.', class(v)) +end +end + +function s = jsonString(s) +% JSON escapes, as Python's json.dumps(ensure_ascii=False). +out = ''; +for c = s + switch c + case '"', out = [out '\"']; %#ok + case '\', out = [out '\\']; %#ok + case char(8), out = [out '\b']; %#ok + case char(9), out = [out '\t']; %#ok + case char(10), out = [out '\n']; %#ok + case char(12), out = [out '\f']; %#ok + case char(13), out = [out '\r']; %#ok + otherwise + if double(c) < 32 + out = [out sprintf('\\u%04x', double(c))]; %#ok + else + out = [out c]; %#ok + end + end +end +s = ['"' out '"']; +end + +function t = formatKcatNumber(x) +if isnan(x) + t = ''; +elseif isinf(x) + t = 'Inf'; + if x < 0 + t = '-Inf'; + end +elseif x == round(x) && abs(x) < 1e15 + t = sprintf('%d', x); +else + for digits = 1:17 + t = sprintf('%.*g', digits, x); + if str2double(t) == x + break + end + end +end +end + +function [key, value] = parseHeaderLine(line) +key = ''; value = ''; +tok = regexp(unquoteField(line), '^\s*#\s*([A-Za-z][A-Za-z0-9_]{0,62})\s*:\s?(.*)$', 'tokens', 'once'); +if ~isempty(tok) + key = lower(tok{1}); + value = strtrim(tok{2}); +end +end + +function t = unquoteField(t) +t = strtrim(t); +if numel(t) >= 2 && t(1) == '"' && t(end) == '"' + t = strtrim(strrep(t(2:end-1), '""', '"')); +end +end + +function checkOptions(args, names, functionName) +% Options are name-value pairs with known names, so a misspelled option is +% an error instead of being taken as a positional argument. +id = [functionName ':arguments']; +given = args(1:2:end); +if mod(numel(args), 2) ~= 0 || ~all(cellfun(@(x) ischar(x) || (isstring(x) && isscalar(x)), given)) + error(id, 'Options of %s are name-value pairs. Options: %s.', functionName, strjoin(names, ', ')) +end +unknown = setdiff(cellfun(@char, given, 'UniformOutput', false), names); +if ~isempty(unknown) + error(id, 'Unknown option of %s: %s. Options: %s.', functionName, strjoin(unknown, ', '), ... + strjoin(names, ', ')) +end +end diff --git a/test/VerificationMatrix.tsv b/test/VerificationMatrix.tsv index 36c28481..f6df8d6b 100644 --- a/test/VerificationMatrix.tsv +++ b/test/VerificationMatrix.tsv @@ -17,7 +17,6 @@ gather_kcats/assignKcatValues.m tc0011 - we did not test all possible parameteri gather_kcats/writeDLKcatInput.m partly by tc0011 gather_kcats/getEcRxnKeys.m tc0070, tc0076 gather_kcats/saveKcatSource.m tc0071, tc0073, tc0076, tc0078 -gather_kcats/readKcatSource.m tc0071, tc0077 gather_kcats/loadKcatSource.m tc0071, tc0072, tc0073, tc0076, tc0077, tc0078 gather_kcats/kcatsFromBrenda.m tc0074, tc0078 gather_kcats/kcatsFromDLKcat.m tc0074, tc0078 - reads an existing DLKcat output; running DLKcat needs Docker diff --git a/test/unit_tests/kcatSourceTests.m b/test/unit_tests/kcatSourceTests.m index 4228cd74..d40829f1 100644 --- a/test/unit_tests/kcatSourceTests.m +++ b/test/unit_tests/kcatSourceTests.m @@ -28,7 +28,7 @@ function testSaveLoadRoundTrip_tc0071(testCase) folder = tempname; file = saveKcatSource(full, list, 'brenda', 'folder', folder, ... 'parameters', struct('aggregation','max'), 'data', 'BRENDA release test'); - [records, meta] = readKcatSource(file); + [records, meta] = rowsOf(full, file); verifyEqual(testCase, meta.source, 'brenda') verifyEqual(testCase, meta.kind, 'gathered') verifyEqual(testCase, meta.data, 'BRENDA release test') @@ -59,7 +59,7 @@ function testFullFileOnLightModel_tc0072(testCase) [loaded, report] = loadKcatSource(light, file); verifyEmpty(testCase, report.unmatched) [baseL, enzL] = getEcRxnKeys(light); - [records] = readKcatSource(file); + records = rowsOf(light, file); for r = 1:numel(records.kcat) hit = strcmp(baseL, records.rxn{r}) & strcmp(enzL, records.enzymes{r}); verifyEqual(testCase, sort(loaded.rxns(loaded.line == records.line(r))), sort(light.ec.rxns(hit))) @@ -106,7 +106,7 @@ function testKcatsFromBrendaAndDLKcat_tc0074(testCase) copyfile(fullfile(adapter.params.path, 'data'), fullfile(project, 'data')); adapter.params.path = project; list = kcatsFromBrenda(full, 'aggregation', 'median', 'modelAdapter', adapter); - [records, meta] = readKcatSource(fullfile(project, 'data', 'kcats', 'brenda.tsv')); + [records, meta] = rowsOf(full, fullfile(project, 'data', 'kcats', 'brenda.tsv')); verifyEqual(testCase, meta.generator(1:15), 'kcatsFromBrenda') verifyEqual(testCase, meta.parameters, '{"aggregation":"median","force_wildcard_level":0}') verifyEqual(testCase, meta.data, 'BRENDA release test-fixture') @@ -121,7 +121,7 @@ function testKcatsFromBrendaAndDLKcat_tc0074(testCase) end fclose(fid); kcatsFromDLKcat(full, 'modelAdapter', adapter); - [records, meta] = readKcatSource(fullfile(project, 'data', 'kcats', 'dlkcat.tsv')); + [records, meta] = rowsOf(full, fullfile(project, 'data', 'kcats', 'dlkcat.tsv')); verifyEqual(testCase, meta.source, 'dlkcat') verifyEqual(testCase, sort(records.kcat), [3.25; 12.5]) verifyEqual(testCase, unique(records.source), {'DLKcat'}) @@ -144,7 +144,7 @@ function testKcatsFromOKPStored_tc0075(testCase) writeLines(fullfile(data, 'OKP_job.txt'), {'jobId: abc123', 'method: CataPro'}); [list, done] = kcatsFromOKP(full, 'method', "CataPro", 'useStored', true, 'modelAdapter', adapter); verifyTrue(testCase, done) - [records, meta] = readKcatSource(fullfile(data, 'kcats', 'okp_catapro.tsv')); + [records, meta] = rowsOf(full, fullfile(data, 'kcats', 'okp_catapro.tsv')); verifyEqual(testCase, meta.data, 'OKP job abc123') verifyEqual(testCase, numel(records.kcat), sum(list.kcats > 0)) verifyEqual(testCase, unique(records.source), {'OKP-CataPro'}) @@ -193,7 +193,7 @@ function testGuardsNumbersAndJson_tc0078(testCase) file = saveKcatSource(full, list, "numbers", 'modelAdapter', adapter, 'parameters', ... struct('z', 1, 'a', struct('y', ['caf' char(233)], 'b', true), 'n', NaN, ... 'f', 1/3, 'g', 1e-7, 'h', 1e15, 'k', 1e20, 'l', {{5}}, 's', 'tab quote"')); - [records, meta] = readKcatSource(file); + [records, meta] = rowsOf(full, file); verifyTrue(testCase, all(isfinite(records.kcat) & records.kcat > 0)) % The same text as geckopy's test_parameters_json_matches_gecko. verifyEqual(testCase, meta.parameters, ['{"a":{"b":true,"y":"caf' char(233) '"},"f":0.3333333333333333,' ... @@ -208,7 +208,7 @@ function testGuardsNumbersAndJson_tc0078(testCase) 'R3 P4 m1 2 b ', 'R3 P4 m1 2 a wc=1', 'R3 P4 5 a ', 'R2_REV P3 m2 1 a '}); l2 = loadKcatSource(full, fullfile(folder, 'unsorted.tsv')); file = saveKcatSource(full, l2, 'sorted', 'modelAdapter', adapter); - records = readKcatSource(file); + records = rowsOf(full, file); verifyEqual(testCase, records.rxn', {'R2_REV','R3','R3','R3'}) verifyEqual(testCase, records.substrate', {'m2','','m1','m1'}) verifyEqual(testCase, records.source', {'a','a','a','b'}) @@ -238,7 +238,8 @@ function testSharedFixture_tc0077(testCase) % a quoted header line and field, comments and blank lines between % rows, a row of tabs, a decimal comma, a negative kcat. file = fullfile(findGECKOroot, 'test', 'unit_tests', 'ecTestGEM', 'data', 'kcatSourceFixture.tsv'); - [records, meta] = readKcatSource(file); + [full] = buildKcatModels(); + [records, meta] = rowsOf(full, file); verifyEqual(testCase, meta.source, 'fixture') verifyEqual(testCase, meta.parameters, '{"aggregation":"median"}') verifyEqual(testCase, meta.kind, 'curated') @@ -248,7 +249,6 @@ function testSharedFixture_tc0077(testCase) verifyEqual(testCase, records.kcat', [12.5 1e-05 NaN 7 3 -1 4]) verifyEqual(testCase, records.log10_sd(1), 0.3) - [full] = buildKcatModels(); [loaded, report] = loadKcatSource(full, file); verifyEqual(testCase, loaded.rxns(loaded.line == 7), {'R2_EXP_1'}) verifyEqual(testCase, sort(loaded.rxns(loaded.line == 12)), {'R2_REV_EXP_1';'R2_REV_EXP_2'}) @@ -259,7 +259,14 @@ function testSharedFixture_tc0077(testCase) % A row with more filled fields than columns is an error. bad = [tempname '.tsv']; writeLines(bad, {'rxn kcat source', 'R2 1 brenda extra'}); - verifyError(testCase, @() readKcatSource(bad), 'readKcatSource:fields') + verifyError(testCase, @() loadKcatSource(full, bad), 'loadKcatSource:fields') +end + +function [records, header] = rowsOf(model, file) + % The rows and header of a kcat source file, as loadKcatSource reads them. + [~, report] = loadKcatSource(model, file); + records = report.records; + header = report.header; end function writeLines(file, lines)