In genehackers/optizyme: Optizyme

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

Vignette Contents:

• List/Explanation of Modelling and Visualization functions available in Optizyme
• List/Explanation of inputs to these functions
• Explanation of Model System
• Demonstration of how to define the inputs to the functions in a way that is easily checked for errors
• Example call for each function

The Functions

The functions EnzymeModel, EnzymeModelPlot, and EnzymeModelTimeCourse all take the same inputs and differ in their outputs.

EnzymeModel outputs the terminal concentration of the species that is supposed to be optimized.

EnzymeModelPlot generates the graph for the time-course of the cell free system.

EnzymeModelTimeCourse outputs the time course for each species as a row in the matrix.

EnzymeModel can be used as the optimization function for MCGD, and the GuideToMCGD vignette explains how to do that.

Inputs

Their inputs will be described below, and then the use of each function will be detailed after.

NumberOfEnzymes is a scalar value that should be equal to the number of enzymes in the system.

InitialState is a vector containing the initial concentrations of each fo the species involved in the system. The last entry in InitialState should be a zero, and won't represent any species in the system. This last entry being zero is used for indexing purposes by the function.

Time is a scalar value that represents how long the simulated time course should represent. The units for time will depend on the units of the inputs to the function. For example, if turnover rates of enzymes are expressed in inverse seconds, the units of time will be seconds.

dt is a scalar value that represents the time step used in the integration of the differential equations representing the enzyme time course. As dt decreases, the computational cost increases, but accuracy of the model also increases.

EnzymeNumberOfSubstrates is a vector of size n, where n is equal to the number of enzymes. For example, if the first enzyme is a single substrate enzyme, then the first entry in EnzymeNumberOfSubstrates should be 1. If the second enzyme in the pathway is a bisubstrate enzyme, then the second entry in EnzymeNumberOfSubstrates should be 2. It is important to note that these functions do not support more than two substrates per enzyme.

EnzymeSubstrates is a matrix of size 2xn, where n is equal to the number of enzymes. Its entries are which species in the InitialState vector are substrates to the nth enzyme in the pathway. For example, if the first enzyme has one substrate and it is compound 1 listed in InitialState, then the first column of the EnzymeSubstrates matrix is 1,0. There must be two entries into the matrix for each column, and if one of the slots isnt used then the value should be 0. The first entry in the column shoul d be filled first. If the second enzyme has two substrates and they are compounds 2 and 3, then the second column of EnzymeSubstrates should be 2,3.

EnzymeProducts is also a matrix of size 2xn, with the same structure and properties as EnzymeSubstrates. However, the entries should represent the products of the enzyme rather than the substrates.

EnzymeConcentrations is a vector of size n. The first entry in the vector should be the concentration of the first enzyme in the pathway. The second entry is the concentration of the second enzyme, and so on.

EnzymeTurnoverRates is a vector of size n. It should be filled in the same manner as EnzymeConcentrations.

EnzymeMichaelisConstants is a matrix of dimensions 2xn. The entries in the matrix are the corresponding Michaelis constants for the enzymatic catalysis of each of the substrates. For example, in our example for EnzymeSubstrates, we said that the first column of EnzymeSubstrates was 1,0. The corresopnding column of EnzymeMichaelisConstants should be Km,1, where Km is the Michaelis constant for the catalysis of compound 1 by enzyme 1 in this example. If a Michaelis Constant does not exist because an enzyme is single substrate, the second entry in the column should be a 1. However, any nonzero constant will work here.

EnzymeS1CompetitiveInhibitors is a matrix of size 10xn, where n is the number of enzymes in the system. Each column of the matrix contains the competitive inhibitors with the first substrate listed in the corresponding EnzymeSubstrates matrix. If there are unused slots, they should be filled in with the final entry of InitialState (Assuming that the final entry is the entry that doesn't represent a species and is just zero). For example, if InitialState is a vector of size 12, then unused slots in the columns of EnzymeS1CompetitiveInhibitors should have the value 12. If compounds 5,6,7 are competitive inhibitors of compound 1 for enzyme 1, then the first column of the matrix is <5,6,7,12,12,12,12,12,12,12>, where 12 is the terminal entry of InitialState.

EnzymeS1CompetitiveInhibitionConstants is a matrix of the same dimensions as EnzymeS1CompetitiveInhibitors, and contains the corresponding competitive inhibition constants for each of the competitive inhibitors with the first substrate binding site. The competitive inhibition constants associated with unused competitive inhibitors should be arbitrary nonzero constants, it is recommended to use zero.

EnzymeS2CompetitiveInhibitors and EnzymeS2CompetitiveInhibitionConstants are structurally the same as their S1 counterparts, but the inputs differ based on whether the compounds are competitive inhibitors with the first substrate or the second substrate listed in the EnzymeSubstrates matrix.

EnzymeNonCompetitiveInhibitors and EnzymeNonCompetitiveInhibitionConstants are also matrices of size 10xn. They should be filled in the same manner as the competitive inhibition matrices.

Model System

In order to best demonstrate the workflow involved in using the Modelling and Visualization capabilities of Optizyme, an example system will be used, and the inputs used for the modelling and visualization of that system will be explained.

The system used in this example is the cell-free system described in the following paper: A combined experimental and modelling approach for the Weimberg pathway optimisation

Code to Define the Inputs

Here is some code that defines the inputs needed to use these functions. Note that many of the inputs are not directly entered into the function call itself. Instead, they are defined as separate variables (for example, InitialState), and the names of those variables will be entered into the function call. The function would work just as well if the vectors were entered directly into the functin call, but doing it this way helps to check for errors.

Another thing to note here is that instead of turnover rate and concentration, the entries to EnzymeTurnoverRates and EnzymeConcentrations are specific activies in U/mg and then concentrations in mg/ml. Within this specific model, there are two reactions that Optizyme does not represent exactly like the model in the paper. There is a spontaneous conversion of compound 2 to compound 3 whose rate equation doesn't follow the structure of enzyme kinetic differential equations, and so it is left out of the model. Furthermore, the NAD+ recycling mechamism described in the paper is approximated as an enzyme reaction, but is not represented perfectly.

InitialState<-c(4.25,0,0,0,0,0,11.5,0,0,10,0)
EnzymeNumberOfSubstrates<-c(2,1,1,1,2,2)
EnzymeSubstratesFilling<-c(1,7,2,0,3,0,4,0,5,7,10,8)
EnzymeSubstrates<-matrix(EnzymeSubstratesFilling,ncol=6)
EnzymeProductsFilling<-c(2,8,3,0,4,0,5,0,6,8,0,11)
EnzymeProducts<-matrix(EnzymeProductsFilling,ncol=6)
EnzymeConcentrations<-c(2,1.5,5.7,9.8,3.6,1)
EnzymeTurnoverRates<-c(.123697,.839404,.043636,.019801,.0583,1000)
EnzymeMichaelisConstantsFilling<-c(0.19850960840785867,0.16123927045886804,0.4452685845967363,1,0.7939547312720485,1,0.20765619156952056,1,0.021750105300746496,0.5962856010942978,10,10)
EnzymeMichaelisConstants<-matrix(EnzymeMichaelisConstantsFilling,ncol=6,nrow=2)
E1S1CI<-c(2,9,9,9,9,9,9,9,9,9,9)
E2S1CI<-c(3,9,9,9,9,9,9,9,9,9,9)
E3S1CI<-c(4,9,9,9,9,9,9,9,9,9,9)
E4S1CI<-c(5,9,9,9,9,9,9,9,9,9,9)
E5S1CI<-c(6,4,9,9,9,9,9,9,9,9,9)
E6S1CI<-c(9,9,9,9,9,9,9,9,9,9,9)
EnzymeS1CompetitiveInhibitors<-cbind(E1S1CI,E2S1CI,E3S1CI,E4S1CI,E5S1CI,E6S1CI)
E1S1C<-c(0.5357339888880199,9,9,9,9,9,9,9,9,9)
E2S1C<-c(0.038100367072300925,9,9,9,9,9,9,9,9,9)
E3S1C<-c(0.8624275767674491,9,9,9,9,9,9,9,9,9)
E4S1C<-c(0.28865820908799766,9,9,9,9,9,9,9,9,9)
E5S1C<-c(0.27905898103297405,0.2136,9,9,9,9,9,9,9,9)
E6S1C<-c(9,9,9,9,9,9,9,9,9,9)
EnzymeS1CompetitiveInhibitionConstants<-cbind(E1S1C,E2S1C,E3S1C,E4S1C,E5S1C,E6S1C)
E1S2CI<-c(8,9,9,9,9,9,9,9,9,9,9)
E2S2CI<-c(9,9,9,9,9,9,9,9,9,9,9)
E3S2CI<-c(9,9,9,9,9,9,9,9,9,9,9)
E4S2CI<-c(9,9,9,9,9,9,9,9,9,9,9)
E5S2CI<-c(8,9,9,9,9,9,9,9,9,9,9)
E6S2CI<-c(9,9,9,9,9,9,9,9,9,9,9)
EnzymeS2CompetitiveInhibitors<-cbind(E1S2CI,E2S2CI,E3S2CI,E4S2CI,E5S2CI,E6S2CI)
E1S2C<-c(0.02970905659171423,9,9,9,9,9,9,9,9,9,9)
E2S2C<-c(9,9,9,9,9,9,9,9,9,9,9)
E3S2C<-c(9,9,9,9,9,9,9,9,9,9,9)
E4S2C<-c(9,9,9,9,9,9,9,9,9,9,9)
E5S2C<-c(0.2674,9,9,9,9,9,9,9,9,9,9)
E6S2C<-c(9,9,9,9,9,9,9,9,9,9,9)
EnzymeS2CompetitiveInhibitionConstants<-cbind(E1S2C,E2S2C,E3S2C,E4S2C,E5S2C,E6S2C)

E1NCI<-c(9,9,9,9,9,9,9,9,9,9)
E2NCI<-c(9,9,9,9,9,9,9,9,9,9)
E3NCI<-c(8,9,9,9,9,9,9,9,9,9)
E4NCI<-c(1,6,10,11,9,9,9,9,9,9)
E5NCI<-c(9,9,9,9,9,9,9,9,9,9)
E6NCI<-c(9,9,9,9,9,9,9,9,9,9)
EnzymeNonCompetitiveInhibitors<-cbind(E1NCI,E2NCI,E3NCI,E4NCI,E5NCI,E6NCI)
E1NC<-c(9,9,9,9,9,9,9,9,9,9)
E2NC<-c(9,9,9,9,9,9,9,9,9,9)
E3NC<-c(10.462860450236784,9,9,9,9,9,9,9,9,9)
E4NC<-c(18.30014950694962,14.791937154740701,17.91684582661257,27.97615617816409,9,9,9,9,9,9)
E5NC<-c(9,9,9,9,9,9,9,9,9,9)
E6NC<-c(9,9,9,9,9,9,9,9,9,9)
EnzymeNonCompetitiveInhibitionConstants<-cbind(E1NC,E2NC,E3NC,E4NC,E5NC,E6NC)

##

Now that the inputs have been defined, we will show examples of using these inputs to call these three functions.

##

EnzymeModel

EnzymeModel outputs the terminal concentration of the species that is supposed to be optimized. It multiplies the terminal concentration because our algorithm is a minimization, so maximizing yield corresponds to minimizing yield*-1.

EnzymeModel(6,InitialState,63.07,.01,EnzymeNumberOfSubstrates,EnzymeSubstrates,EnzymeProducts,EnzymeConcentrations,EnzymeTurnoverRates,EnzymeMichaelisConstants,EnzymeS1CompetitiveInhibitors,EnzymeS1CompetitiveInhibitionConstants,EnzymeS2CompetitiveInhibitors,EnzymeS2CompetitiveInhibitionConstants,EnzymeNonCompetitiveInhibitors,EnzymeNonCompetitiveInhibitionConstants,6)

EnzymeModelPlot

EnzymeModelPlot generates the graph for the time-course of the cell free system.

EnzymeModelPlot(6,InitialState,63.07,.01,EnzymeNumberOfSubstrates,EnzymeSubstrates,EnzymeProducts,EnzymeConcentrations,EnzymeTurnoverRates,EnzymeMichaelisConstants,EnzymeS1CompetitiveInhibitors,EnzymeS1CompetitiveInhibitionConstants,EnzymeS2CompetitiveInhibitors,EnzymeS2CompetitiveInhibitionConstants,EnzymeNonCompetitiveInhibitors,EnzymeNonCompetitiveInhibitionConstants,6)

EnzymeModelTimeCourse

EnzymeModelTimeCourse outputs the time course for each species as a row in a matrix. This function can be useful since it allows one to view the data in a different way than with EnzymeModelPlot.

In order to correctly format the output of EnzymeModelPLot, the InitialState vector must be extended for n entries in each row based on the Time and dt chosen. Below is some code that does that:

Time<-63.07 dt<-.01 NumberOfSteps<-Time/dt

Integers<-(1:(NumberOfSteps)) Time<-dt*Integers rainbow<-rainbow(length(TimeCourseMatrix[,1]))

Once the inputs have been defined, the function can be called:

EnzymeModelTimeCourse(6,InitialState,63.07,.01,EnzymeNumberOfSubstrates,EnzymeSubstrates,EnzymeProducts,EnzymeConcentrations,EnzymeTurnoverRates,EnzymeMichaelisConstants,EnzymeS1CompetitiveInhibitors,EnzymeS1CompetitiveInhibitionConstants,EnzymeS2CompetitiveInhibitors,EnzymeS2CompetitiveInhibitionConstants,EnzymeNonCompetitiveInhibitors,EnzymeNonCompetitiveInhibitionConstants,6)

If you want to plot the output of EnzymeModelTimeCourse, you can do it by defining a variable to be equal to the output of EnzymeModelTimeCourse, like so:

TimeCourseMatrix<-EnzymeModelTimeCourse(6,InitialState,63.07,.01,EnzymeNumberOfSubstrates,EnzymeSubstrates,EnzymeProducts,EnzymeConcentrations,EnzymeTurnoverRates,EnzymeMichaelisConstants,EnzymeS1CompetitiveInhibitors,EnzymeS1CompetitiveInhibitionConstants,EnzymeS2CompetitiveInhibitors,EnzymeS2CompetitiveInhibitionConstants,EnzymeNonCompetitiveInhibitors,EnzymeNonCompetitiveInhibitionConstants,6)

plot(Time,TimeCourseMatrix[1,],type="l",col=rainbow[1],xlab="Time",ylab="Concentration",main="Concentrations Over Time") lines(Time,TimeCourseMatrix[2,],col=rainbow[2]) lines(Time,TimeCourseMatrix[3,],col=rainbow[3]) lines(Time,TimeCourseMatrix[4,],col=rainbow[4]) lines(Time,TimeCourseMatrix[5,],col=rainbow[5]) lines(Time,TimeCourseMatrix[6,],col=rainbow[6]) legend("bottomright",legend=c("S1","S2","S3","S4","S5","S6"),col=rainbow,title="Compounds",lty=1,cex=.5)

More species can be plotted as needed by adding more "lines" commands and labels to the legend.

genehackers/optizyme documentation built on Oct. 24, 2020, 9 a.m.