Merge pull request #1378 from opencobra/develop

Regular merge of develop

src/base/install/prepareTest.m

@@ -9,7 +9,8 @@
 % OPTIONAL INPUTS:
 %    varagin:      `ParameterName` value pairs with the following options:
 %
-%                   - `toolboxes` or `requiredToolboxes`: Names of required toolboxes (the license feature name) (default: {})
+%                   - `toolboxes` or `requiredToolboxes`: Names of required toolboxes (the license feature name).(default: {})
+%                   - `minimalMatlabSolverVersion`: Minimal version of the optimization toolbox required to use the matlab solver.(default: 0)
 %                   - `requiredSolvers`: Names of all solvers that MUST be available. If not empty, the resulting solvers struct will contain cell arrays (default: {})
 %                   - `useSolversIfAvailable`: Names of solvers that should be used if available. If not empty, the resulting solvers struct will contain cell arrays (will not throw an error if not). (default: {})
 %                   - `requireOneSolverOf`: Names of solvers, at least one of which has to be available
@@ -86,6 +87,8 @@
 parser = inputParser();
 parser.addParamValue('toolboxes', {}, @iscell);
 parser.addParamValue('requiredToolboxes', {}, @iscell);
+parser.addParamValue('minimalMatlabSolverVersion',0,@isnumeric);
+
 parser.addParamValue('requiredSolvers', {}, @iscell);
 parser.addParamValue('useSolversIfAvailable', {}, @iscell);
 parser.addParamValue('requireOneSolverOf', {}, @iscell);
@@ -103,6 +106,7 @@
 parser.addParamValue('needsWebAddress', '', @ischar);
 parser.addParamValue('needsWebRead', false, @(x) islogical(x) || x == 1 || x == 0);
 
+
 parser.parse(varargin{:});
 
 useQP = parser.Results.needsQP;
@@ -130,6 +134,8 @@
 useMinimalNumberOfSolvers = parser.Results.useMinimalNumberOfSolvers;
 runtype = getenv('CI_RUNTYPE');
 
+minimalMatlabSolverVersion = parser.Results.minimalMatlabSolverVersion;
+
 errorMessage = {};
 infoMessage = {};
 
@@ -191,6 +197,17 @@
 % restrict the solvers available for this test
 solversForTest = availableSolvers;
 
+if any(ismember(solversForTest.LP,'matlab'))
+    boxes = ver();
+    optBox = find(ismember({boxes.Name},'Optimization Toolbox'));
+    if ~isempty(optBox)
+        optVer = boxes(optBox).Version;
+        if str2double(optVer) < minimalMatlabSolverVersion                
+            excludedSolvers = [columnVector(excludedSolvers);'matlab'];
+        end
+    end
+end
+
 if ~isempty(excludedSolvers)
     solverTypes = fieldnames(availableSolvers);
     for i = 1:numel(solverTypes)

src/base/io/utilities/writeSBML.m

@@ -552,54 +552,58 @@
 
 %% Set the objective sense of the FBC objective according to the osenseStr in
 %the model.
-objectiveSense = 'maximize';
-
-if isfield(model,'osenseStr') && strcmpi(model.osenseStr,'min')
-    objectiveSense = 'minimize';
-end
-
-tmp_fbc_objective=getSBMLDefaultStruct('SBML_FBC_OBJECTIVE',sbmlLevel, sbmlVersion,sbmlPackages, sbmlPackageVersions);
-tmp_fbc_objective.fbc_id = 'obj';
-tmp_fbc_objective.fbc_type = objectiveSense;
-sbmlModel.fbc_activeObjective = 'obj';
-
-tmp_fbc_objective.fbc_fluxObjective=getSBMLDefaultStruct('SBML_FBC_FLUXOBJECTIVE',sbmlLevel, sbmlVersion,sbmlPackages, sbmlPackageVersions);
-
 %%%%% multiple objectives
 if ~isnumeric(model.c)
     model.c=double(cell2mat(model.c)); % convert the variable type to double
 end
-% check=~isempty(model.c(model.c~=0)) the following block is necessry for
-% libSBML library 5.11
-% % % % % % % if isempty(model.c(model.c~=0)) % if no objective function is defined, the first reaction is set as the objective function.
-% % % % % % %     model.c(1)=1;
-% % % % % % % end
-%if ~isempty(model.c(model.c~=0))
-% Construct a default structure of objective reactions and set intial values.
-%     fbc_objective.fbc_fluxObjective=fbc_fluxObjective;
-ind=find(model.c); % Find the index numbers for the objective reactions
-% The fields of a COBRA model are converted into respective fields of a FBCv2 structure.
-if isempty(ind)
-    tmp_fbc_objective.fbc_fluxObjective.fbc_coefficient=0; % no objective function is set
-    sbmlModel.fbc_objective=tmp_fbc_objective;
-else
-    for i=1:length(ind)
-        %     model.c(ind(i));
-        values=model.c(model.c~=0);
-        tmp_fbc_objective.fbc_fluxObjective.fbc_reaction=sbmlModel.reaction(ind(i)).id; % the reaction.id contains the  % model.rxns{ind(i)};
-        tmp_fbc_objective.fbc_fluxObjective.fbc_coefficient=values(i);
-        tmp_fbc_objective.fbc_fluxObjective.isSetfbc_coefficient=1;
-        if i==1
-            sbmlModel.fbc_objective=tmp_fbc_objective;
-        else
-            sbmlModel.fbc_objective.fbc_fluxObjective=[sbmlModel.fbc_objective.fbc_fluxObjective,tmp_fbc_objective.fbc_fluxObjective];
+if ~all(model.c == 0)
+
+
+    objectiveSense = 'maximize';
+
+    if isfield(model,'osenseStr') && strcmpi(model.osenseStr,'min')
+        objectiveSense = 'minimize';
+    end
+
+    tmp_fbc_objective=getSBMLDefaultStruct('SBML_FBC_OBJECTIVE',sbmlLevel, sbmlVersion,sbmlPackages, sbmlPackageVersions);
+    tmp_fbc_objective.fbc_id = 'obj';
+    tmp_fbc_objective.fbc_type = objectiveSense;
+    sbmlModel.fbc_activeObjective = 'obj';
+
+    tmp_fbc_objective.fbc_fluxObjective=getSBMLDefaultStruct('SBML_FBC_FLUXOBJECTIVE',sbmlLevel, sbmlVersion,sbmlPackages, sbmlPackageVersions);
+
+
+    % check=~isempty(model.c(model.c~=0)) the following block is necessry for
+    % libSBML library 5.11
+    % % % % % % % if isempty(model.c(model.c~=0)) % if no objective function is defined, the first reaction is set as the objective function.
+    % % % % % % %     model.c(1)=1;
+    % % % % % % % end
+    %if ~isempty(model.c(model.c~=0))
+    % Construct a default structure of objective reactions and set intial values.
+    %     fbc_objective.fbc_fluxObjective=fbc_fluxObjective;
+    ind=find(model.c); % Find the index numbers for the objective reactions
+    % The fields of a COBRA model are converted into respective fields of a FBCv2 structure.
+    if isempty(ind)
+        tmp_fbc_objective.fbc_fluxObjective.fbc_coefficient=0; % no objective function is set
+        sbmlModel.fbc_objective=tmp_fbc_objective;
+    else
+        for i=1:length(ind)
+            %     model.c(ind(i));
+            values=model.c(model.c~=0);
+            tmp_fbc_objective.fbc_fluxObjective.fbc_reaction=sbmlModel.reaction(ind(i)).id; % the reaction.id contains the  % model.rxns{ind(i)};
+            tmp_fbc_objective.fbc_fluxObjective.fbc_coefficient=values(i);
+            tmp_fbc_objective.fbc_fluxObjective.isSetfbc_coefficient=1;
+            if i==1
+                sbmlModel.fbc_objective=tmp_fbc_objective;
+            else
+                sbmlModel.fbc_objective.fbc_fluxObjective=[sbmlModel.fbc_objective.fbc_fluxObjective,tmp_fbc_objective.fbc_fluxObjective];
+            end
         end
     end
 end
 
 
 
-
 %end
 
 sbmlModel.namespaces=struct('prefix',{'','fbc','groups'},...

src/dataIntegration/metabotools/calculateLODs.m

@@ -14,15 +14,22 @@
 %   lod_mM:            Detection limits in mM
 %
 % .. Author: - Maike K. Aurich 27/05/15
+%            - Modified by Loic Marx, November 2018
 
-lod_gL = lod_ngmL* 0.000001;
+% input checking
+if size(lod_ngmL) == 1
+    lod_ngmL = repmat(lod_ngmL, 1, length(theo_mass)); 
+    end    
 
-lodmM =[];
-
-for i=1:length(lod_gL)
+if size(lod_ngmL) ~= size(theo_mass)
+    error('The number of elements in the input vectors do not match. They have to be either the same size, or lod_ngmL has to be a single value which is used for all elements');
+end
 
-    lod_mM(i,1)= lod_gL(i,1)/theo_mass(i,1)*1000;
+lod_gL = lod_ngmL * 1e-6; 
+lod_mM = zeros(length(lod_gL), 1);
 
+for i = 1 : length(lod_gL)
+    lod_mM(i) = lod_gL(i)/theo_mass(i)*1000;
 end
 
 end

src/dataIntegration/transcriptomics/eFlux/verifyEFluxExpressionStruct.m

@@ -17,7 +17,7 @@
     % checkoptional field preprocessed
     % This leads to -1 for unassociated genes
     if ~isempty(setdiff(model.genes,expression.target))
-        error('The following genes are lacking assignments:\n%s\nAll genes need to be assigned in order to use eFlux. You can remove genes by using the removeGenes function',strjoin(setdiff(model.genes,expression.target),'\n')); 
+        error('The following genes are lacking assignments:\n%s\nAll genes need to be assigned in order to use eFlux. You can remove genes by using the removeGenesFromModel function',strjoin(setdiff(model.genes,expression.target),'\n')); 
     end                
 end
 % check sizes of the value/target fields. They must be of equal size.

src/reconstruction/refinement/removeGenes.m → src/reconstruction/refinement/removeGenesFromModel.m

@@ -1,12 +1,12 @@
-function [model, affectedRxns, originalGPRs, deletedReaction] = removeGenes(model, geneList, varargin)
+function [model, affectedRxns, originalGPRs, deletedReaction] = removeGenesFromModel(model, geneList, varargin)
 % Removes the given genes from the model. GPR rules will be adjusted to reflect the removal.
 % By default, the rules are converted to DNF and all clauses containing any of the given
 % genes are removed. Note, that this function is not supposed to be used to model
 % single gene deletions.
 %
 % USAGE:
 %
-%    [model, affectedRxns, originalGPRs, deletedReactions] = deleteModelGenes(model, geneList)
+%    [model, affectedRxns, originalGPRs, deletedReactions] = removeGenesFromModel(model, geneList)
 %
 % INPUT:
 %    model:               COBRA model with the appropriate constrains for a

test/testAll.m

@@ -25,7 +25,10 @@
     else
         launchTestSuite = false;
     end
-else
+else            
+    % on the CI, always reset the path to make absolutely sure, that we test
+    % the current version
+    restoredefaultpath()
     launchTestSuite = true;
 end
 
@@ -166,7 +169,7 @@
     if launchTestSuite
         % save the userpath
         originalUserPath = path;
-
+        
         % run the tests in the subfolder verifiedTests/ recursively
         [result, resultTable] = runTestSuite();
 

test/verifiedTests/analysis/testFBA/testMinimizeModelFlux.m

@@ -29,21 +29,21 @@
     % qpng does not support this model size, so we don't use quadratic
     % minimzation.
     sol = minimizeModelFlux(currentmodel,'min','one');
-    assert(sol.x(end) == 0);
+    assert(abs(sol.x(end)) < tol);
     sol = minimizeModelFlux(currentmodel);
-    assert(sol.x(end) == 12000); %Since its reversible, exchangers cycle and all rev reactions cycle.
+    assert(abs(sol.x(end) - 12000) <= tol); %Since its reversible, exchangers cycle and all rev reactions cycle.
     sol = minimizeModelFlux(currentmodel,'max','one'); %same as before.
-    assert(sol.x(end) == 12000); %Since its reversible, exchangers cycle and all rev reactions cycle.
+    assert(abs(sol.x(end) - 12000 ) <= tol); %Since its reversible, exchangers cycle and all rev reactions cycle.
     currentmodel.osenseStr = 'min';
     modelChanged = changeRxnBounds(currentmodel,'R3',5,'l');
     sol = minimizeModelFlux(modelChanged);
-    assert(sol.x(end) == 20); %Can only come from the cycle
+    assert(abs(sol.x(end) - 20)  <= tol); %Can only come from the cycle
     sol = minimizeModelFlux(modelChanged,'max');
-    assert(sol.x(end) == 11000); %MAx flux through C-> E and Exchangers, + 3 reactions from the cycle.
+    assert(abs(sol.x(end) - 11000)  <= tol); %MAx flux through C-> E and Exchangers, + 3 reactions from the cycle.
     modelChanged = changeRxnBounds(currentmodel,'EX_E',5,'l'); %Force production of E
     sol = minimizeModelFlux(modelChanged);
     modelChanged = rmfield(modelChanged,'osenseStr');
-    assert(sol.x(end) == 25); %Flux -> A -> B -> C -> E -> 
+    assert(abs(sol.x(end) - 25)  <= tol); %Flux -> A -> B -> C -> E -> 
 end
 
 % change the directory

test/verifiedTests/analysis/testFVA/testFVA.m

@@ -22,7 +22,10 @@
 
 % define the solver packages to be used to run this test, can't use
 % dqq/Minos for the parallel part.
-solverPkgs = prepareTest('needsLP',true,'needsMILP',true,'needsQP',true,'useSolversIfAvailable',{'ibm_cplex'},'excludeSolvers',{'dqqMinos','quadMinos'});
+solverPkgs = prepareTest('needsLP',true,'needsMILP',true,'needsQP',true,...
+                         'useSolversIfAvailable',{'ibm_cplex'},...
+                         'excludeSolvers',{'dqqMinos','quadMinos'},...
+                         'minimalMatlabSolverVersion',8.0);
 
 
 
@@ -143,7 +146,7 @@
 %there
 pttoolboxPath = which('parpool');
 % here, we can use dqq and quadMinos again, because this is not parallel.
-solverPkgs = prepareTest('needsLP',true,'needsMILP',true,'needsQP',true,'useSolversIfAvailable',{'ibm_cplex'});
+solverPkgs = prepareTest('needsLP',true,'needsMILP',true,'needsQP',true,'useSolversIfAvailable',{'ibm_cplex'},'minimalMatlabSolverVersion',8.0);
 if ~isempty(pttoolboxPath)
     %We also have to shut down the parallel pool, as otherwise problems
     %occur with the pool.

test/verifiedTests/analysis/testGeometricFBA/testGeometricFBA.m

@@ -26,7 +26,7 @@
 
 %When detectDeadEnds is changed according to Ronans suggestion, we need to test
 %multiple solvers.
-solverPkgs = prepareTest('needsLP', true);
+solverPkgs = prepareTest('needsLP', true, 'minimalMatlabSolverVersion',8.0);
 tol = 1e-4;
 for k = 1:length(solverPkgs.LP)
     % set the solver

test/verifiedTests/base/testSolvers/testSolveCobraQP.m

@@ -23,7 +23,9 @@
 
 % test solver packages
 useIfAvailable = {'tomlab_cplex','ibm_cplex', 'gurobi','qpng','ibm_cplex','mosek'};
-solverPkgs = prepareTest('needsQP',true,'useSolversIfAvailable', useIfAvailable);
+% pdco is a normalizing solver not a general purpose QP solver, so it will
+% fail the test
+solverPkgs = prepareTest('needsQP',true,'useSolversIfAvailable', useIfAvailable,'excludeSolvers',{'pdco'});
 
 %QP Solver test: http://tomopt.com/docs/quickguide/quickguide005.php
 

test/verifiedTests/dataIntegration/testEFlux/testEFlux.m

@@ -16,7 +16,7 @@
 solverPkgs = prepareTest('needsLP',true,'toolboxes',{'bioinformatics_toolbox'},'needsWebAddress','https://www.ncbi.nlm.nih.gov/geo/query/','excludeSolvers',{'matlab'});
 
 model = getDistributedModel('ecoli_core_model.mat');
-model = removeGenes(model,'s0001','keepReactions',true);
+model = removeGenesFromModel(model,'s0001','keepReactions',true);
 geoaddress = 'https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?form=text&acc=%s&view=full';
 anaerobic = geosoftread(urlread(sprintf(geoaddress,'GSM1010240')));
 aerobic = geosoftread(urlread(sprintf(geoaddress,'GSM1126445')));

test/verifiedTests/dataIntegration/testMetaboTools/testCalculateLODs.m

@@ -0,0 +1,42 @@
+% The COBRA Toolbox: testcalculateLODs
+%
+% Purpose:
+%     - tests the calculateLODs function
+%
+% Author: - Loic Marx, November 2018
+%    
+
+% save the current path
+currentDir = pwd;
+
+% initialize the test
+fileDir = fileparts(which('testcalculateLODs'));
+cd(fileDir);
+
+% define the inputs
+lod_ngmL = 1;
+theo_mass = [1, 2];
+
+% calculate the reference value
+v_ref = 1e-3 * [1; 0.5];
+
+% calculate the actual value by using the function
+calculateLOD = calculateLODs(theo_mass,lod_ngmL);
+
+% compare the calculated value and the reference value
+assert(isequal(calculateLOD, v_ref))
+
+% throw an error if the dimensions of both inputs are not equal
+lod_ngmL_unequal = [1, 2, 3];
+assert(verifyCobraFunctionError('calculateLODs', 'inputs', {theo_mass,lod_ngmL_unequal'}, 'testMessage', 'The number of elements in the input vectors do not match. They have to be either the same size, or lod_ngmL has to be a single value which is used for all elements'));
+
+% lod_ngmL and theo_mass are two vectors of the same size
+lod_ngmL_vector = [1, 2];
+v_ref_vector = 1e-3 * [1; 1];
+
+% compare the values calculated when the inputs have the same size and the reference value
+calculateLOD_vector = calculateLODs(theo_mass,lod_ngmL_vector);
+assert(isequal(calculateLOD_vector, v_ref_vector));
+
+% change the directory
+cd(currentDir)

test/verifiedTests/reconstruction/testModelManipulation/testModelManipulation.m

@@ -455,14 +455,14 @@
 fprintf('>> Done \n\n >> Testing Gene removal...\n');
 
 %Test removal of a gene
-modelMod1 = removeGenes(modelMod,'G1');
+modelMod1 = removeGenesFromModel(modelMod,'G1');
 % now, rules{1} and rules{3} should be equal;
 fp = FormulaParser();
 rule = fp.parseFormula(modelMod1.rules{1});
 rule2 = fp.parseFormula(modelMod1.rules{3});
 assert(rule2.isequal(rule));
 % and now without keeping the clauses
-modelMod2 = removeGenes(modelMod,'G1','keepClauses',false);
+modelMod2 = removeGenesFromModel(modelMod,'G1','keepClauses',false);
 fp = FormulaParser();
 rule = fp.parseFormula(modelMod2.rules{1});
 rule2 = fp.parseFormula(modelMod2.rules{2});