Pairwise Interaction Modeling

calculateInteractionsByTaxon(pairwiseInteractions, taxonInformation)

Part of the Microbiome Modeling Toolbox. This script analyzes computed pairwise interactions on different taxonomical levels to yield the interactions that each taxon (genus, family, order, class, phylum) displayed in total. The computed pairwise interactions and the taxnomical information for each analyzed strain are required as inputs.

USAGE:

[interactionsByTaxon]=calculateinteractionsByTaxon(pairwiseInteractions,taxonInfo)

INPUTS:

pairwiseInteractions: a table with pairwise interactions computed

for a number of analyzed microbe models (output of the function simulatePairwiseInteractions)

taxonInformation: a table with taxonomical information on the

analyzed microbes. Needs to contain at least six columns: 1. The names of the analyzed model structures in the first column 2-6. A column with the appropriate information and the header ‘Genus’,’Family’,’Order’,’Class’,’Phylum’ respectively.

OUTPUT:

interactionsByTaxon: a structure with the outcomes predicted for

all taxa on each taxonomical level.

computeParetoOptimality(model, rxn1, rxn2, varargin)

Performs Pareto optimality analysis for two objective functions by simultaneously optimizing two reactions (e.g., the biomass objective functions of two joined organisms). The result is a depiction of the tradeoff between the two competing objectives.

Multiobjective analysis of two reactions is performed by using the approach described in Oberhardt, M. A., J. B. Goldberg, et al. (2010). “Metabolic network analysis of Pseudomonas aeruginosa during chronic cystic fibrosis lung infection.” J Bacteriol 192(20): 5534-5548. One reaction is fixed at different intervals and the other reaction is optimized.

USAGE:

[ParetoFrontier,fluxSolutions,minFluxes,maxFluxes] = computeParetoOptimality(model,rxn1,rxn2,dinc,FVAflag)

INPUTS:

model: COBRA metabolic reconstruction rxn1: Reaction ID of the first reaction to be optimized rxn2: Reaction ID of the second reaction to be optimized

OPTIONAL INPUTS:

dinc: An index which indicates the distance between steps at

which the flux of each reaction is fixed (default=0.001).

FVAflag: If true, flux variability analysis is performed at each

step.

OUTPUTS:

ParetoFrontier: Lists the objective values for both reactions next to

the interval step. Column 1:the interval step. Column 2 and 3: the flux values of rxn1 in and rxn2, respectively, at the interval step.

OPTIONAL OUTPUTS:

fluxSolutions: Contains the flux solutions every interval step as a

matrix of structures

minFluxes: Reports fastFVA results for every step (minFlux) with

one column for each interval step

maxFluxes: Reports fastFVA results for every step (maxFlux) with

one column for each interval step

computeRescuedGenes(varargin)

Part of the Microbiome Modeling Toolbox. This function determines the effect of the presence of another species on gene deletions in a species. A joint model consisting of two species is entered. For each gene deletion that has an effect in each single model, the function conputes whether the presence of the other species can rescue the flux through the objective function (e.g., biomass).

USAGE:

[OptSolKO,OptSolWT,OptSolRatio,fluxesKO]=computeRescuedGenes(varargin)

INPUTS:

modelJoint: Joint model structure consisting of two COBRA models Rxn1: Reaction in the first joined COBRA model for which the

effect of gene deletions is calculated

Rxn2: Reaction in the second joined COBRA model for the effect

of gene deletions is calculated

nameTag1: Species identifier for reactions of the first model nameTag2: Species identifier for reactions of the second model OriModel1: Original COBRA model structure of first model

(needed to determine reactions to be constrained)

OriModel2: Original COBRA model structure of second model

(needed to determine reactions to be constrained)

OUTPUTS:

OptSolKO: Matlab structure containing the computed optimal solutions

for the gene deletions that had an effect on the two reactions

OptSolRatio: Matlab structure containing all knockout to wildtype

optimal solution ratios for the gene deletions that had an effect on the two reactions

OptSolWT: Matlab structure containing wildtype flux values for both

reactions in the two models

RescuedGenes: List of genes for which deletion is lethal or causes >50%

growth reduction in single but not pairwise model

fluxesKO: Matlab structure containing all solutions for all gene

deletions and for both reactions in the two models. Can be used to identify mechanisms for rescued genes.

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.

Joins one or more COBRA models with or without another COBRA model representing the host. The created setup when a host is entered is depicted schematically in Figures 1 and 2 in PMID 27893703.

Creates a common space u (lumen) through which all cells can feed and exchange metabolites, and separate extracellular spaces for all joined models. If a host model is entered, a separate compartment b (body fluids) which has no connection to extracellular space is created for the host. Metabolites can be transported from the lumen to e, from e to the cytosol and from the cytosol to b, but not from body fluids to the lumen.

USAGE:

[modelJoint] = createMultipleSpeciesModel(models, varargin)

INPUTS:

models: cell array of COBRA models(at least one).

Format

• 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

corresponding to each model. Format

• nameTagsModels{1,1} = ‘name tag 1’

• nameTagsModels{2,1} = ‘name tag 2’…

modelHost: COBRA model for host nameTagHost: string of tag for reaction/metabolite abbreviation of host model mergeGenesFlag: If true, the gene associations in both models are

included in the joined model. If false, empty fields are created instead (default:false). Note: merging genes is time-consuming and may crash certain models.

remCyclesFlag: If true, the function will attempt to remove futile

cycles that appear after merging the models (default: false).

OUTPUT:

modelJoint: model structure for joint model

Note

This function assumes, that exchange reactions are identified by containing ‘EX’ in the reaction name and that no other reactions do have this property!

getMultiSpeciesModelId(modelJoint, nameTagsModels, nameTagHost, metTagRe, rxnTagRe, compCom, compHost)

Get names and IDs for metabolites and exchange reactions in the [u] and [b] space

USAGE:

[names, ids] = getMultiSpeciesModelId(modelJoint, nameTagsModels, nameTagHost, metTagRe, rxnTagRe, compCom, compHost)

INPUTS:

modelJoint: COBRA multi-organism model nameTagsModels: cell array of tags for species to identify the respective

reactions and metabolites from modelJoint.rxns and .mets

OPTIONAL INPUTS:

nameTagHost: string of tag for the host model if exist. Input [] if no host is present metTagRe: a regular expression to identify the tag in nameTagsModels from modelJoint.mets.

Use ‘%s’ for the tag in nameTagsModels for each species. Default ‘^%s’, the beginning of the string. In this case, if the tag for species 1 is ‘SP1’, then mets with id ‘SP1glc-D[e]’ will be identified as belonging to species 1. E.g., ‘met[compartment_SP1]’ where ‘SP1’ is the tag for species 1, input ‘[[^_]+_%s]$’.

rxnTagRe: a regular expression to idenify the tag in nameTagsModels from modelJoint.rxns.

Default ‘^%s’, the beginning of the string

compCom: compartment Id for the community exchange compartment. Default ‘u’ compHost: compartment Id for the exchange compartment specific to the host. Default ‘b’

OUTPUTS:

names: structure of reaction/metabolite name (.rxns/.mets) with the following fields:

• spAbbr - species abbreviation (= nameTags). Host put at the end

• EXcom - #met[u]-by-1 cell. Exchange reactions for community metabolites met[u]

• EXhost - #met[b]-by-1 cell. Exchange reactions for host-specific exchange metabolites met[b]

• EXsp - #met[u]-by-#species cell. EXsp(i,j) is the species-community exchange reactions for the i-th met[u] and the j-th species

• Mcom - #met[u]-by-1 cell. All community metabolites met[u]

• Mhost - #met[u]-by-1 cell. Host-specific exchange metabolites met[b]

• Msp - #met[u]-by-#species cell. Msp(i,j) is the metabolite in met[e] of the j-th species participating in reaction EXsp(i,j)

• rxnSps - #rxns-by-1 cell. rxnSps(i) would be spAbbr(j) if the i-th reaction belongs to the j-th species

• metSps - #mets-by-1 cell. metSps(i) would be spAbbr(j) if the i-th metabolite belongs to the j-th species

ids: structure of reaction/metabolite index with the same fields as names except without spAbbr

joinModelsPairwiseFromList(modelList, modelFolder, varargin)

This function joins a list of microbial genome-scale reconstructions in all combinations. Models are not paired with themselves and pairs are only created once (model1+model2 = model2+model1). The reactions in each joined reconstruction structure receive a suffix (entered as input parameter model list). Coupling constraints are created in which the flux through all reactions of each individual microbe is coupled to the flux through its own biomass objective function. By default, coupling constraints are implemented with a coupling factor c of 400 and a threshold u of 0. The coupling parameters can be changed by entering the input parameters c and u. By default, the gene-protein-reaction associations of the joined models are not merged as the gene associations are not needed for the simulation of pairwise interactions and significantly slow down pairwise model creation. If merging of genes is desired, please set the input parameter mergeGenes to true. If you are not using AGORA/AGORA2 reconstructions, you may need to define the biomass objective functions manually.

Please note: the function takes a cell array with the names of the COBRA models to join as well as the name of the folder where they are located at as the inputs, e.g.: joinModelsPairwiseFromList(modelList,’/Users/almut.heinken/Documents/AGORA’);

USAGE:

joinModelsPairwiseFromList(modelList, modelFolder, varargin)

INPUTS:

modelList: Cell array with names of reconstruction structures to be

joined

modelFolder: Path to folder with COBRA model structures to be joined

OPTIONAL INPUTS:

biomasses: Cell array containing names of biomass objective

functions of models to join. Needs to be the same length as modelList.

c: Coupling factor by which reactions in each joined model are

coupled to its biomass reaction (default: 400)

u: Threshold representing minimal flux allowed when flux

through biomass reaction in zero (default: 0)

mergeGenes: If true, genes in the joined models are merged and included

in the joined model structure (default: false)

remCyclesFlag: If true, futile cycles that may arise from joining the

models are removed (default: false). Recommended if AGORA/AGORA2 reconstructions are used, does not work with other namespaces.

numWorkers: Number of workers in parallel pool if desired pairwiseModelFolder: Folder where pairwise models will be saved

printUptakeBoundCom(model, SpFlag, metFlag)

Prints the uptake bounds of the whole community and individual species for a community COBRA model

USAGE:

rxnIDs = printUptakeBoundCom(model, SpFlag, metFlag)

INPUT:

model: the community model with field .infoCom or .indCom indicating the indicies of

community exchange reactions/metabolites. Can be obtained from getMultiSpeciesModelId.m

OPTIONAL INPUTS:

SpFlag: true to show individual uptake rates though community uptake is not allowed (default false) metFlag: true to print with model.metNames (default false)

OUTPUTS:

printMatrix: matrix of the uptake bounds being printed printMet: column of metabolites whose uptake bounds are printed

simulatePairwiseInteractions(pairwiseModelFolder, varargin)

This function predicts the outcome of pairwise simulations in every combination from a given list of pairwise models. The pairwise models need to be created first with the function joinModelsPairwiseFromList. This script requires the COBRA Toolbox function solveCobraLP. Due to the coupling constraints on the pairwise models, the simulations cannot currently be run with optimizeCbModel.

Below is a description of all possible consequences for the two joined organisms. Please note that the outcomes depend on the two genome-scale reconstructions joined and are highly dependent on the applied constraints.

• Competition: both organisms grow slower in co-growth than separately (same outcome for both).

• Parasitism: one organism grows faster in co-growth than separately, while the other grows slower in co-growth than separately. Outcome for faster-growing organism: Parasitism_Taker, for slower-growing organism: Parasitism_Giver

• Amensalism: one organism’s growth is unaffected by co-growth, while the other grows slower in co-growth than separately. Outcome for unaffected organism: Amensalism_Unaffected, for slower-growing organism: Amensalism_Affected

• Neutralism: both organisms’ growths are unaffected by co-growth (same outcome for both)

• Commensalism: one organism’s growth is unaffected by co-growth, while the other grows fatser in co-growth than separately. Outcome for unaffected organism: Commensalism_Giver, for slower-growing organism: Commensalism_Taker

• Mutualism: both organisms growth faster in co-growth than separately (same outcome for both)

Please note: the function takes the name of the folder containing previously created pairwise models as the input, e.g.: pairwiseInteractions, pairwiseSolutions] = simulatePairwiseInteractions(‘/Users/almut.heinken/Documents/PairwiseModels’);

USAGE:

[pairwiseInteractions, pairwiseSolutions] = simulatePairwiseInteractions(pairwiseModelFolder, varargin)

INPUTS:

pairwiseModelFolder: Folder where pairwise models and pairwise model info file are located

OPTIONAL INPUTS:

inputDiet: Cell array of strings with three columns containing

exchange reaction model abbreviations, lower bounds, and upper bounds. If no diet is input then the input pairwise models will be used with unchanged constraints.

saveSolutionsFlag: If true, flux solutions are stored (may result in

large output files)

numWorkers: Number of workers in parallel pool if desired sigD: Difference in growth rate that counts as significant

change (by default: 10%)

OUTPUTS:

pairwiseInteractions: Table with computed pairwise and single growth

rates for all entered microbe-microbe models

OPTIONAL OUTPUT:

pairwiseSolutions: Table with all computed flux solutions saved

(may result in large output files)

useDiet(modelIn, dietConstraints, printLevel)

Implements diet constraints in a COBRA model structure.

USAGE:

[modelOut] = useDiet(modelIn, dietConstraints)

INPUTS:

modelIn: original model dietConstraints: cell array of three columns containing

exchanges, lower bounds, and upper bounds respectively

printLevel: Verbose level (default: printLevel = 1)

OUTPUT:

modelOut: model with applied constraints