Skip to content

Commit 462fa83

Browse files
rmtflemingrmtfleming
authored
Merge pull request #2219 from farid-zare/debug-generateMetaboScore
Debug generateMetaboScore.m
2 parents 788cd18 + 0e3e00d commit 462fa83

File tree

1 file changed

+45
-25
lines changed

1 file changed

+45
-25
lines changed

‎src/reconstruction/metaboRePort/generateMetaboScore.m

Lines changed: 45 additions & 25 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,26 @@
11
function [modelProp,ScoresOverall] = generateMetaboScore(model,nworkers)
2-
3-
% Ines Thiele June 2022
2+
% Analyzes metabolic model properties and generates scores.
3+
% This function analyzes various aspects of the metabolic model,
4+
% including basic properties, stoichiometric consistency, matrix conditioning,
5+
% annotation of metabolites, reactions, genes, SBO terms, and generates
6+
% overall scores for different aspects.
7+
%
8+
% USAGE:
9+
% [modelProp, ScoresOverall] = generateMetaboScore(model, nworkers)
10+
%
11+
% INPUTS:
12+
% - model: metabolic model structure (COBRA Toolbox format)
13+
% - nworkers: (optional) number of workers for parallel processing
14+
%
15+
% OUTPUT:
16+
% - modelProp: structure containing various properties of the metabolic model
17+
% - ScoresOverall: structure containing overall scores for different aspects
18+
%
19+
% EXAMPLE:
20+
% [modelProp, ScoresOverall] = generateMetaboScore(model)
21+
% [modelProp, ScoresOverall] = generateMetaboScore(model, 4)
22+
%
23+
% .. Author: -Ines Thiele June 2022
424

525
if ~exist('nworkers','var')
626
nworkers = 4;
@@ -95,14 +115,14 @@
95115
modelProp.Details.SinkRxns = [model.rxns(contains(model.rxns,'Sink_'));model.rxns(contains(model.rxns,'sink_'))];
96116
% biomass reactions
97117
modelProp.BiomassRxns = BioR;
98-
modelProp.Details.BiomassRxns = model.rxns(contains(lower(model.rxns,'biomass')));% exclude EX_biomass?
118+
modelProp.Details.BiomassRxns = model.rxns(contains(lower(model.rxns),'biomass'));% exclude EX_biomass?
99119

100120
modelProp.MetabolicRxns = MetR;
101121
modelProp.Details.MetabolicRxns =MetRxns';
102122
modelProp.TransportRxns = TransR;
103123
modelProp.Details.TransportRxns =TransRxns';
104124

105-
% internal reactions without GPR
125+
% internal reactions without GPR
106126
RxnsWOGpr = model.rxns(find(cellfun(@isempty,model.grRules)));
107127
External = [modelProp.Details.ExchangeRxns;modelProp.Details.DemandRxns;modelProp.Details.SinkRxns ;modelProp.Details.BiomassRxns ];
108128
RxnsWOGpr = setdiff(RxnsWOGpr,External);
@@ -257,31 +277,31 @@
257277
% remove met from variable metWOAnno (metabolite without
258278
% annotation)
259279
metWOAnno = (intersect(metWOAnno,missingMet)) ;
260-
280+
261281
% get number of metabolites with annotation that conform with
262282
% database id format
263283
% test format
264284
clear confF* missing
265285
if contains(fields(i,2),'###')
266286
strs = strsplit(fields{i,2},'###');
267-
287+
268288
confFormat1 = model.(fields{i,1})(find(cellfun(@(x)~isempty(x),regexp(model.(fields{i,1}), strs{1}))));
269289
confFormat2 = model.(fields{i,1})(find(cellfun(@(x)~isempty(x),regexp(model.(fields{i,1}), strs{2}))));
270290
confFormat = unique([confFormat1;confFormat2]);
271-
291+
272292
confFormatM1 = model.mets(find(cellfun(@(x)~isempty(x),regexp(model.(fields{i,1}), strs{1}))));
273293
confFormatM2 = model.mets(find(cellfun(@(x)~isempty(x),regexp(model.(fields{i,1}), strs{2}))));
274294
confFormatM = unique([confFormatM1;confFormatM2]);
275-
295+
276296
% remove any potential NaN
277297
confFormat(ismember(confFormat,'NaN'))=[];
278298
confFormatM(ismember(confFormat,'NaN'))=[];
279-
299+
280300
confFormatM = split(confFormatM,'[');
281301
confFormatM = unique(confFormatM(:,1));
282302
modelProp.(strcat('AnnoMetConf', fields{i}))= (length(confFormatM))*100/(modelProp.metUnique - length(missingMet)); % how many have it
283303
p = setdiff(modelProp.Details.metabolites_unique, missingMet);
284-
304+
285305
NonConf = setdiff(p,confFormatM);
286306
modelProp.Details.(strcat('AnnoMetNonConf', fields{i})) = NonConf;
287307
else
@@ -290,12 +310,12 @@
290310
% remove any potential NaN
291311
confFormatM(ismember(confFormat,'NaN'))=[];
292312
confFormat(ismember(confFormat,'NaN'))=[];
293-
313+
294314
confFormatM = split(confFormatM,'[');
295315
confFormatM = unique(confFormatM(:,1));
296316
modelProp.(strcat('AnnoMetConf', fields{i}))= (length(confFormatM))*100/(modelProp.metUnique - length(missingMet)); % how many have it
297317
p = setdiff(modelProp.Details.metabolites_unique, missingMet);
298-
318+
299319
NonConf = setdiff(p,confFormatM);
300320
modelProp.Details.(strcat('AnnoMetNonConf', fields{i})) = NonConf;
301321
end
@@ -304,8 +324,8 @@
304324
modelProp.(strcat('AnnoMet' ,fields{i})) = 0; % how many have it
305325
modelProp.(strcat('AnnoMetConf', fields{i})) = 0; % how many have it
306326
modelProp.Details.(strcat('AnnoMetNonConf', fields{i})) = {};
307-
308-
327+
328+
309329
end
310330
end
311331
% Presence of Metabolite Annotation
@@ -345,13 +365,13 @@
345365
for i = 1 : size(fields)
346366
if isfield(model,fields{i,1})
347367
missingRxn = model.rxns(cellfun('isempty', model.(fields{i,1})));
348-
368+
349369
modelProp.Details.(strcat('missing', fields{i})) = missingRxn;
350370
modelProp.(strcat('AnnoRxn', fields{i}))= (modelProp.n - length(missingRxn))*100/modelProp.n; % how many have it
351371
% remove rxn from variable rxnWOAnno (reaction without
352372
% annotation)
353373
rxnWOAnno = (intersect(rxnWOAnno,missingRxn)) ;
354-
374+
355375
% get number of reactions with annotation that conform with
356376
% database id format
357377
% test format
@@ -361,9 +381,9 @@
361381
confFormat2 = model.rxns(find(cellfun(@(x)~isempty(x),regexp(model.(fields{i,1}), strs{2}))));
362382
confFormat3 = model.rxns(find(cellfun(@(x)~isempty(x),regexp(model.(fields{i,1}), strs{3}))));
363383
confFormat = [confFormat1;confFormat2;confFormat3];
364-
384+
365385
modelProp.(strcat('AnnoRxnConf', fields{i}))= (length(confFormat))*100/(modelProp.n - length(missingRxn)); % how many have it
366-
386+
367387
p = setdiff(model.rxns, missingRxn);
368388
NonConf = setdiff(p,confFormat);
369389
modelProp.Details.(strcat('AnnoRxnNonConf', fields{i})) = NonConf;
@@ -373,14 +393,14 @@
373393
p = setdiff(model.rxns, missingRxn);
374394
NonConf = setdiff(p,confFormat);
375395
modelProp.Details.(strcat('AnnoRxnNonConf', fields{i})) = NonConf;
376-
396+
377397
end
378398
else
379399
modelProp.Details.(strcat('missing' ,fields{i})) = model.rxns;
380400
modelProp.(strcat('AnnoRxn' ,fields{i})) = 0; % how many have it
381401
modelProp.(strcat('AnnoRxnConf', fields{i})) = 0; % how many have it
382402
modelProp.Details.(strcat('AnnoRxnNonConf', fields{i})) = {};
383-
403+
384404
end
385405
%
386406
end
@@ -425,12 +445,12 @@
425445
for i = 1 : size(fields)
426446
if isfield(model,fields{i,1})
427447
missingGene = model.genes(cellfun('isempty', model.(fields{i,1})));
428-
448+
429449
modelProp.Details.(strcat('missing', fields{i})) = missingGene;
430450
modelProp.(strcat('AnnoGene', fields{i}))= (length(model.genes) - length(missingGene))*100/length(model.genes); % how many have it
431451
% remove rxn from variable geneWOAnno (reaction without
432452
% annotation)
433-
453+
434454
% get number of reactions with annotation that conform with
435455
% database id format
436456
% test format
@@ -449,14 +469,14 @@
449469
confFormat = [confFormat1;confFormat2;confFormat3;confFormat4];
450470
end
451471
modelProp.(strcat('AnnoGeneConf', fields{i}))= (length(confFormat))*100/(length(model.genes) - length(missingGene)); % how many have it
452-
472+
453473
p = setdiff(model.genes, missingGene);
454474
NonConf = setdiff(p,confFormat);
455475
modelProp.Details.(strcat('AnnoGeneNonConf', fields{i})) = NonConf;
456476
else
457477
confFormat = model.(fields{i,1})(find(cellfun(@(x)~isempty(x),regexp(model.(fields{i,1}), fields(i,2)))));
458478
modelProp.(strcat('AnnoGeneConf', fields{i}))= (length(confFormat))*100/(length(model.genes) - length(missingGene)); % how many have it
459-
479+
460480
p = setdiff(model.genes, missingGene);
461481
NonConf = setdiff(p,confFormat);
462482
modelProp.Details.(strcat('AnnoGeneNonConf', fields{i})) = NonConf;
@@ -566,4 +586,4 @@
566586

567587
modelProp.Scores.Overall = (3* modelProp.Scores.Consistency + modelProp.Scores.AnnotationMetabolites + modelProp.Scores.AnnotationReactions + ...
568588
modelProp.Scores.AnnotationGenes + 2* modelProp.Scores.AnnotationSBO)*100/(3*100 + 100 +100 +100 + 2*100);
569-
ScoresOverall = modelProp.Scores.Overall;
589+
ScoresOverall = modelProp.Scores.Overall;

0 commit comments

Comments
 (0)