Difference between revisions of "Team:HZAU-China/Model"

Model

The distribution density of bacteria in a culture medium depends on two factors: motion and reproduction of bacteria. For convenience, we will discuss these two factors separately.

Dynamic model of bacteria motion

In this project, we are trying to control the motion of bacteria by light signal to generate a colony of a specific shape. The light signal can affect the bacteria motility by influencing the expression of a protein cheZ. There are two motion forms of bacteria: tumbling and swimming. By tumbling, bacteria can change swimming direction but not position; by swimming, bacteria will move forward and change the position. In fact, the ratio of tumbling to swimming is different from one bacterium to another. At the microscopic level, the angle of tumbling is random, so we have no idea about the movement direction of a bacterium in any time. But from a macroscopic viewpoint, the probability of all directions are equal. This property is similar to the diffusion of a chemical molecule. So, bacteria swimming can be treated as the diffusion of bacteria and the rate of diffusion is related to the expression of protein cheZ.

Figure 1. The 2D plane of bacteria diffusion.

In Figure 1, the simple square represents the area of bacteria in coordinate \((x,y)\) approximately. The number of bacteria in this area is \(S(x,y,t)\) at time t. So, we have

\(S(x,y,t)=\rho(x,y,t)*\Delta x\Delta y\).      \((1)\)

In Equation \((1)\), \(\rho(x,y,t)\) represents the bacteria density in that area. \(\Delta x\) and \(\Delta y\) are the small increments along the \(x\) and \(y\) axes, respectively. Taking the derivative of Equation (1), we have

\(\frac {\partial S(x,y,t)}{\partial t}=\frac {\partial\rho(x,y,t)*\Delta x\Delta y}{\partial t}\).      \((2)\)

The left side of Equation \((2)\) represents the change rate of \(S(x,y,t)\). Ignoring the reproduction of bacteria, the change rate depends only on the rate of bacteria diffusion. As shown in Figure 1, bacteria can move along the \(x\) axis and the \(y\) axis. Suppose the movements along the \(x\) and \(y\) axes are positive, we arrive at

\(\frac {\partial\rho (x,y,t)\Delta x\Delta y}{\partial t}=\phi_x(x,y,t)\Delta y- \phi_x(x+\Delta x,y,t)\Delta y+\phi_y(x,y,t)\Delta x-\phi_y(x,y+\Delta y,t)\Delta x\)      \((3)\)

In Equation \((3)\), \(\phi_x\) and \(\phi_y\) represent the diffusion of bacteria in \(x\) and \(y\) axes, respectively. Dividing \(\Delta x*\Delta y\) on both sides of Equation \((3)\), and we can get

\(\frac {\partial\rho(x,y,t)}{\partial t}=\frac {\phi_x(x,y,t)-\phi_x(x+\Delta x,y,t)}{\Delta x}+\frac {\phi_y(x,y,t)-\phi_y(x,y+\Delta y,t)}{\Delta y}\)      \((4)\)

Let \(\Delta x\) and \(\Delta y\) take a limit to 0, Equation \((4)\) will become

\(\frac {\partial\rho(x,y,t)}{\partial t}=\lim_{\Delta x\rightarrow 0}\frac {\phi_x(x,y,t)-\phi_x(x+\Delta x,y,t)}{\Delta x}+\lim_{\Delta x\rightarrow 0}\frac {\phi_y(x,y,t)-\phi_y(x,y+\Delta y,t)}{\Delta y}\)      \((5)\)

This equation is equivalent to

\(\frac{\partial \rho}{\partial t}=-(\frac{\partial\phi_x}{\partial x}+\frac{\partial\phi_y}{\partial y})\)      \((6)\)

In this equation, \(\phi_x\) and \(\phi_y\) are unknown. In fact, the diffusion rates along \(x\) and \(y\) axes are the same. Refer to the assumption of Fourier Heat Equation:

1. If the temperature is constant within an area, there is no flow of heat.

2. If temperature difference exists in adjacent areas, heat will flow from areas of high temperature to low temperature.

3. For the same kind of material, the bigger the temperature difference between two adjacent areas, the faster the heat flow between them.

We have three similar characters for the diffusion of bacteria colony:

1. If the density of bacteria is constant in an area, the movement of bacteria will not lead to the change of density.

2. If density difference exists in adjacent areas, bacteria will swim from high density to low density areas.

3. The bigger the density difference between two adjacent areas, the faster the diffusion rate between them.

So, the rate of diffusion is related with the difference of density:

\(\phi_x = -k\frac{\partial\rho}{\partial x}\)      \((7)\)

\(\phi_y = -k\frac{\partial\rho}{\partial y}\)      \((8)\)

Combining Equation\((6)\), Equation\((7)\) and Equation\((8)\), we can get

\(\frac{\partial\rho}{\partial t} = k(\frac{\partial^2\rho}{\partial x^2}+\frac{\partial^2\rho}{\partial y^2})\)      \((9)\)

Equation (9) is the dynamic model of bacteria motion.

Dynamic model of bacteria reproduction

There is a famous model about the growth of bacteria under the conditions of limited space and nutrients, called logistic growth model:

\(\frac{d\rho}{dt} = \gamma\rho(1-\frac{\rho}{\rho_s})\)      \((10)\)

In Equation \((10)\), \(\rho\) represents the density of bacteria, and \(\gamma\) is the growth rate constant, and \(\rho_s\) is the saturated density.

Dynamic model of both reproduction and motion (R-M model)

Combining equation \((9)\) and equation \((10)\), we have

\(\frac{\partial\rho}{\partial t}=k(\frac{\partial^2\rho}{\partial x^2}+\frac{\partial^2\rho}{\partial y^2})+\gamma\rho(1-\frac{\rho}{\rho_s})\)      \((11)\)

Referring to the paper by Liu et al [Ref], we have the parameter values:

$$k=200~1000\mu m^2/s,$$ $$\gamma = 3.89e-4s^{-1},$$ $$\rho_s = 1500cell/\mu m^2.$$

We will use the Finite Element Method to solve this PDE. Initial conditions are defined as

$$\rho=matrix.zeros(200,200),$$ $$\rho[100,100]=100,$$

where \(matrix.zeros(200, 200)\) means a zero matrix with \(200*200\) elements and the initial value of the element at the index \([100, 100]\) of this matrix is \(100\). And the boundary conditions are

$$\rho[0:,]=0,$$ $$\rho[100:,0]=0,$$ $$\rho[0,0:]=0,$$ $$\rho[0,100:]=0.$$

With the initial and boundary conditions, Equation \((11)\) was solved by writing python program with Numpy and the result was visualized by the Mayavi python package, as shown in the following video (Movie 1).

Movie 1.

R-M model in a restricted area

In Movie 1, we can see that colony will be a roundness eventually. If we add a restrictive condition about the area that bacteria can move, then what the colony will become? The way to add area limit is using light to control the motility of bacterias. So we will give green light in a specific area and the remainder will be given red light. And bacterias can move only in this green light area. In equation (11), Bacterial motility is mirrored by parameter k. So the model will be the following equations.

\(\frac{\partial\rho}{\partial t}=k(\frac{\partial^2\rho}{\partial x^2}+\frac{\partial^2\rho}{\partial y^2})+\gamma\rho(1-\frac{\rho}{\rho_s})\)      \((11)\)

\(k=f(x,y)\)      \((12)\)

\(k\) in equation \((11)\) means the diffusion rate of colony. \(f(x,y)\) is area limit function, and it’s an image matrix. In this project, we choose Pikachu who is a famous role in a hot AR game Pokimon Go. See this picture in Figure2.

Figure 2.

The black part in Figure 2 means that the parameter \(k\) is normal value. That is \(k=k_{normal}=200~1000\mu m^2/s\). And the white part in Figure2 means that the parameter \(k\) is a small value which will be set to \(0\). But \(k\) is a matrix in equation \((11)\), and it can be solved out using Figure 2.

$$k = img*k_{normal}/255$$

Under the aforementioned initial conditions and boundary conditions, using python program to solve these equations.

Movie 2.

From Movie2 we can see that a specific pattern is basically formed. But such a pattern formation is just equivalent to using mould. In this project, the pattern formation will be adjusted by computer in real time.

R-M model with dynamic regulation

Area limit is static regulation, but dynamic regulation is more suitable in our expectation. In dynamic regulation, we will compare the shape of colony with target picture in real time, and the result of comparing will be converted to light matrix to be irradiate on colony. So, the dynamic model will like the following context.

\(\frac{\partial\rho}{\partial t}=k(\frac{\partial^2\rho}{\partial x^2}+\frac{\partial^2\rho}{\partial y^2})+\gamma\rho(1-\frac{\rho}{\rho_s})\)      \((11)\)

\(k=f(t,\rho,img)\)      \((13)\)

Unlike equation\((12)\), there are two additional parameters in equation\((13)\). The two additional parameters are time \(t\) and density \(\rho\) respectively. That means parameter \(k\) will be change in the process of adjustment. And equation\((13)\) can be interpreted as,

\(\rho_1=\rho*\frac{255}{max(\rho)}\)      \((14)\)

\(\rho_2=threshold(\rho_1,THRES\_BINARY\_INV)\)      \((15)\)

\(edges=Canny(\rho_2)\)      \((16)\)

\(edges_2=bitwise\_and(edges, mask=img)\)      \((17)\)

\(k_1=dilate(edges_2.kernel)\)      \((18)\)

\(k_2=bitwise\_and(k_1,mask=img)\)      \((19)\)

\(k=k_2*\frac{k_{normal}}{255}\)      \((20)\)

\((14)\rightarrow(20)\) is using OpenCV language. \((14)\) is transforming the matrix of colony to image. \((15)\) is image binaryzation. \((16)\) is finding the edge of image. \((17)\) is comparing the image above-mentioned with target image. \((18)\) is expansion the edge. \((19)\) is to present the overflow of the edge. \((20)\) will transform the result of \((19)\) to light message. Under the aforementioned initial conditions and boundary conditions, using python program to solve these equations. We can get it.

Movie 3.

You can download all programs in here