
Approximating Diffusion as a "Flux"From MathBioIntroductionWe are interested in approximating diffusion as a "flux" for the easy analysis of the biological system. Currently, using systems biology techniques, we are only able to apply fundamental subspace analysis to chemical reactions of a particular biological system, and not to the inevitable diffusive aspects. Mathematically, we have a biological model that follows the form:
,
.
Approximation Using the Second DerivativeFick's laws of diffusion provide us with a general basis for the construction of a matrix that includes information on the diffusion of a particular biological system. Namely, we have:
Using Fick's Second Law in particular, we see that the general form of the partial differential equations in a reactiondiffusion model is as follows:
.
,
We see that it becomes logical to approximate diffusion by examining the behavior at a finite number of evenly spaced points, using this limit definition for the second derivative in accordance with Fick's Second Law. Specifically, to approximate the concentration change at a particular point , we utilize information regarding the current concentration at three points: . Of course, in general, the approximation improves as the distance between the discrete points shrinks to zero, by definition. Here, we will use three points in the spatial domain, because this number is appropriate for our C5a production example. An extension to a greater number of approximation points within the spatial domain may be found here. Additionally, it may be helpful to view diffusion as effectively a system of chemical "reactions", where and the behavior of and is determined by the boundary conditions at those points.
Does This Approximation Work?Again, we are approximating diffusion at three discrete, evenly spaced points. Again, by the nature of diffusion and our model, this necessarily requires a point to the left and to the right of each point at which we are interested in approximating diffusion at. Then we have five discrete, evenly spaced points at which we require concentration information to apply Fick's laws to our approximation: 0, 1, 2, 3, 4 (Figure 1). From our second derivative approximation in the previous section, we have, in general, the following system of differential equations:
The question that remains is how to convert this information into stoichiometric matrix form and thereby conduct analysis via SVD decomposition and the like. We know the "reactants" that must be present in this simple diffusion matrix: the three "reactants" really correspond to a single reactant at our chosen three equallyspaced points within the spatial domain. We are concerned with how to define the "fluxes" in such a matrix. It is logical to separate fluxes by diffusive component  for example, and in the first differential equation above, corresponding to the separate diffusion gradients on the "left hand side" and the "right hand side" of the spatial domains relative to point 1. It is worth noting that we are not explicitly specifying whether there is flux "into" or "out of" the system in terms of diffusion by Fick's First Law in our stoichiometric matrix. That is, if (or ), a 1 should be entered into the stoichiometric matrix, but if (or ), a +1 should be entered into the stoichiometric matrix. Still, as these are effectively equal but opposite fluxes, we can simply modify them for whatever biological situation we are interested in. We are now interested in examining the effectiveness of our approximation for various boundary condition types. For our boundary conditions to "work", this approximation must blur the line between point and domain, such that the "boundary" conditions are taken across an entire area. Since our spatial domain we are interested in is from points 1 through 3, our left and right boundary conditions are applied in the space from point 0 to point 1 and point 3 to point 4, respectively. Therefore, in context of our diffusion approximation, we must specify "boundary" conditions as differences between concentrations at the two points specifying the edges of the boundary (not the spatial domain). For example, from point 0 to point 1 in the diffusion approximation for point 1, we define:
We will explore toy examples for different boundary condition scenarios:
by constructing diffusiononly matrices such that and exploring their biological validity via SVD decomposition and interpretation.
If we define Dirchlet boundary conditions at both the left and right boundaries of our spatial domain, we have:
We generate a diffusion matrix after having eliminated these fluxes at the boundaries. Note that we have three "reacting species" (really, the same reactant at three points in our spatial domain) and four fluxes, as two have been canceled by the Dirichlet boundary condition. But first, note that there are two equal and opposite fluxes included in our model; namely:
This eliminates two fluxes from our diffusion matrix. Thus, we can construct a matrix:
Let , corresponding to the number of rows of our matrix, and , corresponding to the number of columns of our matrix. We now perform SVD analysis on this matrix using MATLAB's SVD command, generating three matrices , with , where is an unitary matrix, is an diagonal matrix containing the singular values, and is an unitary matrix. Furthermore, using MATLAB's rank operation for gives a rank . Now, first, note that the dimension of the left null space is , which means that the left null space has one spanning vector. Through our SVD analysis, we see from the matrix that an orthonormal basis for the left null space is composed of the vectors , which in this case is just . This corresponds to a conserved quantity, where is time invariant, which is indeed the case, since in this approximation we have no net flux ath the left and right boundaries. Meanwhile, the dimension of the right null space is , meaning that the right null space has no spanning vectors. Biologically, there are no steadystate flux distributions, besides the trivial state in which the initial reactant concentrations are equal and no diffusive gradients form.
If we define a Dirchlet boundary condition at the left boundary of our spatial domain and a Neumann boundary condition at the right boundary, we have:
The presence of the constant flux at the right boundary gives us an additional flux in our system. Otherwise, this matrix is constructed analogously to the DirichletDirichlet case Thus, we can construct a matrix:
Let , corresponding to the number of rows of our matrix, and , corresponding to the number of columns of our matrix. MATLAB's rank operation for gives a rank . It turns out that we do not need the SVD decomposition: The dimension of the left null space is , which means that the left null space has no spanning vectors. Biologically, there are no conserved quantities. This makes sense, due to the existence of a flux on our right hand boundary. Note that in the case of the Neumann boundary condition that specifies zero flux, the system will look exactly the same as a Dirichlet boundary condition under this approximation, and so will have a conserved quantity of the sum of all three fluxes. Meanwhile, the dimension of the right null space is , meaning that the right null space has no spanning vectors. Biologically, there are no steadystate flux distributions, besides the trivial state in which the initial reactant concentrations are equal and no diffusive gradients form.
If we define Neumann boundary conditions at both the left and right boundaries of our spatial domain, we have::
with the subscripts 1 and 2 to distinguish the two potentially different fluxes at the left and right boundaries. We have added yet another flux to the system. Thus, we can construct a matrix:
Let , corresponding to the number of rows of our matrix, and , corresponding to the number of columns of our matrix. MATLAB's rank operation for gives a rank . The dimension of the left null space is , which means that the left null space has no spanning vectors. Biologically, there are no conserved quantities. This makes sense, due to the existence of a flux on both the left and right hand boundaries. Again, if the Dirichlet boundary condition specifies zero flux at both of these ends, the DirichletDirichlet conserved quantity arises. Meanwhile, the dimension of the right null space is , meaning that the right null space has one spanning vector. Indeed, through our SVD analysis, we see from taking the transpose of the matrix (not reproduced here) that an orthonormal basis is composed of the vectors , or in this case, . We reproduce this vector here using MATLAB's svd function: Remember that the Neumann boundary conditions were specified as +1 in the diffusive matrix, which corresponds to fluxes into the spatial domain. Biologically, our 'V' vector corresponds to the steady state distribution: This makes biological sense, as it corresponds to the situation in which each point on the spatial domain is balanced by equal and opposite diffusion fluxes. Specifically, this steadystate case really corresponds to the reactant in question diffusing across the entirety of the spatial domain at a uniform rate. 