Simulation of the mycelium development

Fungal growth is initiated by $n$ number of spores, and a branch will start to extend from each of the spores added to space. Following along one of these branches originating from a single spore, the hyphae will grow in a direction $\theta$ with a tip extension $r_{tip, i}$. A branching event, in which a new branch is formed from the first branch, can occur with a probability $q$. Tip extension rate is calculated by using the equation below, which considers growth kinetics for the fungus and the amount of substrate available. It essentially outputs the accelerated growth dependent on the amount of substrate available, where the accelerated growth equation depends on fungal kinetics and branch lengths (1, 3). \begin{equation} r_{tip, i} = \bigg(k_{tip,1} + k_{tip,2}\cdot \frac{l_{br,i}}{l_{br, i} + K_t}\bigg)\cdot\bigg(\frac{S}{S+K_s}\bigg) \end{equation} $k_{tip,1}$ is the initial tip extension rate of the branch and $k_{tip,2}$ is the difference between the maximum extension rate and $k_{tip,1}$ (1, 3). $K_t$ symbolizes the time it takes to reach half of the maximum extension rate. The length of branch $i$ described by $l_{br,i}$. $S$ is the substrate concentration and $K_S$ corresponds to the substrate concentration to reach half of the maximum growth level (4). As the starting coordinates of this simulation are $(x_0, y_0)$ and end coordinates $(x, y)$, the length of branch $i$ can be calculated as the distance between two points: \begin{equation} l_{br, i} = \sqrt{(x-x_0)^2 + (y-y_0)^2} \end{equation} By dividing the growth area into a grid of $w\cdot h$ areas, it is possible to investigate the uptake of the substrate, hyphal movement and the development of biomass through the simulation.

At the start of the simulation, it can be assumed that the initial substrate concentration $S_0$ will be distributed evenly across the grid. When the spores are placed randomly, the branches start to develop by using the substrate available around the hyphal tip. So for each branch, it is checked which tip is in which grid and whether there is any substrate available. That means that one can track the substrate depletion in each grid, and when there is no more substrate in one of the areas of the grid, the hyphal tips located there can no longer grow and can be considered as inactive hyphae. If the total area around the mycelium no longer contains substrate, the growth will simply stop.

In the same way as tracking substrate depletion, biomass production can also be viewed in a grid. But instead of calculating a loss of substrate, each new branch extension is considered a gain in biomass or the density $d$. It should show the same mechanics as the substrate depletion, as these two are directly related to each other by the equation $l_{br,i}$.

As the fungus grows in size, the computational cost also increases. We were limited in running the simulation due to the power of our computers, thus resulting in us introducing two restrictions: the number of hyphae $M$ and the number of steps in the simulation $N$.

Substrate dependent growth

The fungus can grow and fill out the space since the substrate concentration is sufficient to support this pattern. If there is not enough substrate available, it is impossible for the fungus to fill out the space. This can be illustrated in the simulation if there is not enough substrate available. For example in fig. 3, the substrate concentration is a fifth of the level in fig. 6. As can be observed, most of the hyphae use up all the substrate available in their grid area, meaning they do not grow into another grid area as they are locked in.

The Model

The model presented below is based on Falconer (5), but greatly simplified because of the non-complex morphological behavior of our fungus. It should be noted that the full Falconer model was implemented, but all relevant effects can be recreated with the simplified model. \begin{equation} \frac{\partial b_i}{\partial t} = \zeta D_b \nabla^2 b_n \end{equation} \begin{equation} \frac{\partial b_n}{\partial t} = (1-\zeta) D_b \nabla^2 b_n \end{equation} \begin{equation} \frac{\partial n}{\partial t} = D_n \nabla^2 n +s(\lambda_1b_n+\lambda_2b_i) \end{equation} \begin{equation} \frac{\partial s}{\partial t} =-(\lambda_1 b_n+\lambda_2 b_i)s \end{equation} The four equations describe how the fungal biomass and the substrate change over time. The fungus is modeled as consisting of mobile and immobile biomass, the immobile biomass can be either non-insular ($b_n$), corresponding to the hyphal tips of the colony that are capable of large-scale extraction of external resources. The mobile biomass ($n$) represents the uptake of nutrients from the available substrate, and is used to regulate the diffusion rate $D_n$, which is lowered to $10^{-7}D_b$ if $n>n_0$, ie. the biomass exceeds a critical value: The parameter $\zeta$ is the conversion rate of non-insulated biomass to insulated biomass. The uptake efficiency of the mobile and immobile biomass are denoted respectively $\lambda_1$ and $\lambda_2$.

The system was spatially discretized on a 150x150 grid, which created a set of ordinary differential equations for each grid point (for a more thorough introduction see here). These were numerically integrated using a numerical $2^{nd}$ order strong stability preserving Runge Kutta method, using the package DifferentialEquations.jl in Julia (6). No flux boundary conditions were imposed, which is also known as a Neumann boundary condition, these are deemed realistic since we will be growing our fungus in molds. Alternatively, periodic boundary conditions could have been interesting if it was decided to make a spatial grid of inoculum. All calculations were vectorized and parallelized, which allowed us to simulate large systems on personal laptops.

Substrate Limited Growth

As described in the micro-model, the meso-model of course also limits the growth of the fungus based on the available substrate. This naturally occurs since any increase in fungal biomass is caused by a conversion of substrate to biomass. The model predicts a growth similar to the real measured growth of fungi (7), although the scaling of time and biomass are system dependent the overall trend is the same. These results are similar to those obtained in the micro-model, and shows that the model can also simulate realistic growth curves in the simple case. Note the asymmetry of the final growth, caused by an initial asymmetry of the initial spores grown in the micro-model.

Future Work

Expanding the model to 3-spatial dimensions should be fairly trivial, and would offer a better representation of how the density of mycelium changes as a function of time in the final molds. But the primary challenge is verifying that the density distribution is as predicted by the model, this could perhaps be achieved by measuring the absorbance of an agar plate with fungal growth. Furthermore, the parameters of the model have not been experimentally measured.