diff --git a/.gitattributes b/.gitattributes index a1f0cebc..012b630f 100644 --- a/.gitattributes +++ b/.gitattributes @@ -3,7 +3,8 @@ # Use MATLAB style diff *.m text diff=matlab # Mark as binary to prevent EOL changes. Also no need to track -mafft-linux*/* binary +mafft.bat binary +mafft/ binary # Files/folders that should not be included when downloading repo as ZIP .gitattributes export-ignore .gitignore export-ignore diff --git a/.github/workflows/testing.yml b/.github/workflows/testing.yml index 888c2aba..cad0bcf5 100644 --- a/.github/workflows/testing.yml +++ b/.github/workflows/testing.yml @@ -7,22 +7,24 @@ jobs: runs-on: self-hosted steps: + - name: Checkout + uses: actions/checkout@v2 - - name: Checkout - uses: actions/checkout@v2 + - name: Run tests + id: matlab-test + run: | + TEST_RESULTS=$(/usr/local/bin/matlab -nodisplay -nosplash -nodesktop -r "addpath(genpath('.')); cd('testing/unit_tests'); runtests(struct2table(dir('*.m')).name); exit;") + PARSED_RESULTS=$(echo $TEST_RESULTS | awk -F'com\\.' '{ n = split($2, v, "(___+)|(===+)|(---+)"); for (i = 0; ++i <= n;) { print v[i] } }') + PARSED_RESULTS="${PARSED_RESULTS//'%'/'%25'}" + PARSED_RESULTS="${PARSED_RESULTS//$'\n'/'%0A'}" + PARSED_RESULTS="${PARSED_RESULTS//$'\r'/'%0D'}" + echo ::set-output name=results::$PARSED_RESULTS - - name: Run tests - id: matlab-test - run: | - TEST_RESULTS=$(/usr/local/bin/matlab -nodisplay -nosplash -nodesktop -r "addpath(genpath('.')); cd('testing/unit_tests'); runtests(struct2table(dir('*.m')).name); exit;") - PARSED_RESULTS=$(echo $TEST_RESULTS | awk -F'.com.' '{ n = split($2, v, "__________"); for (i = 0; ++i <= n;) { print v[i] } }') - echo ::set-output name=results::$PARSED_RESULTS - - - name: Post comment - uses: NejcZdovc/comment-pr@v1 - with: - file: "testing-comment.md" - env: - GITHUB_TOKEN: ${{secrets.GITHUB_TOKEN}} - TEST_RESULTS: ${{steps.matlab-test.outputs.results}} - GH_ACTION_RUN: ${{github.run_id}} \ No newline at end of file + - name: Post comment + uses: NejcZdovc/comment-pr@v1 + with: + file: "testing-comment.md" + env: + GITHUB_TOKEN: ${{secrets.GITHUB_TOKEN}} + TEST_RESULTS: ${{steps.matlab-test.outputs.results}} + GH_ACTION_RUN: ${{github.run_id}} diff --git a/INIT/getINITModel.m b/INIT/getINITModel.m index 1a944403..c4f54a13 100755 --- a/INIT/getINITModel.m +++ b/INIT/getINITModel.m @@ -109,6 +109,8 @@ if nargin<3 celltype=[]; +else + celltype=char(celltype); end if nargin<4 hpaData=[]; @@ -121,6 +123,8 @@ end if nargin<7 taskFile=[]; +else + taskFile=char(taskFile); end if nargin<8 || isempty(useScoresForTasks) useScoresForTasks=true; diff --git a/INIT/runINIT.m b/INIT/runINIT.m index 3f71bf03..b35b180e 100755 --- a/INIT/runINIT.m +++ b/INIT/runINIT.m @@ -67,24 +67,17 @@ if isempty(rxnScores) rxnScores=zeros(numel(model.rxns),1); end -if nargin<3 +if nargin<3 || isempty(presentMets) presentMets={}; +else + presentMets=convertCharArray(presentMets); end -if isempty(presentMets) - presentMets={}; -end -presentMets=presentMets(:); -if nargin<4 - essentialRxns={}; -end -if isempty(essentialRxns) +if nargin<4 || isempty(essentialRxns) essentialRxns={}; +else + essentialRxns=convertCharArray(essentialRxns); end -essentialRxns=essentialRxns(:); -if nargin<5 - prodWeight=0.5; -end -if isempty(prodWeight) +if nargin<5 || isempty(prodWeight) prodWeight=0.5; end if nargin<6 diff --git a/core/FSEOF.m b/core/FSEOF.m index fc5abbc3..cc575e7e 100755 --- a/core/FSEOF.m +++ b/core/FSEOF.m @@ -23,6 +23,9 @@ % % Usage: targets=FSEOF(model,biomassRxn,targetRxn,iterations,coefficient,outputFile) +biomassRxn=char(biomassRxn); +targetRxn=char(targetRxn); + if nargin<4 iterations=10; coefficient=0.9; @@ -96,6 +99,7 @@ %Generating output formatSpec='%s\t%s\t%s\t%s\t%s\t%s\t%s\n'; if output == 1 %Output to a file + outputFile=char(outputFile); fid=fopen(outputFile,'w'); fprintf(fid,formatSpec,'Slope','rowID','Enzyme ID','Enzyme Name','Subsystems','Direction','Gr Rule'); else %Output to screen diff --git a/core/addExchangeRxns.m b/core/addExchangeRxns.m index 9fc40c2e..660b08bf 100755 --- a/core/addExchangeRxns.m +++ b/core/addExchangeRxns.m @@ -23,11 +23,14 @@ if nargin<3 mets=model.mets; +elseif ~islogical(mets) && ~isnumeric(mets) + mets=convertCharArray(mets); end -reactionType=upper(reactionType); J=getIndexes(model,mets,'mets',false); mets=model.mets(J); +reactionType=char(upper(reactionType)); + %Production is positive for OUT and BOTH if strcmp(reactionType,'IN') I=speye(numel(model.mets)); diff --git a/core/addGenesRaven.m b/core/addGenesRaven.m index 26059125..b853ce36 100755 --- a/core/addGenesRaven.m +++ b/core/addGenesRaven.m @@ -32,11 +32,8 @@ if ~isfield(genesToAdd,'genes') EM='genes is a required field in genesToAdd'; dispEM(EM); -end - -if ~iscellstr(genesToAdd.genes) - EM='genesToAdd.genes must be a cell array of strings'; - dispEM(EM); +else + genesToAdd.genes=convertCharArray(genesToAdd.genes); end %Number of genes @@ -68,14 +65,11 @@ %Some more checks and if they pass then add each field to the structure if isfield(genesToAdd,'geneShortNames') + genesToAdd.geneShortNames=convertCharArray(genesToAdd.geneShortNames); if numel(genesToAdd.geneShortNames)~=nGenes EM='genesToAdd.geneShortNames must have the same number of elements as genesToAdd.genes'; dispEM(EM); end - if ~iscellstr(genesToAdd.geneShortNames) - EM='genesToAdd.geneShortNames must be a cell array of strings'; - dispEM(EM); - end %Add empty field if it doesn't exist if ~isfield(newModel,'geneShortNames') newModel.geneShortNames=largeFiller; diff --git a/core/addMets.m b/core/addMets.m index 303ddd23..223776f8 100755 --- a/core/addMets.m +++ b/core/addMets.m @@ -58,6 +58,8 @@ end if nargin<4 prefix='m_'; +else + prefix=char(prefix); end newModel=model; @@ -69,35 +71,19 @@ %Check some stuff regarding the required fields if ~isfield(metsToAdd,'mets') metsToAdd.mets=generateNewIds(newModel,'mets',prefix,numel(metsToAdd.metNames)); +else + metsToAdd.mets=convertCharArray(metsToAdd.mets); end if ~isfield(metsToAdd,'metNames') metsToAdd.metNames=metsToAdd.mets; +else + metsToAdd.metNames=convertCharArray(metsToAdd.metNames); end if ~isfield(metsToAdd,'compartments') EM='compartments is a required field in metsToAdd'; dispEM(EM); -end -if ischar(metsToAdd.mets) - metsToAdd.mets={metsToAdd.mets}; -elseif ~iscellstr(metsToAdd.mets) - EM='metsToAdd.mets must be a cell array of strings'; - dispEM(EM); -end -if ischar(metsToAdd.metNames) - metsToAdd.metNames={metsToAdd.metNames}; -elseif ~iscellstr(metsToAdd.metNames) - EM='metsToAdd.metNames must be a cell array of strings'; - dispEM(EM); -end -if ~iscellstr(metsToAdd.compartments) - if ischar(metsToAdd.compartments) - temp=cell(numel(metsToAdd.mets),1); - temp(:)={metsToAdd.compartments}; - metsToAdd.compartments=temp; - else - EM='metsToAdd.compartments must be a cell array of strings'; - dispEM(EM); - end +else + metsToAdd.compartments=convertCharArray(metsToAdd.compartments); end %Number of metabolites @@ -193,14 +179,11 @@ end if isfield(metsToAdd,'inchis') + metsToAdd.inchis=convertCharArray(metsToAdd.inchis); if numel(metsToAdd.inchis)~=nMets EM='metsToAdd.inchis must have the same number of elements as metsToAdd.mets'; dispEM(EM); end - if ~iscellstr(metsToAdd.inchis) - EM='metsToAdd.inchis must be a cell array of strings'; - dispEM(EM); - end %Add empty field if it doesn't exist if ~isfield(newModel,'inchis') newModel.inchis=largeFiller; @@ -214,14 +197,11 @@ end if isfield(metsToAdd,'metFormulas') + metsToAdd.metFormulas=convertCharArray(metsToAdd.metFormulas); if numel(metsToAdd.metFormulas)~=nMets EM='metsToAdd.metFormulas must have the same number of elements as metsToAdd.mets'; dispEM(EM); end - if ~iscellstr(metsToAdd.metFormulas) - EM='metsToAdd.metFormulas must be a cell array of strings'; - dispEM(EM); - end %Add empty field if it doesn't exist if ~isfield(newModel,'metFormulas') newModel.metFormulas=largeFiller; @@ -243,6 +223,9 @@ EM='metsToAdd.metCharges must be of type "double"'; dispEM(EM); end + if ~isfield(newModel,'metCharges') + newModel.metCharges=NaN(numel(largeFiller),1); + end newModel.metCharges=[newModel.metCharges;metsToAdd.metCharges(:)]; else %Add default @@ -252,19 +235,16 @@ end if isfield(metsToAdd,'metNotes') + metsToAdd.metNotes=convertCharArray(metsToAdd.metNotes); if numel(metsToAdd.metNotes)~=nMets EM='metsToAdd.metNotes must have the same number of elements as metsToAdd.mets'; dispEM(EM); end - if ~iscellstr(metsToAdd.metNotes) - EM='metsToAdd.metNotes must be a cell array of strings'; - dispEM(EM); - end %Add empty field if it doesn't exist if ~isfield(newModel,'metNotes') newModel.metNotes=largeFiller; end - newModel.metNotes=[newModel.inchis;metsToAdd.metNotes(:)]; + newModel.metNotes=[newModel.metNotes;metsToAdd.metNotes(:)]; else %Add empty strings if structure is in model if isfield(newModel,'metNotes') @@ -326,16 +306,3 @@ end end end - -%For getting the numerical form of metabolite ids on the form "m1". -function I=getInteger(s) -%Checks if a string is on the form "m1" and if so returns the value of the -%integer -I=0; -if strcmpi(s(1),'m') - t=str2double(s(2:end)); - if ~isnan(t) && ~isempty(t) - I=t; - end -end -end diff --git a/core/addRxns.m b/core/addRxns.m index 6240fe3b..05172dfa 100755 --- a/core/addRxns.m +++ b/core/addRxns.m @@ -103,16 +103,20 @@ if nargin<4 compartment=[]; +else + compartment=char(compartment); end if nargin<5 allowNewMets=false; +elseif ~islogical(allowNewMets) + allowNewMets=char(allowNewMets); end if nargin<6 allowNewGenes=false; end if allowNewGenes & isfield(rxnsToAdd,'grRules') - genesToAdd.genes = strjoin(rxnsToAdd.grRules); + genesToAdd.genes = strjoin(convertCharArray(rxnsToAdd.grRules)); genesToAdd.genes = regexp(genesToAdd.genes,' |)|(|and|or','split'); % Remove all grRule punctuation genesToAdd.genes = genesToAdd.genes(~cellfun(@isempty,genesToAdd.genes)); % Remove spaces and empty genes genesToAdd.genes = setdiff(unique(genesToAdd.genes),model.genes); % Only keep new genes @@ -147,10 +151,7 @@ end if eqnType==2 || (eqnType==1 && allowNewMets==true) - if ~ischar(compartment) - EM='compartment must be a string'; - dispEM(EM); - end + compartment=char(compartment); if ~ismember(compartment,model.comps) EM='compartment must match one of the compartments in model.comps'; dispEM(EM); @@ -160,7 +161,8 @@ if ~isfield(rxnsToAdd,'rxns') EM='rxns is a required field in rxnsToAdd'; dispEM(EM); -elseif iscell(rxnsToAdd.rxns) +else + rxnsToAdd.rxns=convertCharArray(rxnsToAdd.rxns); %To fit with some later printing rxnsToAdd.rxns=rxnsToAdd.rxns(:); end @@ -174,30 +176,16 @@ dispEM(EM); end -if ~iscellstr(rxnsToAdd.rxns) && ~ischar(rxnsToAdd.rxns) - %It could also be a string, but it's not encouraged - EM='rxnsToAdd.rxns must be a cell array of strings'; - dispEM(EM); -else - rxnsToAdd.rxns=cellstr(rxnsToAdd.rxns); -end - %Normal case: equations provided if isfield(rxnsToAdd,'equations') - if ~iscellstr(rxnsToAdd.equations) && ~ischar(rxnsToAdd.equations) - %It could also be a string, but it's not encouraged - EM='rxnsToAdd.equations must be a cell array of strings'; - dispEM(EM); - else - rxnsToAdd.equations=cellstr(rxnsToAdd.equations); - end - + rxnsToAdd.equations=convertCharArray(rxnsToAdd.equations); + %Alternative case: mets+stoichiometry provided else %In the case of 1 rxn added (array of strings + vector), transform to %cells of length=1: if iscellstr(rxnsToAdd.mets) - rxnsToAdd.mets = {rxnsToAdd.mets}; + rxnsToAdd.mets={rxnsToAdd.mets}; end if isnumeric(rxnsToAdd.stoichCoeffs) rxnsToAdd.stoichCoeffs = {rxnsToAdd.stoichCoeffs}; @@ -257,6 +245,7 @@ newModel.rxns=[newModel.rxns;rxnsToAdd.rxns(:)]; if isfield(rxnsToAdd,'rxnNames') + rxnsToAdd.rxnNames=convertCharArray(rxnsToAdd.rxnNames); if numel(rxnsToAdd.rxnNames)~=nRxns EM='rxnsToAdd.rxnNames must have the same number of elements as rxnsToAdd.rxns'; dispEM(EM); @@ -312,7 +301,7 @@ end %Fill with standard if it doesn't exist if ~isfield(newModel,'ub') - newModel.ub=repmat(newUb,nOldrxns,1); + newModel.ub=repmat(newUb,nOldRxns,1); end newModel.ub=[newModel.ub;rxnsToAdd.ub(:)]; else @@ -340,6 +329,7 @@ end if isfield(rxnsToAdd,'eccodes') + rxnsToAdd.eccodes=convertCharArray(rxnsToAdd.eccodes); if numel(rxnsToAdd.eccodes)~=nRxns EM='rxnsToAdd.eccodes must have the same number of elements as rxnsToAdd.rxns'; dispEM(EM); @@ -414,6 +404,7 @@ end end if isfield(rxnsToAdd,'grRules') + rxnsToAdd.grRules=convertCharArray(rxnsToAdd.grRules); if numel(rxnsToAdd.grRules)~=nRxns EM='rxnsToAdd.grRules must have the same number of elements as rxnsToAdd.rxns'; dispEM(EM); @@ -431,6 +422,7 @@ end if isfield(rxnsToAdd,'rxnFrom') + rxnsToAdd.rxnFrom=convertCharArray(rxnsToAdd.rxnFrom); if numel(rxnsToAdd.rxnFrom)~=nRxns EM='rxnsToAdd.rxnFrom must have the same number of elements as rxnsToAdd.rxns'; dispEM(EM); @@ -448,6 +440,7 @@ end if isfield(rxnsToAdd,'rxnNotes') + rxnsToAdd.rxnNotes=convertCharArray(rxnsToAdd.rxnNotes); if numel(rxnsToAdd.rxnNotes)~=nRxns EM='rxnsToAdd.rxnNotes must have the same number of elements as rxnsToAdd.rxns'; dispEM(EM); @@ -465,6 +458,7 @@ end if isfield(rxnsToAdd,'rxnReferences') + rxnsToAdd.rxnReferences=convertCharArray(rxnsToAdd.rxnReferences); if numel(rxnsToAdd.rxnReferences)~=nRxns EM='rxnsToAdd.rxnReferences must have the same number of elements as rxnsToAdd.rxns'; dispEM(EM); @@ -482,6 +476,7 @@ end if isfield(rxnsToAdd,'pwys') + rxnsToAdd.pwys=convertCharArray(rxnsToAdd.pwys); if numel(rxnsToAdd.pwys)~=nRxns EM='rxnsToAdd.pwys must have the same number of elements as rxnsToAdd.rxns'; dispEM(EM); diff --git a/core/addRxnsGenesMets.m b/core/addRxnsGenesMets.m index cdada989..ff452e1a 100755 --- a/core/addRxnsGenesMets.m +++ b/core/addRxnsGenesMets.m @@ -5,8 +5,8 @@ % % model draft model where reactions should be copied to % sourceModel model where reactions and metabolites are sourced from -% rxns cell array with reaction IDs (from source model) -% string allowed if only one reaction is added +% rxns cell array with reaction IDs (from source model). Can also +% be string if only one reaction is added % addGene three options: % false no genes are annotated to the new reactions % true grRules ared copied from the sourceModel and @@ -16,10 +16,12 @@ % array, and any new genes are added when % required % (opt, default false) -% rxnNote string explaining why reactions were copied to model, -% is included as newModel.rxnNotes (opt, default +% rxnNote cell array with strings explaining why reactions were copied +% to the model, to be included as newModel.rxnNotes. Can also +% be string if same rxnNotes should be added for each new +% reaction, or only one reaction is to be added (opt, default % 'Added via addRxnsAndMets()') -% confidence double specifying confidence score for all reactions. +% confidence integer specifying confidence score for all reactions. % 4: biochemical data: direct evidence from enzymes % assays % 3: genetic data: knockout/-in or overexpression @@ -45,19 +47,21 @@ if nargin<6 confidence=0; end +rxns=convertCharArray(rxns); if nargin<5 - rxnNote='Added via addRxnsGenesMets()'; + rxnNote={'Added via addRxnsGenesMets()'}; +else + rxnNote=convertCharArray(rxnNote); + if numel(rxnNote)==1 && numel(rxns)>1 + rxnNoteArray=cell(numel(rxns),1); + rxnNoteArray{:}=rxnNote; + rxnNote=rxnNoteArray; + end end if nargin<4 addGene=false; end -%If the supplied object is a character array, then convert it to a cell -%array -if ischar(rxns) - rxns={rxns}; -end - % Obtain indexes of reactions in source model [notNewRxn,oldRxn]=ismember(rxns,model.rxns); rxns=rxns(~notNewRxn); @@ -150,8 +154,7 @@ rxnToAdd.rxns=sourceModel.rxns(rxnIdx); rxnToAdd.lb=sourceModel.lb(rxnIdx); rxnToAdd.ub=sourceModel.ub(rxnIdx); -rxnToAdd.rxnNotes=cell(1,numel(rxnToAdd.rxns)); -rxnToAdd.rxnNotes(:)={rxnNote}; +rxnToAdd.rxnNotes(:)=rxnNote(~notNewRxn); rxnToAdd.rxnConfidenceScores=NaN(1,numel(rxnToAdd.rxns)); if ~isnumeric(confidence) EM='confidence score must be numeric'; diff --git a/core/addTransport.m b/core/addTransport.m index dcf9f69a..a97aa836 100755 --- a/core/addTransport.m +++ b/core/addTransport.m @@ -25,18 +25,14 @@ % Usage: [model, addedRxns]=addTransport(model,fromComp,toComps,metNames,... % isRev,onlyToExisting,prefix) -if iscell(fromComp) - fromComp=fromComp{1}; -end +fromComp=char(fromComp); [I, fromID]=ismember(model.comps,fromComp); fromID=find(fromID); if sum(I)~=1 EM='fromComps must have exactly one match in model.comps'; dispEM(EM); end -if ischar(toComps) - toComps={toComps}; -end +toComps=convertCharArray(toComps); [I, toIDs]=ismember(toComps,model.comps); if ~all(I) EM='All compartments in toComps must have a match in model.comps'; @@ -45,12 +41,11 @@ if nargin<4 %Find all metabolites in fromComp metNames=model.metNames(model.metComps==fromID); -end - -%If an empty set was given -if isempty(metNames) +elseif isempty(metNames) %Find all metabolites in fromComp metNames=model.metNames(ismember(model.metComps,model.comps(fromID))); +else + metNames=convertCharArray(metNames); end if nargin<5 @@ -61,12 +56,11 @@ end if nargin<7 prefix='tr_'; +else + prefix=char(prefix); end %Check that the names are unique -if ischar(metNames) - metNames={metNames}; -end if numel(unique(metNames))~=numel(metNames) dispEM('Not all metabolite names are unique'); end diff --git a/core/buildEquation.m b/core/buildEquation.m index e3741da7..a4a20563 100755 --- a/core/buildEquation.m +++ b/core/buildEquation.m @@ -2,7 +2,7 @@ % buildEquation % Construct single equation string for a given reaction % -% mets string array with metabolites involved in the reaction. +% mets cell array with metabolites involved in the reaction. % stoichCoeffs vector with corresponding stoichiometric coeffs. % isrev logical indicating if the reaction is or not % reversible. @@ -11,12 +11,7 @@ % % Usage: equationString=buildEquation(mets,stoichCoeffs,isrev) -if ~iscellstr(mets) && ~ischar(mets) - EM = 'mets must be a cell array of strings'; - dispEM(EM); -else - mets = cellstr(mets); -end +mets=convertCharArray(mets); if ~isnumeric(stoichCoeffs) EM = 'stoichCoeffs must be a numeric vector'; dispEM(EM); diff --git a/core/canConsume.m b/core/canConsume.m index a8da0f2b..931d0e21 100755 --- a/core/canConsume.m +++ b/core/canConsume.m @@ -15,6 +15,8 @@ if nargin<2 mets=model.mets; +elseif ~islogical(mets) && ~isnumeric(mets) + mets=convertCharArray(mets); end [model, rxns]=addExchangeRxns(model,'in',mets); diff --git a/core/canProduce.m b/core/canProduce.m index dda218f3..1472f5b1 100755 --- a/core/canProduce.m +++ b/core/canProduce.m @@ -16,6 +16,8 @@ if nargin<2 mets=model.mets; +elseif ~islogical(mets) && ~isnumeric(mets) + mets=convertCharArray(mets); end [model, rxns]=addExchangeRxns(model,'out',mets); diff --git a/core/changeGrRules.m b/core/changeGrRules.m index 33545e98..15ceea91 100755 --- a/core/changeGrRules.m +++ b/core/changeGrRules.m @@ -20,13 +20,8 @@ replace=true; end -if isstr(rxns) - rxns={rxns}; -end - -if isstr(grRules) - grRules={grRules}; -end +rxns=convertCharArray(rxns); +grRules=convertCharArray(grRules); if ~(numel(grRules)==numel(rxns)) error('Number of rxns and grRules should be identical') diff --git a/core/changeRxns.m b/core/changeRxns.m index d40b0ec0..e4b068cc 100755 --- a/core/changeRxns.m +++ b/core/changeRxns.m @@ -64,12 +64,9 @@ allowNewMets=false; end -if ischar(rxns) - rxns={rxns}; -end -if ischar(equations) - equations={equations}; -end +rxns=convertCharArray(rxns); +equations=convertCharArray(equations); +compartment=char(compartment); %Find the indexes of the reactions and throw an error if they aren't all %found diff --git a/core/checkProduction.m b/core/checkProduction.m index f8416868..1f49002c 100755 --- a/core/checkProduction.m +++ b/core/checkProduction.m @@ -48,6 +48,8 @@ if nargin<3 excretionFromCompartments=model.comps; +else + excretionFromCompartments=convertCharArray(excretionFromCompartments); end if nargin<4 diff --git a/core/checkRxn.m b/core/checkRxn.m index 3503ab3c..4632e400 100755 --- a/core/checkRxn.m +++ b/core/checkRxn.m @@ -5,7 +5,7 @@ % reactions which cannot have flux % % model a model structure -% rxn the id of the reaction to check +% rxn the id of one reaction to check % cutoff minimal flux for successful production/consumption (opt, % default 10^-7) % revDir true if the reaction should be reversed (opt, default @@ -22,10 +22,7 @@ % % Usage: report=checkRxn(model,rxn,cutoff,revDir,printReport) -%Convert to cell string -if ischar(rxn) - rxn={rxn}; -end +rxn=char(rxn); if nargin<3 cutoff=10^-7; end diff --git a/core/compareMultipleModels.m b/core/compareMultipleModels.m index 3e6ea8b5..a13d84f1 100755 --- a/core/compareMultipleModels.m +++ b/core/compareMultipleModels.m @@ -41,6 +41,11 @@ % Usage: compStruct=compareMultipleModels(models,printResults,... % plotResults,groupVector,funcCompare,taskFile); +%% Stats toolbox required +if ~(exist('mdscale.m','file') && exist('pdist.m','file') && exist('squareform.m','file') && exist('tsne.m','file')) + error('The MATLAB Statistics and Machine Learning Toolbox is required for this function') +end + %% Set up input defaults if nargin < 2 || isempty(printResults) printResults=false; @@ -62,6 +67,8 @@ end if nargin < 6 taskFile = []; +else + taskFile=char(taskFile); end if numel(models) <= 1 EM = 'Cannot compare only one model. Use printModelStats if you want a summary of a model'; @@ -182,16 +189,14 @@ proj_coords = tsne(double(binary_matrix'),'Distance','hamming','NumDimensions',3); % 3D compStruct.structCompMap = proj_coords; axis_labels = {'tSNE 1';'tSNE 2';'tSNE 3'}; - fprintf('*** Done \n\n') -else - fprintf('\nWARNING: Could not complete full structural comparison because the function \n') - fprintf(' "tsne" does not exist in your Matlab version. \n') - fprintf(' Using MDS to project data instead of tSNE. \n') - fprintf(' Please upgrade to Matlab 2017b or higher for full functionality. \n\n') +else % Seems odd to use mdscale if tsne is not found, as both are + % distributed with stats toolbox. If tsne is not present, then also + % mdscale is missing. Anyway, will leave this for now. [proj_coords,stress,disparities] = mdscale(pdist(double(binary_matrix'),'hamming'),3); compStruct.structCompMap = proj_coords; axis_labels = {'MDS 1';'MDS 2';'MDS 3'}; end +fprintf('*** Done \n\n') % plot structure comparison results if plotResults == true diff --git a/core/constructEquations.m b/core/constructEquations.m index 8acdb53a..8672ac8c 100755 --- a/core/constructEquations.m +++ b/core/constructEquations.m @@ -34,8 +34,10 @@ % Usage: equationStrings=constructEquations(model,rxns,useComps,... % sortRevRxns,sortMetNames,useMetID,useFormula,useRevField) -if nargin<2 +if nargin<2 || isempty(rxns) rxns=model.rxns; +elseif ~islogical(rxns) && ~isnumeric(rxns) + rxns=convertCharArray(rxns); end if nargin<3 useComps=true; @@ -55,9 +57,6 @@ if nargin<8 useRevField=true; end -if isempty(rxns) && nargin>2 - rxns=model.rxns; -end %Sort reversible equations if sortRevRxns==true diff --git a/core/constructS.m b/core/constructS.m index d6b5699c..6cf79943 100755 --- a/core/constructS.m +++ b/core/constructS.m @@ -22,6 +22,14 @@ % % Usage: [S, mets, badRxns, reversible]=constructS(equations,mets) +equations=convertCharArray(equations); +switch nargin + case 2 + mets=convertCharArray(mets); + case 3 + rxns=convertCharArray(rxns); +end + badRxns=false(numel(equations),1); %Check that no equations are too short to have reversibility data diff --git a/core/consumeSomething.m b/core/consumeSomething.m index e363b2bf..b7c4d181 100755 --- a/core/consumeSomething.m +++ b/core/consumeSomething.m @@ -42,6 +42,8 @@ if nargin<2 ignoreMets=[]; +elseif ~islogical(ignoreMets) && ~isnumeric(ignoreMets) + ignoreMets=convertCharArray(ignoreMets); end if nargin<3 isNames=false; diff --git a/core/convertCharArray.m b/core/convertCharArray.m new file mode 100644 index 00000000..40d97877 --- /dev/null +++ b/core/convertCharArray.m @@ -0,0 +1,27 @@ +function inputConverted = convertCharArray(funcInput) +%convertCharArray +% Converts input to make sure it is a cell array of character vectors. +% String arrays are transformed into character vectors, and if only one +% character vector is given. Output is always a cell array, also if only +% one character vector is given as input. +% +% Input: +% funcInput function input that should be checked +% +% Output: +% inputConverted cell array of character vectors +% + +if ~isempty(funcInput) + if ischar(funcInput) || isstring(funcInput) + inputConverted = {char(funcInput)}; + elseif ~iscell(funcInput) && ~(all(cellfun(@ischar, funcInput)) || all(cellfun(@isstring, funcInput))) + error([inputname(1) ' should be a (cell array of) character vector(s)']) + else + inputConverted = cellstr(funcInput); + end +else + inputConverted=[]; +end +end + diff --git a/core/convertToIrrev.m b/core/convertToIrrev.m index aae1e62f..5cb353d7 100755 --- a/core/convertToIrrev.m +++ b/core/convertToIrrev.m @@ -16,6 +16,7 @@ if nargin<2 I=true(numel(model.rxns),1); else + rxns=convertCharArray(rxns); I=getIndexes(model,rxns,'rxns',true); end diff --git a/core/copyToComps.m b/core/copyToComps.m index c5ae846d..fb54ecc4 100755 --- a/core/copyToComps.m +++ b/core/copyToComps.m @@ -27,16 +27,22 @@ if nargin<3 rxns=model.rxns; +elseif ~islogical(rxns) && ~isnumeric(rxns) + rxns=convertCharArray(rxns); end if nargin<4 deleteOriginal=false; end if nargin<5 compNames=toComps; +else + compNames=convertCharArray(compNames); end if nargin<6 compOutside=cell(numel(toComps),1); compOutside(:)={''}; +else + compOutside=convertCharArray(compOutside); end originalID=model.id; @@ -46,7 +52,7 @@ for i=1:numel(toComps) %Check if the compartment exists, otherwise add it - [I J]=ismember(toComps(i),model.comps); + [I,J]=ismember(toComps(i),model.comps); if I==false model.comps=[model.comps;toComps(i)]; model.compNames=[model.compNames;compNames(i)]; @@ -78,6 +84,10 @@ model=mergeModels({model;modelToAdd},'metNames'); end +model=rmfield(model,'rxnFrom'); +model=rmfield(model,'metFrom'); +model=rmfield(model,'geneFrom'); + if deleteOriginal==true model=removeReactions(model,rxns,true,true,true); %Also delete unused compartments end diff --git a/core/dispEM.m b/core/dispEM.m index 36f297d3..43d28183 100755 --- a/core/dispEM.m +++ b/core/dispEM.m @@ -21,8 +21,8 @@ function dispEM(string,throwErrors,toList,trimWarnings) toList=[]; elseif isempty(toList) return; -elseif ~iscell(toList) - toList={toList}; +else + toList=convertCharArray(toList); end if nargin<4 trimWarnings=true; diff --git a/core/findGeneDeletions.m b/core/findGeneDeletions.m index 96f6fbe8..8d049542 100755 --- a/core/findGeneDeletions.m +++ b/core/findGeneDeletions.m @@ -50,6 +50,8 @@ if nargin<5 oeFactor=10; end +testType=char(testType); +analysisType=char(analysisType); %Check that the test type is correct if ~strcmpi(testType,'sgd') && ~strcmpi(testType,'dgd') && ~strcmpi(testType,'sgo') && ~strcmpi(testType,'dgo') diff --git a/core/findRAVENroot.m b/core/findRAVENroot.m new file mode 100644 index 00000000..b8f6f3c0 --- /dev/null +++ b/core/findRAVENroot.m @@ -0,0 +1,27 @@ +function [ravenPath, prevDir] = findRAVENroot() +% findRAVENroot +% Finds the root of the RAVEN directory, by searching for the path to +% RAVEN2.png. Can also record the current directory, in case a function will +% use the ravenPath to navigate to a precise folder, and it should return to +% the previous directory afterwards. See e.g. optimizeProb calling glpk. + +ST=dbstack('-completenames'); +prevDir = pwd(); +if length(ST)>1 + ravenPath=ST(2).file; % In case findRAVENroot is run via another function +else + ravenPath=ST(1).file; +end +rootFound = 0; +while rootFound == 0 + isRoot = exist(fullfile(ravenPath,'RAVEN2.png'),'file'); + if isRoot == 2 + rootFound = 1; + else + ravenPathOld = ravenPath; + ravenPath = fileparts(ravenPath); + if strcmp(ravenPathOld,ravenPath) + error('Cannot find the RAVEN root directory. Make sure you have not removed the RAVEN2.png file from your RAVEN installation.') + end + end +end \ No newline at end of file diff --git a/core/fitParameters.m b/core/fitParameters.m index b25a203a..a72334cd 100755 --- a/core/fitParameters.m +++ b/core/fitParameters.m @@ -3,7 +3,8 @@ % Fits parameters such as maintenance ATP by quadratic programming % % model a model structure -% xRxns cell array with the IDs of the reactions that will be fixed for each data point +% xRxns cell array with the IDs of the reactions that will be +% fixed for each data point % xValues matrix with the corresponding values for each % xRxns (columns are reactions) % rxnsToFit cell array with the IDs of reactions that will be fitted to @@ -47,6 +48,9 @@ plotFitting=false; end +xRxns=convertCharArray(xRxns); +rxnsToFit=convertCharArray(rxnsToFit); + %Find the indexes of reactions that will be fitted [I, rxnsToFitIndexes]=ismember(rxnsToFit,model.rxns); diff --git a/core/followChanged.m b/core/followChanged.m index a8c9b4d5..fd64f1ac 100755 --- a/core/followChanged.m +++ b/core/followChanged.m @@ -32,6 +32,8 @@ function followChanged(model,fluxesA,fluxesB, cutOffChange, cutOffFlux, cutOffDi end if nargin<7 metaboliteList=[]; +else + metaboliteList=convertCharArray(metaboliteList); end %If a metabolite list is to be used, then find all the reactions involving diff --git a/core/generateNewIds.m b/core/generateNewIds.m index 58e4911f..ea58afd5 100755 --- a/core/generateNewIds.m +++ b/core/generateNewIds.m @@ -16,6 +16,8 @@ % % Usage: newIds=generateNewIds(model,type,prefix,quantity,numLength) % +type=char(type); +prefix=char(prefix); if type=='rxns' existingIds=model.rxns; diff --git a/core/getAllRxnsFromGenes.m b/core/getAllRxnsFromGenes.m index 4e2bb35d..10471e74 100755 --- a/core/getAllRxnsFromGenes.m +++ b/core/getAllRxnsFromGenes.m @@ -15,10 +15,8 @@ % % Usage: allRxns=getAllRxnsFromGenes(model,rxns) -%If the supplied object is a character array, then convert it to a cell -%array -if ischar(rxns) - rxns={rxns}; +if ~islogical(rxns) && ~isnumeric(rxns) + rxns=convertCharArray(rxns); end rxnIdx=getIndexes(model,rxns,'rxns'); diff --git a/core/getAllowedBounds.m b/core/getAllowedBounds.m index b837d1da..de72a5b0 100755 --- a/core/getAllowedBounds.m +++ b/core/getAllowedBounds.m @@ -18,7 +18,8 @@ if nargin<2 rxns=1:numel(model.rxns); -else +elseif ~islogical(rxns) && ~isnumeric(rxns) + rxns=convertCharArray(rxns); rxns=getIndexes(model,rxns, 'rxns'); end diff --git a/core/getElementalBalance.m b/core/getElementalBalance.m index 84acd258..abe3832e 100755 --- a/core/getElementalBalance.m +++ b/core/getElementalBalance.m @@ -28,6 +28,8 @@ if nargin<2 rxns=[]; +elseif ~islogical(rxns) && ~isnumeric(rxns) + rxns=convertCharArray(rxns); end if nargin<3 diff --git a/core/getEssentialRxns.m b/core/getEssentialRxns.m index 631a7af4..db478802 100755 --- a/core/getEssentialRxns.m +++ b/core/getEssentialRxns.m @@ -16,6 +16,8 @@ if nargin<2 ignoreRxns={}; +else + ignoreRxns=convertCharArray(ignoreRxns); end %Too make sure that it doesn't try to optimize for something @@ -31,7 +33,7 @@ %Check which reactions have flux. Only those can be essential. This is not %the smallest list of reactions, but it's a fast way -rxnsToCheck=setdiff(model.rxns(abs(sol.x)>10^-8),ignoreRxns); +rxnsToCheck=setdiff(model.rxns(abs(sol.x)>10^-12),ignoreRxns); nToCheck=numel(rxnsToCheck); minimize=true; while 1 diff --git a/core/getExchangeRxns.m b/core/getExchangeRxns.m index 86b01492..5d47b60d 100755 --- a/core/getExchangeRxns.m +++ b/core/getExchangeRxns.m @@ -18,6 +18,8 @@ if nargin<2 reactionType='both'; +else + reactionType=char(reactionType); end hasNoProducts=sparse(numel(model.rxns),1); diff --git a/core/getIndexes.m b/core/getIndexes.m index fb2c418f..53b31606 100755 --- a/core/getIndexes.m +++ b/core/getIndexes.m @@ -24,11 +24,10 @@ returnLogical=false; end -%If the supplied object is a character array, then convert it to a cell -%array -if ischar(objects) - objects={objects}; +if ~islogical(objects) && ~isnumeric(objects) + objects=convertCharArray(objects); end +type=char(type); indexes=[]; @@ -79,7 +78,7 @@ elseif ~isempty(index) indexes(i)=index; else - error(['Could not find object ' objects{i} ' in the model']); + error(['Could not find object ''' objects{i} ''' in the model']); end end else diff --git a/core/getMetsInComp.m b/core/getMetsInComp.m index 9ef1ec61..ae6580ed 100755 --- a/core/getMetsInComp.m +++ b/core/getMetsInComp.m @@ -10,9 +10,7 @@ % % Usage: [I, metNames]=getMetsInComp(model,comp) -if ischar(comp) - comp={comp}; -end +comp=char(comp); J=find(ismember(upper(model.comps),upper(comp))); diff --git a/core/getMinNrFluxes.m b/core/getMinNrFluxes.m index bc920701..25deba8c 100755 --- a/core/getMinNrFluxes.m +++ b/core/getMinNrFluxes.m @@ -27,20 +27,14 @@ % % Usage: [x,I,exitFlag]=getMinNrFluxes(model, toMinimize, params, scores) -% glpk solver as implemented by COBRA does not work well for MILP. -global CBT_MILP_SOLVER -if strcmp(getpref('RAVEN','solver'),'cobra') && strcmp(CBT_MILP_SOLVER,'glpk') - dispEM('The current solver is set to ''cobra'', while in COBRA the MILP solver has been set to ''glpk''. The COBRA implementation of glpk is not well suitable for solving MILPs. Please install the Gurobi or an alternative MILP solver.',true); -end - exitFlag=1; if nargin<2 toMinimize=model.rxns; +elseif ~islogical(toMinimize) && ~isnumeric(toMinimize) + toMinimize=convertCharArray(toMinimize); else - if ~iscell(toMinimize) - toMinimize=model.rxns(toMinimize); - end + toMinimize=model.rxns(toMinimize); end %For passing parameters to the solver @@ -127,8 +121,7 @@ prob.lb = [prob.blx; prob.blc]; prob.ub = [prob.bux; prob.buc]; prob.osense=1; -prob.csense=char(zeros(size(prob.a,1),1)); -prob.csense(:)='E'; +prob.csense=repmat('E', 1, size(prob.a,1),1); prob.b=zeros(size(prob.a,1), 1); %Use the output from the linear solution as starting point. Only the values @@ -136,8 +129,9 @@ prob.sol.int.xx=zeros(numel(prob.c),1); prob.sol.int.xx(prob.ints.sub(sol.x(indexes)>10^-12))=1; prob.x0=[]; -prob.vartype=repmat('C', 1, size(prob.A,2)); -prob.vartype(prob.ints.sub) = 'B'; +prob.vartype=repmat('C', size(prob.A,2), 1); +prob.vartype(prob.ints.sub) = 'I'; % with .lb = 0 and .ub = 1, they are binary +% integers (glpk in octave only allows 'continuous' or '', not 'binary') prob=rmfield(prob,{'blx','bux','blc','buc'}); % Optimize the problem diff --git a/core/getModelFromHomology.m b/core/getModelFromHomology.m index 76119130..56c098fe 100755 --- a/core/getModelFromHomology.m +++ b/core/getModelFromHomology.m @@ -69,8 +69,12 @@ hitGenes.oldGenes = []; % collect the old genes from the template model (organism) hitGenes.newGenes = []; % collect the new genes of the draft model (target organism) +getModelFor=char(getModelFor); + if nargin<4 preferredOrder=[]; +else + preferredOrder=convertCharArray(preferredOrder); end if nargin<5 strictness=1; @@ -91,8 +95,6 @@ mapNewGenesToOld=true; end -preferredOrder=preferredOrder(:); - if isfield(models,'S') models={models}; end diff --git a/core/getRxnsInComp.m b/core/getRxnsInComp.m index 2cbd903f..98da6aea 100755 --- a/core/getRxnsInComp.m +++ b/core/getRxnsInComp.m @@ -13,9 +13,7 @@ % % Usage: [I, rxnNames]=getRxnsInComp(model,comp,includePartial) -if ischar(comp) - comp={comp}; -end +comp=char(comp); if nargin<3 includePartial=false; end diff --git a/core/haveFlux.m b/core/haveFlux.m index 57de806d..4b6477e3 100755 --- a/core/haveFlux.m +++ b/core/haveFlux.m @@ -27,6 +27,8 @@ end if nargin<3 rxns=model.rxns; +elseif ~islogical(rxns) && ~isnumeric(rxns) + rxns=convertCharArray(rxns); end %This is since we're maximizing for the sum of fluxes, which isn't possible @@ -49,7 +51,9 @@ %Maximize for all fluxes first in order to get fewer rxns to test smallModel.c=ones(numel(smallModel.c),1); sol=solveLP(smallModel); -J(abs(sol.x(mixIndexes))>cutOff)=true; +if ~isempty(sol.x) + J(abs(sol.x(mixIndexes))>cutOff)=true; +end %Loop through and maximize then minimize each rxn if it doesn't already %have a flux diff --git a/core/makeSomething.m b/core/makeSomething.m index 3c1cb6ab..71ce9eb4 100755 --- a/core/makeSomething.m +++ b/core/makeSomething.m @@ -44,6 +44,8 @@ if nargin<2 ignoreMets=[]; +elseif ~islogical(ignoreMets) && ~isnumeric(ignoreMets) + ignoreMets=convertCharArray(ignoreMets); end if nargin<3 isNames=false; diff --git a/core/mergeModels.m b/core/mergeModels.m index b118e820..fece8065 100755 --- a/core/mergeModels.m +++ b/core/mergeModels.m @@ -1,6 +1,9 @@ function model=mergeModels(models,metParam,supressWarnings) % mergeModels -% Merges models into one model structure. +% Merges models into one model structure. Reactions are added without any +% checks, so duplicate reactions might appear. Metabolites are matched by +% their name and compartment (metaboliteName[comp]), while genes are +% matched by their name. % % models a cell array with model structures % metParam string specifying whether to refer to metabolite name @@ -23,6 +26,8 @@ if nargin<2 metParam='metNames'; +else + metParam=char(metParam); end if nargin<3 @@ -43,16 +48,6 @@ model.geneFrom(:)={models{1}.id}; end -if isfield(model,'subSystems') - hasDeletedSubSystem=false; -else - if supressWarnings==false - EM='Cannot add subsystems since the existing model has no subsystems info. All reactions must have a subsystem for this to be included'; - dispEM(EM,false); - end - hasDeletedSubSystem=true; -end - if isfield(model,'equations') model=rmfield(model,'equations'); end @@ -93,16 +88,21 @@ model.c=[model.c;models{i}.c]; model.rev=[model.rev;models{i}.rev]; - if hasDeletedSubSystem==false - if isfield(models{i},'subSystems') + if isfield(models{i},'subSystems') + if isfield(model,'subSystems') model.subSystems=[model.subSystems;models{i}.subSystems]; else - if supressWarnings==false - EM='Cannot add subsystems since the existing model has no subsystems info. All reactions must have a subsystem for this to be included. Deleting subSystems field'; - dispEM(EM,false); - end - hasDeletedSubSystem=true; - model=rmfield(model,'subSystems'); + emptySubSystem=cell(numel(model.rxns)-numel(models{i}.rxns),1); + emptySubSystem(:)={''}; + emptySubSystem=cellfun(@(x) cell(0,0),emptySubSystem,'UniformOutput',false); + model.subSystems=[emptySubSystem;models{i}.subSystems]; + end + else + if isfield(model,'subSystems') + emptySubSystem=cell(numel(models{i}.rxns),1); + emptySubSystem(:)={''}; + emptySubSystem=cellfun(@(x) cell(0,0),emptySubSystem,'UniformOutput',false); + model.subSystems=[model.subSystems;emptySubSystem]; end end diff --git a/core/permuteModel.m b/core/permuteModel.m index 23a06cea..fc672bec 100755 --- a/core/permuteModel.m +++ b/core/permuteModel.m @@ -15,7 +15,7 @@ % Usage: newModel=permuteModel(model, indexes, type) newModel=model; -indexes=indexes(:); +type=char(type); switch type case 'rxns' diff --git a/core/predictLocalization.m b/core/predictLocalization.m index 283dad71..cd28bd7d 100755 --- a/core/predictLocalization.m +++ b/core/predictLocalization.m @@ -92,6 +92,7 @@ dispEM(EM,false); end +defaultCompartment=char(defaultCompartment); I=ismember(defaultCompartment,GSS.compartments); if I==false EM='defaultCompartment not found in GSS'; diff --git a/core/printFluxes.m b/core/printFluxes.m index cbd2bddc..e0929655 100755 --- a/core/printFluxes.m +++ b/core/printFluxes.m @@ -32,7 +32,7 @@ function printFluxes(model, fluxes, onlyExchange, cutOffFlux, outputFile,outputS % hand sides are lumped % % Usage: printFluxes(model, fluxes, onlyExchange, cutOffFlux, -% outputFile,outputString) +% outputFile,outputString,metaboliteList) if nargin<3 onlyExchange=true; @@ -47,19 +47,21 @@ function printFluxes(model, fluxes, onlyExchange, cutOffFlux, outputFile,outputS fid=1; else if ~isempty(outputFile) + outputFile=char(outputFile); fid=fopen(outputFile,'w'); else fid=1; end end -if nargin<6 - outputString='%rxnID\t(%rxnName):\t%flux\n'; -end -if isempty(outputString) +if nargin<6 || isempty(outputString) outputString='%rxnID\t(%rxnName):\t%flux\n'; +else + outputString=char(outputString); end if nargin<7 metaboliteList={}; +else + metaboliteList=convertCharArray(metaboliteList); end if numel(fluxes)~=numel(model.rxns) EM='The number of fluxes and the number of reactions must be the same'; diff --git a/core/printModel.m b/core/printModel.m index 3080ceff..c05fe4eb 100755 --- a/core/printModel.m +++ b/core/printModel.m @@ -35,23 +35,25 @@ function printModel(model,rxnList,outputString,outputFile,metaboliteList) % % Usage: printModel(model,rxnList,outputString,outputFile,metaboliteList) -if nargin<2 +if nargin<2 || isempty(rxnList) rxnList=model.rxns; +elseif ~islogical(rxnList) && ~isnumeric(rxnList) + rxnList=convertCharArray(rxnList); end -if isempty(rxnList) - rxnList=model.rxns; -end -if nargin<3 - outputString='%rxnID (%rxnName)\n\t%eqn [%lower %upper]\n'; -end -if isempty(outputString) +if nargin<3 || isempty(outputString) outputString='%rxnID (%rxnName)\n\t%eqn [%lower %upper]\n'; +else + outputString=char(outputString); end if nargin<4 outputFile=[]; +else + outputFile=char(outputFile); end if nargin<5 metaboliteList=[]; +else + metaboliteList=convertCharArray(metaboliteList); end I=getIndexes(model,rxnList,'rxns',true)*1.00; %To convert it to "fluxes" diff --git a/core/removeBadRxns.m b/core/removeBadRxns.m index 81ff6f59..0aaceed9 100755 --- a/core/removeBadRxns.m +++ b/core/removeBadRxns.m @@ -74,12 +74,16 @@ end if nargin<3 ignoreMets=[]; +elseif ~islogical(ignoreMets) && ~isnumeric(ignoreMets) + ignoreMets=convertCharArray(ignoreMets); end if nargin<4 isNames=false; end if nargin<5 balanceElements={'C';'P';'S';'N';'O'}; +else + balanceElements=convertCharArray(balanceElements); end if nargin<6 refModel=[]; @@ -167,7 +171,7 @@ %If there are unbalanced rxns then delete one of them and iterate if any(I) - rxnToRemove=I(randsample(numel(I),1)); + rxnToRemove=I(randperm(numel(I),1)); else %If there are no unbalanced rxns in the solution if rxnRules==1 @@ -187,7 +191,7 @@ %If there are any such reactions, remove one of them and %iterate if any(I) - rxnToRemove=I(randsample(numel(I),1)); + rxnToRemove=I(randperm(numel(I),1)); else if rxnRules==2 %This happens when all reactions used are balanced @@ -218,7 +222,7 @@ end end I=find(abs(solution)>10^-8); - rxnToRemove=I(randsample(numel(I),1)); + rxnToRemove=I(randperm(numel(I),1)); end end end diff --git a/core/removeGenes.m b/core/removeGenes.m index d797b6b6..4071e040 100755 --- a/core/removeGenes.m +++ b/core/removeGenes.m @@ -37,8 +37,8 @@ end reducedModel = model; %Only remove genes that are actually in the model -try - ischar(genesToRemove{1}); +if ~islogical(genesToRemove) || ~isnumeric(genesToRemove) + genesToRemove=convertCharArray(genesToRemove); genesToRemove=genesToRemove(ismember(genesToRemove,model.genes)); end if ~isempty(genesToRemove) diff --git a/core/removeMets.m b/core/removeMets.m index 67d92819..1961e493 100755 --- a/core/removeMets.m +++ b/core/removeMets.m @@ -22,9 +22,8 @@ % % Usage: reducedModel=removeMets(model,metsToRemove,isNames,... % removeUnusedRxns,removeUnusedGenes,removeUnusedComps) - -if ischar(metsToRemove) - metsToRemove={metsToRemove}; +if ~islogical(metsToRemove) && ~isnumeric(metsToRemove) + metsToRemove=convertCharArray(metsToRemove); end if nargin<3 @@ -43,16 +42,9 @@ removeUnusedComps=false; end -if isNames==true - %Check that metsToRemove is a cell array - if iscellstr(metsToRemove)==false - if ischar(metsToRemove) - metsToRemove={metsToRemove}; - else - EM='Must supply a cell array of strings if isNames=true'; - dispEM(EM); - end - end +%Check that metsToRemove is a cell array +if isNames==true && ~iscell(metsToRemove) + error('Must supply a cell array of strings if isNames=true'); end reducedModel=model; diff --git a/core/removeReactions.m b/core/removeReactions.m index c0d8615e..59b6c6e5 100755 --- a/core/removeReactions.m +++ b/core/removeReactions.m @@ -27,9 +27,8 @@ if nargin<5 removeUnusedComps=false; end - -if ischar(rxnsToRemove) - rxnsToRemove={rxnsToRemove}; +if ~islogical(rxnsToRemove) && ~isnumeric(rxnsToRemove) + rxnsToRemove=convertCharArray(rxnsToRemove); end reducedModel=model; diff --git a/core/replaceMets.m b/core/replaceMets.m index a50444d2..00f003a3 100755 --- a/core/replaceMets.m +++ b/core/replaceMets.m @@ -17,6 +17,9 @@ % % Usage: model=replaceMets(model,metabolite,replacement,verbose) +metabolite=char(metabolite); +replacement=char(replacement); + if nargin<4 verbose=false; end @@ -26,15 +29,13 @@ % possible. repIdx = find(strcmp(replacement,model.metNames)); if isempty(repIdx) - EM='The replacement metabolite name cannot be found in the model.' - dispEM(EM,true); + error('The replacement metabolite name cannot be found in the model.'); end % Change name and information from metabolite to replacement metabolite metIdx = find(strcmp(metabolite,model.metNames)); if isempty(metIdx) - EM='The to-be-replaced metabolite name cannot be found in the model.' - dispEM(EM,true); + error('The to-be-replaced metabolite name cannot be found in the model.'); end if verbose==true fprintf('\n\nThe following reactions contain the replaced metabolite as reactant:\n') @@ -68,9 +69,11 @@ idxDelete=[]; for i = 1:length(repIdx) metCompsNidx=find(strcmp(metCompsN(repIdx(i)), metCompsN)); - if gt(length(metCompsNidx),1) % If more than 1 metabolite matches - model.S(metCompsNidx(1),:) = model.S(metCompsNidx(1),:) + model.S(metCompsNidx(2:end),:); - idxDelete=[idxDelete; metCompsNidx(2:end)]; % Make list of metabolite IDs to delete + if length(metCompsNidx)>1 + for j = 2:length(metCompsNidx) + model.S(metCompsNidx(1),:) = model.S(metCompsNidx(1),:) + model.S(metCompsNidx(j),:); + idxDelete=[idxDelete; metCompsNidx(j)]; % Make list of metabolite IDs to delete + end end end diff --git a/core/reporterMetabolites.m b/core/reporterMetabolites.m index f43199b4..9e4b08ea 100755 --- a/core/reporterMetabolites.m +++ b/core/reporterMetabolites.m @@ -39,6 +39,7 @@ % Usage: repMets=reporterMetabolites(model,genes,genePValues,printResults,... % outputFile,geneFoldChanges) +genes=convertCharArray(genes); if nargin<4 printResults=false; end diff --git a/core/setExchangeBounds.m b/core/setExchangeBounds.m index 4613f1be..c2631da1 100755 --- a/core/setExchangeBounds.m +++ b/core/setExchangeBounds.m @@ -50,8 +50,8 @@ % handle input arguments if nargin < 2 mets = []; -elseif ischar(mets) - mets = {mets}; % in case only one metabolite is provided as a string +elseif ~islogical(mets) || ~isnumeric(mets) + mets=convertCharArray(mets); end if nargin < 3 || isempty(lb) @@ -186,8 +186,8 @@ [metInd,rxnInd] = find(model_temp.S(exchMetInd,exchRxnInd) ~= 0); % check for any metabolites that are exchanged in more than one reaction -tbl = tabulate(metInd); -repeatedInds = tbl(:,2) > 1; +[tbl,i,~] = unique(metInd,'first'); +repeatedInds = find(not(ismember(1:numel(tbl),i))); multiMetInd = exchMetInd(metInd(repeatedInds)); if ~isempty(multiMetInd) fprintf('WARNING: The following metabolites are involved in more than one exchange reaction:\n'); diff --git a/core/setParam.m b/core/setParam.m index 532c71d2..b7816c23 100755 --- a/core/setParam.m +++ b/core/setParam.m @@ -23,29 +23,22 @@ % % Usage: model=setParam(model, paramType, rxnList, params) +paramType=convertCharArray(paramType); if ~isempty(setdiff(paramType,{'lb';'ub';'eq';'obj';'rev';'var'})) EM=['Incorrect parameter type: "' paramType '"']; dispEM(EM); end -%Allow to set several parameters to the same value -if numel(rxnList)~=numel(params) && numel(params)~=1 - EM='The number of parameter values and the number of reactions must be the same'; - dispEM(EM); -end - if isnumeric(rxnList) || islogical(rxnList) rxnList=model.rxns(rxnList); -elseif ischar(rxnList) - rxnList={rxnList}; +else + rxnList=convertCharArray(rxnList); end -if ischar(paramType) - paramType={paramType}; -end - -if isnumeric(params) - params=[params]; +%Allow to set several parameters to the same value +if numel(rxnList)~=numel(params) && numel(params)~=1 + EM='The number of parameter values and the number of reactions must be the same'; + dispEM(EM); end if length(rxnList)>1 diff --git a/core/simplifyModel.m b/core/simplifyModel.m index be3d97b8..bf970416 100755 --- a/core/simplifyModel.m +++ b/core/simplifyModel.m @@ -59,6 +59,8 @@ end if nargin<9 reservedRxns=[]; +else + reservedRxns=convertCharArray(reservedRxns); end if nargin<10 suppressWarnings=false; diff --git a/core/sortModel.m b/core/sortModel.m index d8d21e5d..bc6b7548 100755 --- a/core/sortModel.m +++ b/core/sortModel.m @@ -73,7 +73,7 @@ end end subsystemsUnique=unique(subsystemsUnique); - for i=1:numel(subsystems) + for i=1:numel(subsystemsUnique) %Get all reactions for that subsystem rxns=find(~cellfun(@isempty,regexp(subsystemsConcatenated,subsystemsUnique(i)))); diff --git a/doc/INIT/getINITModel.html b/doc/INIT/getINITModel.html index f0c1ac98..15549864 100644 --- a/doc/INIT/getINITModel.html +++ b/doc/INIT/getINITModel.html @@ -261,354 +261,358 @@