In our cohort study we found that modern pharmaceutical treatment strategies against IBD target the tumor necrosis factor-alpha (TNFα). But the only drugs available so far act on the human body in a systemic manner. (1) We propose a new treatment that combats inflammation specifically in the gastrointestinal tract. We engineered E. coli Nissle 1917 (2) to sense nitric oxide (NO) as an inflammation signal and consequently release nanobodies which target TNFα. This allows the treatment to be localised to where it is needed.
The human gastrointestinal tract and its microbiota are very complicated systems that are difficult to simulate in laboratory conditions (3). Moreover, clinical trials exceed the scope of an iGEM project, so we opted to use a computational model to give us insights into the necessary conditions for an effective treatment. With our model, we want to test the effectiveness of our treatment depending on how many nanobodies the bacteria can produce and how severe the inflammation is.
The gut environment of IBD patients is characterized by eroded mucus with several localized inflammation sites on the gut surface (4). The lumen is continuously flushed, and anything that is not attached to the gut wall can largely be disregarded (5). In particular, our bacteria would be able to hold onto the gut only when they can penetrate the protective mucus and temporarily colonize(6). Therefore, we propose a 2D simulation of the gut. With the bacteria and the inflammation sites as the two base layers of our construction, we simulate the respective emission of molecules. The constant emission of NO levels at inflammation sites is relevant for out simulation. The larger and less mobile cytokine, TNFα, is also emitted at these sites and is the target for our bacteria.
The bacteria emit nanobodies (Nb) when sensing adequate levels of NO. We assume no division of bacteria and immediate sensing of NO to Nb production. This assumption is made as the time scales are negligible and because we are interested in the effect of the treatment after colonisation and NO detection. Over time we let the emitted particles diffuse, decay, and reemit, in accord with our research.
We model interactions between Nb’s and TNFα based on our lab data, as well as the relationship between reduced inflammation and NO production, to gather insights into the overall reduction of inflammation in the gut environment. As many variables have not yet been uncovered, we simulate a broad range of possibilities from worst to best-case scenarios to gather as much information as possible. In the end, we want to be able to test how overall TNFα levels on the grid surface can be reduced based on varying values for the number of inflammation sites and bacterial nanobody production. This knowledge can be used to give insights into the engineering successes achieved by the wet lab team, also showing which parts need to be optimized in order to bring the treatment closer to an in vivo application.
These are some samples for different settings in our simulation. Three of the parameters can be adjusted with the sliders, to produce different gifs of the processes. In the actual code there are many more parameters to adjust and any possible value can be inserted.
Parameter | Low | Medium | High |
---|---|---|---|
Inflammation sites | 15 | 25 | 100 |
Nb production | 0.1 * Standard | Standard | 10 * Standard |
Sensing Threshold | 0.25 * Standard | Standard | 4 * Standard |
Low
-
Inflammation sites
Nanobody production
Sensing Threshold
High
+
A section of the gut surface is represented as a square matrix with N rows and N columns, where N is a natural number. Each element represents an area of 1 µm². Thus, the whole grid represents an area of N² µm². Grid areas with varying sizes are being randomly assigned the status “inflamed”.
Grid elements with the inflamed condition emit NO and the cytokine TNFα. Bacteria are also randomly placed over the grid with a probability that can be varied in the input parameters. They can sense a certain threshold of NO, resulting in the production of nanobodies. Since E. coli has a length of roughly 2 µm, and a diameter of about 1 µm(9), we assume that on our grid, a single element can contain, at maximum, one bacterium.
Suppose the grid element containing a bacterium also contains an NO concentration above the sensing threshold, the bacteria produce nanobodies, leading to a rise of Nb concentration in the grid element. NO, TNFα and nanobodies are all referred to by us as particles. Each particle is affected by physical processes like diffusion and will decay over time. Nanobodies additionally bind to the receptors of TNFα if they are close enough to each other.
TNFα with nanobodies bound to its binding sites cannot attach to the receptors on the gut cells, which we simulate by removing TNFα from the system. In our model, the maximal distance to make this bond possible is determined by the location on the grid. If there is a nonzero concentration of Nb and TNFα inside the same grid element, TNFα and Nb cancel each other out inside that grid element according to their concentration, respectively. Because nanobodies diffuse relatively fast (used diffusion coefficient: 40 µm²/s), the resolution should be good enough to simulate proximity.
Additionally, the Nb concentration is reduced by the same number of deleted TNFα. TNFα has multiple binding sites, and we do not know how many of them have to be blocked in order to consider the particle to be inactivated. Thus we can also change the ratio in which TNFα gets canceled compared to the present nanobodies. To indicate the number of particles inside a grid element, we need a unit that shows the relation of particles per area because we work with a 2D model.
We chose mol/µm² as our standard unit to indicate concentrations. In this way, we can put the amount of TNFα and nanobodies into a relation that includes their functional effects on each other. In addition, the unit has the same spatial dimension as our model. Concentrations inside a grid element can move into neighboring grid elements within a timestep of one second via diffusion.
In every model that has multiple processes that change over time, the order of operations is essential. There are some operations that can be interchanged without changing the outcome, while some others have profound impacts.
Diffusion and enforcement of half-life are both processes that affect all soluble particles. The diffusion and the half-life are processes that work in proportions, so they are interchangeable in time, with results remaining the same. However, as the molecule production is additive and based on a fixed value, changing its place in the operation pipeline will alter the results.
The same goes for the binding and cancellation of nanobodies and TNFα. In an ideal case, the durations of the processes approach zero and repeat infinitely often. In reality, this is not feasible, as the computational needs are constant for each distinct process, regardless of the time interval. Increasing the number of cycles, even if the simulation period remains the same, thus increases the computational needs. The result is a need to balance computational cost and accuracy. To compensate for the losses in accuracy, we deduced a balancing mechanism for particle diffusion, decay, and emission.
The decay is grouped with the diffusion not only to compensate for a loss in accuracy, but also because it is the more conservative way to simulate. Nb's and TNFα are placed between the processes of emission and cancellation, as the binding process needs time and will not happen the moment the Nb's are created. The decay of Nb is faster than TNFα, so simulating the decay first will make our model less efficient in disposing of TNFα than the other way around.
The production of new particles is the first process in each timestep. TNFα and NO are produced by inflamed tissue. Therefore, the concentration of these two particles will be increased in every grid element with the status “inflamed”. The increase in concentration is based on literature values and is fixed for every particle.
Nanobody production is bound to two conditions: if the grid element contains bacteria and the concentration of NO is above the sensing threshold of 8.36 * 10-11 mol/µm², the grid element’s Nb concentration increases by a fixed value. This value comes from a similar sensing system (10).
After the particles are produced, they can diffuse away from the location at which they were produced. For this, we use the Heat Equation . Through discretization of this partial differential equation (F is equal to the diffusion coefficient multiplied with the difference in time divided by difference in distance), we can simulate diffusion in 1D with the backwards euler scheme. (11)
In 1D, the diffusion matrix has the same dimensions as our grid. The column and row indeces of the diffusion matrix can be viewed as the coordinates on our grid over which substances diffuse. A random walker can distribute a concentration inside a vector to the neighboring elements. Given the concentration within a grid element, we can estimate the local and neighbouring concentrations after simulating one timestep. This estimation is conducted by the mesh Fourier number (F). Each substance has its own diffusion coefficient and needs its own diffusion matrix. The main diagonal of the diffusion matrix contains the value S, which stands for “stay” and is equal to 1+2F. The other diagonal contains the value M, which stands for “move”, and is equal to -F. This is where the particles are going to diffuse to.
To get from this one-dimensional diffusion to a two-dimensional one, we multiply our grid with the inverse of this diffusion matrix. This gives us a diffusion along one axis. Then we multiply the inversed diffusion matrix again with the diffused grid to get diffusion along the other axis. A more accurate method to simulate 2D diffusion uses the Crank-Nicolson scheme paired with the Runge-Kutta scheme(12). Whereas diffusion in the first method happens through two separate matrix multiplications, this method only needs one single matrix-vector multiplication. This allows simultaneous diffusion in all directly neighboring grid elements (left, right, up, down). However, the method requires a diffusion matrix of size N²*N². For large grid sizes, diffusion matrices become very large and computationally demanding. To compare both methods, we simulated 30 timesteps and compared the resulting diffusion pattern. For this purpose, we used a grid size equal to 100µm² and an arbitrary diffusion coefficient of 20 µm²/s, and a starting concentration of 1 mol/µm² at index x=50 and y=50.
The approximation results in indistinguishable diffusion patterns and enables us to reduce computational needs significantly.
When implementing the diffusion, we need to consider that we have discrete time steps. First particles are emitted, then diffused. If we scale this operation up by multiplying the strength of emission and diffusion, we need to compensate for the effects. Let's consider going from the time step of n seconds to n*10 seconds. In the n-scenario, if we simulate 10 time steps, we introduce particles 10 times, and diffuse 10 times. This leads to the first emission being diffused 10 times, the second 9 times and so on until the last emission is diffused exactly once. In the n * 10-scenario however, we emit all particles at once, then do the diffusion 10 times, which results in all particles being diffused 10 times. To counteract this, we model the diffusion process in the following way:
A represents our particle matrix. X replaces z as the emission parameter as described below. M is our diffusion matrix and d our decay parameter. The number of repeats is described by n. This process should follow the same process as the n-scenario, where 1/10th of the particles diffuses ten times, 1/10th diffuses nine times, and so forth. When implemented, this process allows us to calculate the sum of the diffusion matrix and decay once and use it for all further calculations. This way, the computational cost for a single time step is invariant to the time span covered.
This parameter indicates how many countable inflammation sites with varying sizes are being generated. Since the amount and size of inflammation sites can vary strongly between patients and because there is no data available on human IBD patients this parameter could take a broad range of values. We decided to assign this parameter a random value between 20 and 100 inflammation sites. Inflammation sites are randomly generated and can overlap.
We chose a diffusion coefficient of 2000 µm² per second for NO based on current literature(13). For TNFα we chose a diffusion coefficient of 0.16 µm² per second(14). For nanobodies we chose a diffusion coefficient of 40 µm² per second based on diffusion of antibodies(15) because diffusion coefficients on nanobodies were hard to find and we are also constructing bivalent nanobodies which are larger.
The emission of NO- and TNFα particles should be high enough for the density of particles to constantly be at the literature values by compensating for the losses. We can simply subtract the losses from the desired values and add the difference. However, if we scale up, we need to make sure that we introduce enough particles at once to bridge the period during a process in which the concentration will be significantly lower than average. If the particles deplete entirely, going from a time scale of n to n*10 could mean introducing 10 times as many particles at once. To stay faithful to the processes happening, we wanted to find a better estimation.
Let us compare two trials of taking water from a bucket. In the first trial, the bucket contains z liters of water. We take a proportion of k liters from the bucket and transfer it to a second bucket. We fill the first bucket back up to the levels of z and repeat this process for n times. At the end, the liters in our second bucket equal:
In the second trial, we want to fill up our first bucket to a value of x liters, so that we reach the same value as in the first trial, if we take the proportion k, n times, without refilling. The first time we take water, we get k ∗ x liters. The second time, only k ∗ (1 − k) ∗ x, as there is less water remaining. If we repeat this process, our total amount of water equals:
By combining 1 and 2 and solving for x, we get:
Assuming that our model has constant levels of NO, we can think of diffusion and refilling to happen at infinitely small time steps. To translate this into our model, we want to calculate the value x that our model would need to assume for the NO values at the start of diffusion. We can do so based on the loss through diffusion and refilling at very small time scales.
Each grid element that is part of an inflammation site emits a fixed amount of NO and TNFα. Each grid element that contains an NO concentration above the detection threshold and contains a bacterium emits a specified concentration of nanobodies.
Since there is no data on NO concentrations around inflammation sites in the gut, we derived this value from the medium NO concentration in blood serum samples of IBD patients. This NO concentration should lie between 14.54 μmol/L and 15.25 μmol/L (16). We used 15 μmol/L and changed it into our standard unit to get a value of 1.2349 *10-10mol/µm².
Similarly, we have no data on TNFα concentrations in the gut. Therefore, we derived a value for TNFα from blood serum samples based on the TNFα concentrations in UC patients, which lie around 8.3 ± 2.5 pg/ml and 5.4 ± 1.7 pg/ml (17). We chose 5.4 pg/ml and changed it into our standard unit for better comparison, which resulted in a value of 3.84*10-17 mol/µm².
Every particle that leaves its generative environment through diffusion will eventually decay. To simulate this, we used the respective half-lives of each particle. For NO, we found the literature value of 0.09 to 2 seconds half-life in extravascular tissue(18).
For TNFα, we found parameters exploring the half-life of TNFα from intravenous injections in rats. The researchers found near dose-independent decay of around 30 minutes half-life in the high-dosage conditions, which would be similar to the elevated values observed in the gut of IBD patients (19).
No experimental value was available for the specific nanobodies we used in our experiments. However, as nanobodies are structurally similar, we think that the general decay times of nanobodies should simulate the real conditions reasonably well, leading to a half-life of about 12 minutes (20).
By doing simulations with our model, we want to determine whether varying amounts of inflammation sites influence their detection. If so, in which range does inflammation site amount significantly affect the treatment course? To test this, we used a grid size of 1000 µm², a probability for a single element to contain a bacterium of 0.0001 and an inflammation size range of 20 to 100 µm². We chose 1000 seconds for each simulation. We let the simulation run ten times, each time increasing the inflammation site amount by 30 inflammation sites and starting the first simulation with 30 inflammation sites. Because there are some randomized mechanisms in the simulation, like placement and size of inflammation sites and bacterial placement, we decided to repeat this experiment 3 times. Next, we want to find out whether it makes a difference in how many nanobodies a single bacterium can produce and secret. For that, we conducted a similar experiment where we now fixed the number of inflammation sites to 50. We are going to increase the parameter denseNB, which defines the nanobody concentration a single bacterium is able to produce when sensing adequate levels of NO. We started with a concentration of 1.5*10-20mol/ µm² and increased it by the same value each time.
Figure 7 shows TNFα concentrations plotted against time. With more inflammation sites being created, the total TNFα concentration on the surface at timestep zero increases. All TNFα concentrations decrease with time. We did an ANOVA to test whether the simulations with different amounts of inflammation sizes differ from each other (p-value = 0.0035). To see how NO detection worked in the different runs, we plotted the nanobody concentration over time in figure 8. The curve for the inflammation site amount of 30 showed a lower nanobody production. We performed another ANOVA here to see if nanobody production is significantly different depending on the number of inflammation sites (p-value = 4.039*10-208). When dealing with higher numbers of inflammation sites, further increasing the inflammation sites does not make a large difference, as seen from a t.test between inflammation sites = 210 and inflammation sites = 300 (p-value = 0.34). To see if bacterial nanobody productivity influences the TNFα cancellation, we plotted the TNFα concentration against time from the experiment with a fixed inflammation site amount. The course of TNFα reduction looks very clustered and similar on the plot in figure 3. An ANOVA reveals that the simulations with different nanobody productions are not significantly different from each other (p-value = 0.8).
The results show that in our model, the amount of nanobody production of a single bacterium does not affect the TNFα reduction significantly. The number of inflammation sites seems to have a much stronger effect on the TNFα reduction. This could be because inflammation site sizes can vary strongly. Like that, the total area of inflamed tissue can theoretically be higher for a lower amount of inflammation sites by chance. It could also be that our grid is too small for the effect to be visible.
Experiment two did not match our expectations. We would have expected for high production of nanobodies to get a significantly higher reduction in TNFα. This could be because an initial intense burst in Nbs could cause the local inflammation to decrease quickly, reducing the amount of Nbs produced in the long term, as they will not get sufficient concentrations of NO to trigger the production.
When fewer Nbs are produced, there is more time where the local NO threshold enables Nb production. This result indicates that the Nb production levels are less relevant than initially assumed, and further engineering procedures should focus on increasing the survival of bacteria and the chances of adhering to the mucus. In this way a larger surface area of the gastrointestinal tract can be covered and treated.
When more information is collected about IBD in the gut, the parameters can be updated accordingly. The model we have built can also be easily adapted to work with any system that requires bacterial sensing and emission. We took special care to make our systems as close as possible to continuous time, writing functions for diffusion and emissions that could be used by anyone who wants to simulate these processes, even outside of synthetic biology.
Even though we have focused on simulating as many mechanisms as possible, we could have improved some with more time and resources. The first thing we want to point out is the cancellation of TNFα by nanobodies, which we simulate with ratios. This process would probably be more of a sequential nature because a single TNFα molecule with only one blocked binding site would still partly fulfill its function in the cell tissue until the other binding sites are blocked.
TNFα molecules with bound nanobodies would also still diffuse and affect the diffusion of the active TNFα molecules. In a more advanced diffusion simulation, particles could also affect the movement of different kinds of particles. For example, TNFα concentrations could then affect how NO and nanobodies diffuse. A last big improvement of the model could be adding bacterial processes to the model. The movement and attachment of the bacteria, as well as their proliferation and death, could be modeled. However, the temporal orders of magnitudes on which these processes happen greatly exceed our simulation times.
We simulated the gut environment and all relevant particles. The backbone of our model is the emission-diffusion system that propagates particles over space and time. We can manufacture a wide range of scenarios by implementing various parameters such as decay times, diffusion constants, particle concentrations, and severeness of inflammation. The model lets us navigate the very hard-to-explore and still vastly unknown disease, IBD. By keeping the model general and source code available, we allow for adaptation to different illnesses, and certain aspects like diffusion, decay and emission could be used in various models in the future.