Team:DTU-Denmark/Model

• Summary
• The mycelium
• Goal and design
• Assumptions
• Morphologizer
• Mycemulator
• Results
• Back to top

• Summary
• The mycelium
• Goal and design
• Assumptions
• Morphologizer
• Mycemulator
• Results

Summary

To aid our parts characterization, we developed two models for studying the morphological features of Aspergillus niger, both based on our experimental data.

The first model, Morphologizer, automatically analyzes microscopic images of mycelia that we captured with our novel microscopy protocol. By treating the mycelium as a graph, we can extract parameters that describe features of interest such as branching and chemotropism.

The second model, Mycemulator, is a stochastic model that simulates mycelial growth based on morphological parameters. Some of these parameters were estimated by Morphologizer model, some were from our fermentation experiments, and a few were based on literature references. To better reflect industrial conditions, Mycemulator also models substrate diffusion towards the center of the mycelium.

Several of the features we have introduced were not found in literature, and can be considered novel in the field of mycelium modelling. For example the experimentally estimated distribution of hyphal curvatures used to model chemotropic behaviour, or the polar grid used for substrate concentrations.

Together with our novel microscopy protocol, the models present an entirely new, streamlined method for estimating mycelial morphology, and has been well-received by academic and industry scientists alike. To better illustrate our work, we present a showcase of each of our models below.

Model showcase: Morphologizer

Showcase of our first model, Morphologizer. The model analyzes microscope images of mycelia, converts them into graph objects, and extracts morphological features of interest from the graph, such as the number of branches, branching angles and hyphal curvature.

Model showcase: Mycemulator

Showcase of our second model, Mycemulator. The model simulates mycelial growth, and is able to use parameters and distributions estimated by Morphologizer. The green field is substrate availability and is modelled in polar coordinates to accommodate the spherical shape of fungal mycelia.

The mycelium

As we describe in our project description, the morphology of filamentous fungi has a large impact on its productivity in industrial biotechnology. When we talk about morphology, we are more specifically talking about the morphology of the vegetative mycelium, a branching network of filaments or hyphae. While most filamentous fungi have various interesting behaviors and multicellular structures, such as their spore-forming structures, the mycelia was our focus, as that is what we mainly find inside the fermentation tanks. One of our mycelium images of Aspergillus niger can be found in the figure below.

Click the image to zoom in!

An example of a vegetative mycelium from the filamentous fungi Aspergillus niger (click to zoom in). Picture was taken by staining cell walls with Calcofluor White and measuring fluorescence at 410 nm. Colours are artificial, original is grayscale. Captured using our novel microscopy protocol.

Basic mycelium behaviour

While the mycelium is a complex multicellular structure that can have varying behaviors from species to species, some basic patterns are observed across all species, including A. niger. The mycelium begins with a spore, from where several hyphae (filaments) start growing. These hyphae extend at the tips, and can form branches occasionally. If you are not familiar with the details of how the vegetative mycelium grows, check out the figure below.

Illustration of basic behaviour in a mycelium, the branching network of hyphae (filaments) in fungi. These are the behaviours that our modelling efforts focus on, with the exception of how the spore itself is generated.

Modelling goal

The main goal of our modelling work is to analyze and predict the growth and morphology of fungal mycelia. We accomplish this through a pipeline of two separate models, Morphologizer and Mycemulator. Simply put, we capture images with our microscopy protocol, analyze the images with Morphologizer which outputs morphological parameters, and then we simulate mycelial growth with Mycemulator based on these parameters (and a few others). This allows us to understand the morphological differences between the strains that we created and visualize their growth dynamics. The general idea, and some of the important features of each model, can be seen in the figure below.

An example of a vegetative mycelium from the filamentous fungi Aspergillus niger. Picture was taken by staining cell walls with Calcofluor White and measuring fluorescence at 410 nm. Colours are artificial, original is grayscale. Captured using our novel microscopy protocol.

Modelling design

The development of computer models of mycelial growth has been ongoing since the 1970’s (Prosser & Trinci, 1979). Our inspiration for Mycemulator came from the DTU BioBuilders team from 2018 who made a simple growth model. We wanted to take the general concept of their model and refashion it to capture more of the features of mycelial growth that we actually observe under the microscope.

When observing mycelial growth through a microscope, it becomes apparent that many relevant features cannot be extracted by the naked eye. In other words, to accurately assess parameters such as branching frequencies, branching angles, and curvature of extended hyphae, computational analysis is needed

We found that one of the biggest limitations was actually a lack of accurate morphological parameters, especially for non-wildtype strains. We then explored fungal literature for image analysis tools, but found none that could quite fit our needs for flexibility, automation and ease-of-use. So we made our own! With our Morphologizer, we could link our simulation with parameters from our images, thereby connecting our wetlab and drylab work. We present some of the tools we considered in the table below.

Their approach	Limitations for our use case	Our approach
Characterising entire fungal pellets from brightfield images (Cairns et al., 2019)	Model only looks at mycelia as whole clumps, not individual hyphae in mycelia	Our model is more fine-grained, analyzing behaviour at the hyphae-level in images.
Characterising initial spore germination from brightfield images (Brunk et al., 2018)	Model focuses on early development from spore to single hyphae, not branching.	Our model focuses on branching behaviour, as that is a main driver of morphological differences.
Three-dimensional analysis of fungal pellets from X-ray microtomography images (Barthel et al., 2019)	Replicates are time-intensive, and the necessary equipment is not widespread. Model requires 3D image of pellet.	Our model is based on 2D-images from our protocol, which are cheap and easy to reproduce on most microscopes

Related models for morphological parameter estimation of fungi that we considered using, before deciding to build our own. Many of these models were made for other purposes or species. As fungi are quite diverse, it is hard to create a “one-size-fits-all” solution.

Model scope and assumptions

Before we wrote any lines of code, we wanted to clearly define the behaviours that our models should revolve around. We selected these behaviours based on our conversations with industry stakeholders in our Human Practices work, and what we saw in our early image data. Our thoughts are summarized in the table below.

Name of behaviour	Description of behaviour	Examples of parameters	Justification for including in models
Tip extension	A mycelium grows by extending cells at the tips of its filaments. This is substrate dependent.	Maximum extension rate	A good growth rate is crucial for strains in industrial biotechnology, and substrate availability can vary wildly.
Branching	In a mycelium, branching happens both at and away from hyphal tips (apical and lateral branching, respectively)	Lateral branching rate, branching angle	The frequency of branch formation has a large impact on the overall structure of the mycelium
Chemotropism	Hyphae can curve towards substrate gradients, though their rigidity prevent them from taking sharp turns.	Maximum curvature	Our image data indicated that chemotropic behaviour occurred more often than we suspected.
Specific substrate consumption	Different part of the mycelium consume different amounts of substrate.	Tip substrate consumption, Non-tip substrate consumption	The growth regions consume more substrate than just maintaining existing hyphae.
Substrate diffusion	Substrate will commonly diffuse to areas of lower substrate concentration	Diffusion constant, substrate consumption, biomass	Most media will allow some level of diffusion and substrate concentration in different areas is crucial for curvature and chemotropism.

Mycelium behaviours that our models include, what parameters these behaviours translate into, and why we chose to include these behaviours.

Events outside the model scope

To further narrow down the model structure and limit the complexity, we discussed behaviors that we consciously wanted to exclude from our models. These were excluded due to adding too much complexity compared to increase in accuracy, being irrelevant to the modeled situation, or being difficult to study. Some of the major decisions can be found in the table below.

Behaviour	Description	Justification for NOT including in models
Spore creation and germination	When fungi create spores (“seeds”) they often create large, multicellular fruiting bodies to facilitate their production.	Under industrial conditions, the mycelium itself is the majority of the biomass. Spores are still important in the initial biomass production, however.
Cell death	In case of starvation, hyphae will eventually die and degrade.	Under industrial conditions, the morphology of a mycelium arises much earlier than cell death occurs.
Three-dimensional growth	Mycelia are three-dimensional structures. Our image analysis and growth model work with two-dimensional data.	This is a large simplification, but simplifies image analysis and makes it more scalable as we can work with 2D image data. For the simulation, some calculations, e.g. substrate movement, become complex and resource-intensive in 3D.
Hyphal fusion (anastomosis)	A hyphal tip can fuse to adjacent hyphae, effectively making it a branch.	We only observed this very rarely in our image data, and chose not to prioritize it.

Mycelium behaviours we did not factor into our models, and the reason we chose not to include them.

Structure of Morphologizer

Now that we have dealt with the model goals and assumptions, the actual descriptions are fairly straightforward. We begin with Morphologizer, our model for estimating morphological parameters of our strains from microscopic image data. By modelling the mycelium as a graph, Morphologizer can estimate various morphological parameters of interest, as seen in figure 1. It is important to note that only hyphal tips and branch points are graph nodes, whereas the hyphal parts connecting them are considered graph edges.

The main model object of Morphologizer is the mycelium graph. Hyphal tips and branch points are graph nodes, while the hyphal parts connecting them are graph edges. By performing various calculations on the graph, morphological parameters can be estimated, such as the branching frequency.

How do we actually convert a noisy image of a mycelium to a graph that we can work with? That is described on the Software page. Below you can see an example of such a conversion from image to graph.

Flow chart illustrating how the Morphologizer analyzes an image. The model works on images of whole mycelia, but only a small part of the image is shown for illustration. (a) Morphologizer requires microscope images of mycelia, where the hyphae are sufficiently pigmented (e.g. by a fluorescent marker) to be filtered out from the background with a threshold. The mycelium should be relatively planar for best results. We refer to our novel microscopy protocol. (b) The mycelium graph treats hyphal branches and tips as nodes, while hyphal segments connecting them are edges. Graph edges have the same path as hyphae, but lose their width. For more details on how the graph is made, check out the Software page. Red nodes are branches in the model, orange nodes are tips. (c) The mycelium graph can be used for several purposes. For example, we can mark branches and tips in the original image. (d) Branching angles are determined by first detecting the two hyphal segments that form the stem, then assuming that the third segment is the branch. The angle is then calculated between the slope of the stem, and the slope of the branch. Orange colour indicates the branch according to the model, while red indicates the “stem”. Estimated angles shown in white. (e) Hyphal curvature is often an indication of chemotropic behaviour. We calculate curvature at each point of an edge by picking two points close by, and calculating the radius of the circle passing through these three points. The hyphal curvature is the inverse of the radius. On the image, brighter spots indicate higher curvature. Curvature is not calculated close to tips or branch points. (f) Space filling, a novel parameter we introduce, describes how densely the hyphae pack. It is defined as the total length of hyphal edges divided by the convex hull of the mycelium graph. For more information on our image preprocessing, check out the Software page.

Morphologizer Parameters

Morphologizer estimates four parameters. Branching frequency, branching angle, curvature and space filling percentage. The parameters and their calculation are explained in the section below. Figure 1 explains some calculations to help.

Branching frequency
We measure branching frequency $f_{\text{branch}}$ in branch points per micrometer of hyphae per hour. We calculate this from our mycelium graph simply by dividing the number of branches $N_b$ (nodes with three neighbors) with the sum of hyphal edge lengths $|e|$ for all $N_e$ edges. We also divide by the duration of mycelium growth $t_{\text{exp}}$. \begin{equation} f_{\text{branch}} = \frac{N_b}{\sum_{i=1}^{N_{e}} |e_i| \cdot t_{\text{exp}}} \end{equation} The branching frequency is later scaled to fit with the segment length of the growth model.

Branching angle
We measure the branching angle in degrees. For a given branch point, we pick the $n$ closest points from each of the three connecting hyphae. These $3n$ points can be thought of as the hyphal structure right around the branch. We then pick the set of two hyphae that minimize the sum-of-squares error, effectively picking the two hyphae that form the straightest line, deeming these two to be the stem from which the branch is created.

Using least-squares regression, we fit a line $l_1$ through the $2n$ points of the stem, then another line $l_2$ through the $n$ points of the branch. Again, not all points from the edges are used, just the $n$ closest ones for each edge. We extract the slopes $m_1$ and $m_2$ for line $l_1$ and $l_2$ respectively, and calculate the branching angle $\theta_b$ with the equation for angles between linear polynomials: \begin{equation} \theta_b = \tan^{-1}(m_1) - \tan^{-1}(m_2) \end{equation} If the angle is negative we reverse the sign as that is not relevant. Note that this method works best for lateral branches, as apical branching does not result in a straight stem. We picked this method because apical branching is relatively rare, occurring only approximately 25% as often (Du et al, 2016)

Curvature
We measure the curvature $\kappa$ in $\mu m^{-1}$. Curvature is estimated at each point $p_i$ of each edge $e$ of the graph, except near branch nodes due to the risk of sharp turns. First we select two points on the same edge as $p_i$, each $k$ points away from $p$ on each side. With the three points $p_{i-k}$, $p_i$ and $p_{i+k}$, we can perfectly calculate the center $x,y$ and radius $R$ of the circle that fits through them. The curvature for a point $i$ on edge $e$ is then simply \begin{equation} \kappa_{e,i} = \frac{1}{R} \qquad \text{where} \qquad R = \text{radius of circle made from points } p_{i-k}, p_i, p_{i+k} \end{equation} If the points are on a straight line, the curvature is undefined. Another way of estimating curvature for $p$ would be to fit a high degree polynomial to $k$ points on each side, then using the derivatives to calculate it. This is more resource intensive, however, and we have a lot of points in our mycelium edges.

Space filling
To differentiate between strains that explore their far surroundings before filling up space close to their starting point, we introduced the space filling parameter $S$. This is defined as the sum of hyphal edge lengths $|e|$ for all $N_e$ edges, divided by the convex hull of the graph $\text{Hull}(g)$. The convex hull is the “area covered” by the fungi, so to speak, so a higher ratio means that the mycelium is filling out more of the area it covers before expanding. \begin{equation} S = \frac{\sum_{i=1}^{N_{e}} |e_i|}{\text{Hull}(g)} \end{equation}

Structure of Mycemulator

Mycemulator is our model for simulating growth of mycelia in conditions that resemble a fermentation. The aim of Mycemulator is to show the differences between the morphological strains that we established during this project, visualize their growth dynamics over time, and do scale-up predictions of different strains. By analyzing images of the different morphological strains grown in the lab and feeding both the analyzed parameters and experimentally estimated growth parameters to the simulation, we take advantage of all information gained through our project and combine it in Mycemulator. A figure illustrating the basic concepts of Mycemulator can be found in the figure below.

Mycemulator simulates the mycelium by dividing it into small linear segments. Many of the model parameters, e.g. the distribution of branching angles, are experimentally determined by us.

Model steps

Initializing the algorithm
Mycemulator is a stochastic model with discrete time steps, similar to random walk. Each simulation round is intended to approximate 1 minute of hyphal growth and each septated hyphal element is given length 1. The simulation is initiated by hyphal elements growing from the center in several directions (see figure), emulating spore germination.

Hyphal tip extension
In any given simulation round every hyphal tip has a probability of being extended based on the local substrate concentration $S$, the growth rate of the given strain $\mu_{max}$ and the half velocity constant $K_s$. For computing the tip extension rate, we applied the Monod equation: \begin{equation} \text{extension rate} = \mu_{max} \cdot \frac {S} {S+K_s} \end{equation} The maximum growth rates for our project were obtained with plate-scale fermentations of the different strains.

Branching
As mentioned earlier, hyphal elements may experience both apical and lateral branching. Apical branching happens when a tip is split into two tips, whereas lateral branching happens from the middle of a hyphal element where the new hyphae branches off the side. For each simulation round, all hyphal tips thus have probability $p_{apical}$ of branching apically and all non-tip hyphal elements have probability $p_{lateral}$ of branching laterally. The overall branching frequency is estimated by Morphologizer and is then converted into a probability of apical and lateral branching, respectively, based on an estimated ratio of the frequency of the two made by Du and colleagues (Du et al., 2018). Hyphal branching events consume substrate and are thus dependent on substrate concentration. This has been implemented in two ways. First, a minimum substrate concentration for a hyphal element to have a probability of branching. Second, the branching probability of a given hyphal element is scaled dependent on the local substrate concentration.

Position calculations

To increase efficiency of the simulation, we store as little information as possible in each simulation round. Therefore, only four parameters were stored: x_mid, y_mid, angle, and tip status. These represent the Cartesian coordinates of the middle of the given hyphal element, the angle with respect to the first axis, and whether the given septum is a tip or not, respectively. This structure means that the simulation runs relatively fast but requires some additional calculations in order to find hyphal element ends and angles.

To determine the midpoint of element $n+1$ ($P_{n+1}$) based on the midpoint of element $n$ ($P_n$) growing at angle $\alpha$, when the change in angle is $\Delta \alpha$: \begin{equation} \overrightarrow{P_nP_{n+1}} = \frac{1}{2} \cdot \begin{pmatrix} \text{cos}(\alpha_n)+\text{cos}(\alpha_n + \Delta \alpha)\\ \text{sin}(\alpha_n)+\text{sin}(\alpha_n + \Delta \alpha) \end{pmatrix} \end{equation} To determine the beginning ($P_n - \frac{1}{2}$) and end ($P_n+\frac{1}{2}$) point of each hyphal element for plotting: \begin{equation} \overrightarrow{P_nP_{n+\frac{1}{2}}} = \frac{1}{2} \cdot \begin{pmatrix} cos(\alpha_n)\\ sin(\alpha_n) \end{pmatrix} \end{equation}
Chemotropism
The direction at which each tip is extended is determined by the direction of the previous element and a small stochastic change. The change is determined by sampling from a beta distribution of parameters which were derived from the image analysis tool. Furthermore, the angle at which a given hyphal element arcs strongly depends on the substrate gradient as the hyphae will seek towards areas of higher substrate concentration. This phenomenon is called tropism and is implemented in Mycemulator by calculating different possible changes in angle and choosing the one which leads the new hyphal tip to an area of higher substrate concentration.

Substrate concentration
At the outset, we modeled substrate being depleted evenly across the simulation area, not accounting for diffusion of substrate. This was extended to use a diffusion model to create a substrate gradient. For both cases, the initial substrate concentration $S_0$ was uniformly distributed across the simulation area and the rate of substrate depletion was calculated based on biomass with two different rates for tip and non-tip hyphal elements. For the uniform field model, substrate concentration $S$ was a function of time $t$ given by: \begin{equation} \frac{\partial S}{\partial t} = \frac {- c_{tip} \cdot b_{tip} - c_{nontip} + b_{nontip}} {r} \end{equation} Where $c$ is rate of substrate consumption per unit biomass, $b$ is biomass, and $r$ is number of compartments in the field.

For the gradient field model, the mycelium was assumed to take up substrate only at the center of the simulation area. For this version, the substrate concentration $S$ was a function of time $t$ and distance from the centre $r$ and was modeled using the numerical solution to a diffusion equation with a point sink: \begin{equation} \frac{\partial S}{\partial t} = D \cdot \frac{\partial^2 S}{\partial r^2}-c_{tip}\cdot b_{tip} - c_{nontip} \cdot b_{nontip} \end{equation} $D$ represents the diffusion constant, $c$ the rate of substrate consumption per unit biomass and $b$ the biomass. For most purposes, the gradient model will be more appropriate as most media will allow some level of diffusion. To be able to model both a mycelium growing under continuous submerged fermentation and solid-state growth on a plate, the model allows the choice of running the simulation with a source at the edge of the field, i.e. a boundary condition $S(t,r_{max})=k$, or without a source.

Uniform substrate

Simulation of ATCC 1015 using uniformly depleting substrate levels and no substrate source.

Substrate gradient

Simulation of ATCC 1015 with substrate diffusion and no substrate source.

Substrate gradient with source

Simulation of ATCC 1015 with substrate diffusion and infinite substrate source at the edge.

Mycemulator parameters

Below you can find an overview of the parameters used in the simulation, the default values of these and how the values were estimated. All parameters in the table below can be tuned by the user to fit the modelling purpose. The default values are based on ATCC 1015.

Parameter ID	Description	Default value	Source
n_septa	Number of initial growing hyphae	4	Based on simulation purpose
S0	Initial overall substrate concentration	100	Estimated
time_step	Time per simulation round (hours)	1/60	Estimated based on growth of ATCC 1015
branch_angle_beta_params	Parameters of a beta distribution fitted to	(4.288, 0.983)	Output from image analysis of microscopic images (default: ATCC 1015)
curvature_gamma_params	Parameters for a gamma distribution fitted to the curvature	(1.611, 2.023)	Output from image analysis of microscopic images (default: ATCC 1015)
q	Branching frequency ()	0.0169	Output from image analysis of microscopic images (default: ATCC 1015)
p_lateral	Probability of lateral branching	0.775*q	A lattice-based system for modeling fungal mycelial growth in complex environments
p_apical	Probability of apical branching	0.225*q	A lattice-based system for modeling fungal mycelial growth in complex environments
branching_substrate_min	Minimum substrate concentration for lateral branching to occur	10	Estimated based on ATCC 1015
mu_max	Maximum growth level	0.275	Experimentally measured (biolector)
extension_substrate_min	Minimum substrate concentration for tip extension to occur	5	Estimated based on ATCC 1015
D	Diffusion constant	0.99	Estimated based on ATCC 1015
S_tip	Substrate use of a single tip in a single round of simulation	10	Estimated based on ATCC 1015
S_nontip	Substrate use of a single non-tip septum in a single simulation round	1	Estimated based on ATCC 1015

Results - Morphologizer

Below, we show the results of Morphologizer analyzing images of our eight strains. A few interesting patterns are worth pointing. While the distribution of branching angles is very stable across all strains, we see a drastic change in the branching frequency and space fillingness of some strains, in particular in our racA knockout and our double knockouts.

Plots of the four parameters estimated by the Morphologizer across eight of our strains. aplD did not grow well, and we were unable to capture images of it to analyze. A point on the first two plots indicates a single image of a single mycelium. Violin plot explanation: the central, thick black line is the two middle quartiles. The white dot is the mean. The coloured area is an estimated density function of the data, created by kernel density estimation with bandwidth 0.5.

Results - Mycemulator

Here, we demonstrate our final simulation of ATCC 1015 based on parameters obtained from both experimental measurements and microscopic image analysis. The substrate gradient is shown in green, i.e. the greener the higher substrate concentration and the mycelium is shown in purple with newer hyphal elements in lighter colors.


The simulation was run for all of the strains that we were able to do good microscopy of. The simulations of the individual mutants can be seen in the Results section.

Future directions

In the future, we would like to work on a number of improvements to our models. For Morhologizer, we would like to add a capacity analyzing 3D images and images of more densely grown mycelia. For Mycemulator, we want to further differentiate substrate use in different regions of the simulation area to fine tune the substrate distribution. Further, we would like to expand Mycemulator to model several mycelia and their interactions, possibly expand into 3D, and incorporate fluid dynamic effects, but those would be rather distant future perspectives

References