Test printConstraint.m

src/analysis/FVA/fastFVA/mex/generateMexFastFVA.m

@@ -117,7 +117,7 @@ function generateMexFastFVA(rootPathCPLEX, printLevel)
     eval(['mex -setup c' compVerboseMode]);
 
     % define the name of the source code
-    filename = [CBTDIR filesep 'external' filesep 'analysis' filesep 'fastFVAmex' filesep 'cplexFVA.c'];
+    filename = [CBTDIR filesep 'src' filesep 'analysis' filesep 'FVA' filesep 'fastFVA' filesep 'mex' filesep 'cplexFVA.c'];
 
     % generation of MEX string with compiler options
     if isunix

src/analysis/calculateFluxShifts/calculateFluxShifts.m

@@ -1,6 +1,6 @@
 function calculateFluxShifts(source, target)
 % Compares reaction flux values (min, max) between two models or more models
-% in resultsPostOptimization\contextSpecificModels folder and writes the
+% in resultsPostOptimization/contextSpecificModels folder and writes the
 % comparison results ('up', 'down', 1) in a new column in the target model
 %
 % USAGE:
@@ -16,22 +16,23 @@ function calculateFluxShifts(source, target)
 % .. Authors:
 %       - Kristina Grausa 05/16/2022
 %       - Kristina Grausa 08/23/2022 - standard header and formatting
+%       - Farid Zare      20/11/2023 - Repository addresses are corrected, standard variable names
 
     % Get target model/-s
     destTarget = {};
     if string(target) == "All"
-        destTarget = string(strcat('resultsPostOptimization\contextSpecificModels\*.xls'));
+        destTarget = string(strcat('resultsPostOptimization/contextSpecificModels/*.xls'));
     else
-        destTarget = string(strcat('resultsPostOptimization\contextSpecificModels\*',char(target),'*.xls'));
+        destTarget = string(strcat('resultsPostOptimization/contextSpecificModels/*',char(target),'*.xls'));
     end
-    T = dir(destTarget);
+    xlsFilesInfo = dir(destTarget);
 
     % Source model
-    destSource = string(strcat('resultsPostOptimization\contextSpecificModels\*',char(source),'*.xls'));
+    destSource = string(strcat('resultsPostOptimization/contextSpecificModels/*',char(source),'*.xls'));
     S = dir(destSource);
     sheetName = 'Reaction List';
-    filename = strcat(S.folder, '\', S.name);
-    dataSource=readtable(filename,'Sheet',sheetName);
+    filename = strcat(S.folder, '/', S.name);
+    dataSource = readtable(filename,'Sheet',sheetName);
 
     % Get minFlux and MaxFlux columns from source model
     minFluxSource = dataSource.MinFlux;
@@ -47,59 +48,59 @@ function calculateFluxShifts(source, target)
     end
 
     % Get min and max flux data from target model file
-    for i=1:1:length(T)
-        if ~contains(T(i).name,char(source))
+    for i=1:length(xlsFilesInfo)
 
-            filename = strcat(T(i).folder, '\', T(i).name);
+        if ~contains(xlsFilesInfo(i).name,char(source))
+            filename = strcat(xlsFilesInfo(i).folder, '/', xlsFilesInfo(i).name);
             data=readtable(filename,'Sheet',sheetName);
 
             try
             % Get minFlux and MaxFlux columns from target model file
-            minFlux = data.MinFlux;
-            maxFlux = data.MaxFlux;
+            minFluxTarget = data.MinFlux;
+            maxFluxTarget = data.MaxFlux;
 
             % Create and fill cell arrays (min, max) with rounded flux
             % values from target model
-            minFlux_r = cell(height(minFlux),1);
-            maxFlux_r = cell(height(maxFlux),1);
+            minFluxTarget_r = cell(height(minFluxTarget),1);
+            maxFluxTarget_r = cell(height(maxFluxTarget),1);
 
             for k=1:1:height(data)
-                minFlux_r{k} = round(str2double(minFlux{k}),6);
-                maxFlux_r{k} = round(str2double(maxFlux{k}),6);
+                minFluxTarget_r{k} = round(str2double(minFluxTarget{k}),6);
+                maxFluxTarget_r{k} = round(str2double(maxFluxTarget{k}),6);
             end
 
             % Create and fill flux ratio cell arrays (target/source)
-            MinFluxRatio = cell(height(minFlux),1);
-            MaxFluxRatio = cell(height(minFlux),1);
+            MinFluxRatio = cell(height(minFluxTarget),1);
+            MaxFluxRatio = cell(height(minFluxTarget),1);
 
             % Calculate ratios for all reaction fluxes
-            for n=1:1:length(minFlux)
-                minRatio = minFlux_r{n}/minFluxSource_r{n};
-                maxRatio = maxFlux_r{n}/maxFluxSource_r{n};
+            for n=1:1:length(minFluxTarget)
+                minRatio = minFluxTarget_r{n}/minFluxSource_r{n};
+                maxRatio = maxFluxTarget_r{n}/maxFluxSource_r{n};
                 minRatio = round(minRatio,6);
                 maxRatio = round(maxRatio,6);
 
                 % Source flux is 0
-                if minFluxSource_r{n} == 0 && minFlux_r{n} ~= 0
+                if minFluxSource_r{n} == 0 && minFluxTarget_r{n} ~= 0
                     minRatio = 'up';
                 end
-                if maxFluxSource_r{n} == 0 && maxFlux_r{n} ~= 0
+                if maxFluxSource_r{n} == 0 && maxFluxTarget_r{n} ~= 0
                     maxRatio = 'up';
                 end
 
                 % Target flux is 0
-                if minFlux_r{n} == 0 &&  minFluxSource_r{n} ~= 0
+                if minFluxTarget_r{n} == 0 &&  minFluxSource_r{n} ~= 0
                     minRatio = 'down';
                 end
-                if maxFlux_r{n} == 0 && maxFluxSource_r{n} ~= 0
+                if maxFluxTarget_r{n} == 0 && maxFluxSource_r{n} ~= 0
                     maxRatio = 'down';
                 end
 
                 % Both are 0
-                if minFlux_r{n} == 0 &&  minFluxSource_r{n} == 0
+                if minFluxTarget_r{n} == 0 &&  minFluxSource_r{n} == 0
                     minRatio = 1;
                 end
-                if maxFlux_r{n} == 0 && maxFluxSource_r{n} == 0
+                if maxFluxTarget_r{n} == 0 && maxFluxSource_r{n} == 0
                     maxRatio = 1;
                 end
 
@@ -108,7 +109,7 @@ function calculateFluxShifts(source, target)
             end
 
             % Save results
-            folderName = 'resultsPostOptimization\fluxShifts\';
+            folderName = 'resultsPostOptimization/fluxShifts/';
             if ~exist(folderName, 'dir')
                mkdir(folderName)
             end
@@ -117,7 +118,7 @@ function calculateFluxShifts(source, target)
             ratioData.Properties.VariableNames = {strcat('MinFlux_',source{1}), strcat('MaxFlux_',source{1}), 'MinFluxRatio', 'MaxFluxRatio'};
             newData = [data ratioData];
 
-            temp = split(T(i).name,".");
+            temp = split(xlsFilesInfo(i).name,".");
             resultFileName = strcat(temp(1),'_flux_shifts.xls');
             fullResultFilePath = strcat(folderName, resultFileName);
             writetable(newData,fullResultFilePath{1},'AutoFitWidth',false);

src/analysis/filterNonFluxReactions/filterNonFluxReactions.m

@@ -16,8 +16,9 @@ function filterNonFluxReactions(phenotype)
 % .. Authors:
 %       - Kristina Grausa 05/16/2022
 %       - Kristina Grausa 08/23/2022 - standard header and formatting
+%       - Farid Zare      20/11/2023 - Repository addresses are correcteds
 
-    dest = string(strcat('resultsPostOptimization\contextSpecificModels\*', phenotype, '*.*'));
+    dest = string(strcat('resultsPostOptimization/contextSpecificModels/*', phenotype, '*.*'));
     S = dir(dest);
     count = length(S);
 
@@ -27,7 +28,7 @@ function filterNonFluxReactions(phenotype)
 
     for i=1:1:count
         sheetName = 'Reaction List';
-        filename = strcat(S(i).folder, '\', S(i).name);
+        filename = strcat(S(i).folder, '/', S(i).name);
         data=readtable(filename,'Sheet',sheetName);
 
         % Copy headers
@@ -75,7 +76,7 @@ function filterNonFluxReactions(phenotype)
     end
 
     % Save results
-    folderName = 'resultsPostOptimization\nonFluxReactions\';
+    folderName = 'resultsPostOptimization/nonFluxReactions/';
     if ~exist(folderName, 'dir')
        mkdir(folderName)
     end

src/analysis/filterRateLimittingReactions/filterRateLimittingReactions.m

@@ -7,17 +7,19 @@ function filterRateLimittingReactions(phenotype)
 %
 % USAGE:
 %
-%   filterRateLimittingReactions(phenotype)
+%    filterRateLimittingReactions(phenotype)
 %
 % INPUTS:
-%   phenotype:              char representing the phenotype name provided
-%                           in each model name of the same phanotype
+%    phenotype:              char representing the phenotype name provided
+%                            in each model name of the same phanotype
 %
 % .. Authors:
 %       - Kristina Grausa 05/16/2022
 %       - Kristina Grausa 08/23/2022 - standard header and formatting
+%       - Farid Zare      11/21/2023 - Correction of repository address format 
+%
 
-    dest = string(strcat('resultsPostOptimization\contextSpecificModels\*', phenotype, '*.xls'));
+    dest = string(strcat('resultsPostOptimization/contextSpecificModels/*', phenotype, '*.xls'));
     S = dir(dest);
     count = length(S);
 
@@ -27,7 +29,7 @@ function filterRateLimittingReactions(phenotype)
 
     for i=1:1:count
         sheetName = 'Reaction List';
-        filename = strcat(S(i).folder, '\', S(i).name);
+        filename = strcat(S(i).folder, '/', S(i).name);
         data=readtable(filename,'Sheet',sheetName);
 
         % Copy headers
@@ -79,7 +81,7 @@ function filterRateLimittingReactions(phenotype)
     end
 
     % Save results
-    folderName = 'resultsPostOptimization\rateLimittingReactions\';
+    folderName = 'resultsPostOptimization/rateLimittingReactions/';
     if ~exist(folderName, 'dir')
        mkdir(folderName)
     end

src/analysis/multiSpecies/SteadyCom/subroutines/SteadyComSubroutines.m

@@ -234,7 +234,7 @@
         spName = spAbbr;
     end
     if isfield(infoCom, 'spBm')
-        indCom.spBm = modelCom.rxns(infoCom.spBm);
+        indCom.spBm = modelCom.rxns(infoCom.spBm); %modelCom.infoCom.spBm
     end
     if isfield(infoCom,'spATPM')
         indCom.spATPM = modelCom.rxns(infoCom.spATPM);

src/analysis/multiSpecies/microbiomeModelingToolbox/pairwiseInteractionModeling/createMultipleSpeciesModel.m

@@ -1,4 +1,4 @@
-function [modelJoint] = createMultipleSpeciesModel(models, varargin)
+function [modelJoint] = createMultipleSpeciesModel(models, biomasses, varargin)
 % Based on the implementation from *Klitgord and Segre 2010, PMID 21124952*.
 % The present implementation has been used in *PMID 23022739*, *PMID 25841013*,
 % *PMID 25901891*, *PMID 27893703*.
@@ -24,6 +24,9 @@
 %
 %                         * models{1,1} = model 1
 %                         * models{2,1} = model 2...
+%    biomasses:           Cell array containing names of biomass objective
+%                         functions of models to join. Needs to be the same
+%                         length as models.
 %
 % OPTIONAL INPUTS:
 %    nameTagsModels:    cell array of tags for reaction/metabolite abbreviation
@@ -144,7 +147,6 @@
 modelHost.mets(metIndices) = strrep(modelHost.mets(metIndices), '[e0]', '[e]');
 end
 %% Remove futile cycles that may appear after joining (optional)
-
 if remCyclesFlag
     database = loadVMHDatabase;
 
@@ -156,8 +158,7 @@
     unionRxns=unique(unionRxns);
     for i = 1:modelNumber
         model=models{i, 1};
-        biomassReaction=model.rxns(find(strncmp(model.rxns, 'bio', 3)));
-        model = removeFutileCycles(model, biomassReaction, database,unionRxns);
+        model = removeFutileCycles(model, biomasses{i}, database,unionRxns);
         models{i, 1} = model;
     end
 end

src/analysis/multiSpecies/microbiomeModelingToolbox/pairwiseInteractionModeling/joinModelsPairwiseFromList.m

@@ -185,12 +185,16 @@ function joinModelsPairwiseFromList(modelList, modelFolder, varargin)
             model1
             model2
             };
+        modelBiomasses = {
+            biomasses{i}
+            biomasses{k}
+            };
         nameTagsModels = {
             strcat(modelList{i}, '_')
             strcat(modelList{k}, '_')
             };
         if isempty(find(contains(existingModels,['pairedModel', '_', modelList{i}, '_', modelList{k}, '.mat'])))
-            [pairedModel] = createMultipleSpeciesModel(models, 'nameTagsModels', nameTagsModels,'mergeGenesFlag', mergeGenesFlag, 'remCyclesFlag', remCyclesFlag);
+            [pairedModel] = createMultipleSpeciesModel(models, modelBiomasses, 'nameTagsModels', nameTagsModels,'mergeGenesFlag', mergeGenesFlag, 'remCyclesFlag', remCyclesFlag);
             [pairedModel] = coupleRxnList2Rxn(pairedModel, pairedModel.rxns(strmatch(nameTagsModels{1, 1}, pairedModel.rxns)), strcat(nameTagsModels{1, 1}, biomasses{i}), c, u);
             [pairedModel] = coupleRxnList2Rxn(pairedModel, pairedModel.rxns(strmatch(nameTagsModels{2, 1}, pairedModel.rxns)), strcat(nameTagsModels{2, 1}, biomasses{k}), c, u);
             pairedModelsTemp{k} = pairedModel;