Sok Kin Cheng

Sok Kin Cheng

he/him

I am interested in mathematical modeling and data science.

Modeling Nature’s Recyclers: Fungal Decomposition Dynamics

11 minute read

Published:

This is a blog post about the 2021 MCM Problem A on modeling fungal decomposition of woody fibers.

Problem Background

Fungi play a crucial role in the carbon cycle by decomposing plant material and woody fibers. Recent research has identified key fungal traits that determine decomposition rates, particularly growth rate and moisture tolerance. Our goal is to model the decomposition of woody fibers in the presence of multiple fungal species with varying traits.

Requirements

  1. Build a mathematical model that describes the breakdown of ground litter and woody fibers through fungal activity in the presence of multiple species of fungi.
  2. In your model, incorporate the interactions between different species of fungi, which have different growth rates and different moisture tolerances as shown in Figures 1 and 2.
  3. Provide an analysis of the model and describe the interactions between the different types of fungi. The dynamics of the interactions should be characterized and described including both short- and long-term trends. Your analysis should examine the sensitivity to rapid fluctuations in the environment, and you should determine the overall impact of changing atmospheric trends to assess the impact of variation of local weather patterns.
  4. Include predictions about the relative advantages and disadvantages for each species and combinations of species likely to persist, and do so for different environments including arid, semi-arid, temperate, arboreal, and tropical rain forests.
  5. Describe how the diversity of fungal communities of a system impacts the overall efficiency of a system with respect to the breakdown of ground litter. Predict the importance and role of biodiversity in the presence of different degrees of variability in the local environment.

Analysis

To solve this problem, we developed a mathematical model that incorporates:

  1. Fungal growth and interaction
  2. Decomposition of woody fibers
  3. Environmental factors (moisture levels)

Let’s break down the model into its key components:

Fungal Growth Model

We model the growth of each fungal species using a logistic differential equation:

\[\frac{dN_i}{dt} = r_i N_i \left(1 - \frac{N_i}{K_i}\right) - \sum_{j \neq i} \beta_{ij} N_i N_j\]

Where:

  • $N_i$ is the population of species $i$.
  • $r_i$ is the intrinsic growth rate of species $i$.
  • $K_i$ is the carrying capacity of species $i$.
  • $\beta_{ij}$ is the interaction coefficient between species $i$ and $j$.
  • The first term represents logistic growth, and the second term represents interspecies interactions.

In implementation of the model, we define the growth_rate function to represent each fungal species:

def growth_rate(self, species, moisture):
    return species.growth_rate * np.exp(-0.5 * ((moisture - species.moisture_tolerance) / 0.1)**2)

dBdt = [self.growth_rate(species, moisture[moisture_index]) * biomass[i] * (1 - sum(biomass) / species.carrying_capacity)
        for i, species in enumerate(self.fungal_species)]

This equation captures:

  • Species-specific growth rates
  • Moisture dependence
  • Competition between species

Decomposition Model

The decomposition of woody fibers is modeled as:

\[\frac{dL}{dt} = -\sum_i \lambda_i \cdot \text{biomass}_i \cdot L\]

Where:

  • $L$ is the amount of woody fiber litter.
  • $\lambda_i$ is the decomposition rate of species $i$.
  • $\text{biomass}_i$ is the biomass of species $i$.
  • The equation captures the combined decomposition rate of all fungal species.

Where:

In the model implementation, we define the decomposition_rate function to calculate the decomposition rate for each species:

dLdt = -sum(self.decomposition_rate(species, moisture[moisture_index]) * biomass[i] for i, species in enumerate(self.fungal_species)) * litter

This equation relates the decomposition rate to:

  • Fungal biomass
  • Species-specific decomposition rates
  • Moisture levels

Environmental Fluctuations

We simulate different environments (arid, semi-arid, temperate, tropical) by varying moisture levels:

def generate_moisture_series(env_type, num_days, fluctuations=False):
    if env_type == 'arid':
        moisture_mean, moisture_std = 0.2, 0.1
    elif env_type == 'semi-arid':
        moisture_mean, moisture_std = 0.4, 0.15
    elif env_type == 'temperate':
        moisture_mean, moisture_std = 0.6, 0.2
    elif env_type == 'tropical':
        moisture_mean, moisture_std = 0.8, 0.1
    else:
        raise ValueError('Invalid environment type')
    
    if fluctuations:
        moisture = np.random.normal(moisture_mean, moisture_std, num_days)
    else:
        moisture = np.full(num_days, moisture_mean)
    
    return np.clip(moisture, 0, 1)

Solution

Our model reveals several key insights:

  • Fast-growing fungi dominate in stable, resource-rich environments but are more vulnerable to fluctuations.
  • Moisture-tolerant fungi show greater resilience to changing conditions.
  • Fungal diversity generally leads to more stable and efficient decomposition, especially in variable environments.
  • The optimal fungal community composition varies significantly between ecosystem types.

The full solution implementation is revealed in the following code snippet:

# Import necessary libraries
import numpy as np
from scipy.integrate import odeint
from scipy.optimize import minimize
import matplotlib.pyplot as plt
import itertools

Then we define a class called “FungalSpecies” to represent each fungal species:

class FungalSpecies:
    def __init__(self, name, growth_rate, moisture_tolerance, decomposition_rate, carrying_capacity, interaction_coefficients):
        self.name = name
        self.growth_rate = growth_rate
        self.moisture_tolerance = moisture_tolerance
        self.decomposition_rate = decomposition_rate
        self.carrying_capacity = carrying_capacity
        self.interaction_coefficients = interaction_coefficients

In the class “FungalSpecies”, we may see that there are several attributes that define each species:

  • Growth rate: how quickly the species grows under optimal conditions.
  • Moisture tolerance: the range of moisture levels the species can thrive in.
  • Decomposition rate: how quickly the species decomposes woody fibers.
  • Carrying capacity: the maximum population size the species can sustain.
  • Interaction coefficients: the strength of interactions with other species.

We then define a class called “LitterDecompositionModel” to represent the overall model:

class LitterDecompositionModel:
    def __init__(self, fungal_species, initial_biomass, initial_litter):
        self.fungal_species = fungal_species
        self.initial_biomass = initial_biomass
        self.initial_litter = initial_litter

    def growth_rate(self, species, moisture):
        return species.growth_rate * np.exp(-0.15 * ((moisture - species.moisture_tolerance) / 0.1)**2)

    def decomposition_rate(self, species, moisture):
        return species.decomposition_rate * np.exp(-0.15 * ((moisture - species.moisture_tolerance) / 0.1)**2)

    def model_equations(self, y, t, moisture):
        biomass = y[:-1]
        litter = y[-1]

        moisture_index = min(int(t), len(moisture) - 1)

        dBdt = []
        for i, species in enumerate(self.fungal_species):
            growth_term = self.growth_rate(species, moisture[moisture_index]) * biomass[i]
            interaction_term = sum(species.interaction_coefficients[j] * biomass[j] for j in range(len(self.fungal_species)))
            dBdt_i = growth_term * (1 - interaction_term / species.carrying_capacity)
            dBdt.append(dBdt_i)

        dLdt = -sum(self.decomposition_rate(species, moisture[moisture_index]) * biomass[i] for i, species in enumerate(self.fungal_species)) * litter

        return dBdt + [dLdt]

    def simulate(self, time_steps, moisture_series):
        initial_conditions = self.initial_biomass + [self.initial_litter]
        solution = odeint(self.model_equations, initial_conditions, time_steps, args=(moisture_series,))
        return solution

    def optimize_parameters(self, target_biomass, time_steps, moisture_series):
        def objective(params):
            for i, species in enumerate(self.fungal_species):
                species.growth_rate, species.decomposition_rate = params[i*2], params[i*2+1]

            solution = self.simulate(time_steps, moisture_series)
            final_biomass = solution[-1, :-1]
            error = np.sum((final_biomass - target_biomass)**2)
            return error

        initial_params = [species.growth_rate for species in self.fungal_species] + [species.decomposition_rate for species in self.fungal_species]
        bounds = [(0.01, 0.2)] * len(self.fungal_species) + [(0.005, 0.05)] * len(self.fungal_species)

        result = minimize(objective, initial_params, method='L-BFGS-B', bounds=bounds)
        optimized_params = result.x

        for i, species in enumerate(self.fungal_species):
            species.growth_rate, species.decomposition_rate = optimized_params[i*2], optimized_params[i*2+1]

    def plot_results(self, time_steps, solution):
        fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(8, 8))
        
        for i, species in enumerate(self.fungal_species):
            ax1.plot(time_steps, solution[:, i], label=species.name)
        ax1.set_xlabel('Time (days)')
        ax1.set_ylabel('Fungal Biomass')
        ax1.legend()
        
        ax2.plot(time_steps, solution[:, -1], label='Litter')
        ax2.set_xlabel('Time (days)')
        ax2.set_ylabel('Litter Mass')
        ax2.legend()
        
        plt.tight_layout()
        plt.savefig('results.png')
        plt.show()
    
def sensitivity_analysis(model, parameter, range_, time_steps, moisture_series):
        original_value = getattr(model.fungal_species[0], parameter)
        decomposition_rates = []

        for value in range_:
            for species in model.fungal_species:
                setattr(species, parameter, value)
            solution = model.simulate(time_steps, moisture_series)
            decomposition_rate = (model.initial_litter - solution[-1, -1]) / len(time_steps)
            decomposition_rates.append(decomposition_rate)

        setattr(model.fungal_species[0], parameter, original_value)

        plt.figure(figsize=(8, 6))
        plt.plot(range_, decomposition_rates, 'o-')
        plt.xlabel(f'{parameter}')
        plt.ylabel('Decomposition Rate')
        plt.title(f'Sensitivity Analysis: {parameter}')
        plt.grid(True)
        plt.savefig(f'{parameter} Sensitivity Analysis.png')
        plt.show()


def generate_moisture_series(env_type, num_days, fluctuations=False):
    if env_type == 'arid':
        moisture_mean, moisture_std = 0.2, 0.1
    elif env_type == 'semi-arid':
        moisture_mean, moisture_std = 0.4, 0.15
    elif env_type == 'temperate':
        moisture_mean, moisture_std = 0.6, 0.2
    elif env_type == 'tropical':
        moisture_mean, moisture_std = 0.8, 0.1
    else:
        raise ValueError('Invalid environment type')
    
    if fluctuations:
        moisture = np.random.normal(moisture_mean, moisture_std, num_days)
    else:
        moisture = np.full(num_days, moisture_mean)
    
    return np.clip(moisture, 0, 1)

def simulate_environments(model, time_steps, num_days):
    env_types = ['arid', 'semi-arid', 'temperate', 'tropical']
    
    fig, axs = plt.subplots(2, 2, figsize=(12, 8))
    axs = axs.ravel()
    
    for i, env_type in enumerate(env_types):
        moisture_series = generate_moisture_series(env_type, num_days)
        solution = model.simulate(time_steps, moisture_series)
        
        for j, species in enumerate(model.fungal_species):
            axs[i].plot(time_steps, solution[:, j], label=species.name)
        axs[i].set_xlabel('Time (days)')
        axs[i].set_ylabel('Fungal Biomass')
        axs[i].set_title(env_type.capitalize())
        axs[i].legend()
    plt.tight_layout()
    plt.savefig('environments.png')
    plt.show()

def simulate_diversity(model, time_steps, num_days, env_type):
    num_species_range = range(1, len(model.fungal_species) + 1)
    species_names = [species.name for species in model.fungal_species]
    decomposition_rates = []
    
    for num_species in num_species_range:
        species_combinations = list(itertools.combinations(model.fungal_species, num_species))
        
        for combination in species_combinations:
            model.fungal_species = list(combination)
            model.initial_biomass = [0.1] * num_species  # Update initial biomass based on the number of species
            moisture_series = generate_moisture_series(env_type, num_days)
            solution = model.simulate(time_steps, moisture_series)
            decomposition_rate = (model.initial_litter - solution[-1, -1]) / num_days
            decomposition_rates.append(decomposition_rate)
    
    plt.figure(figsize=(8, 6))
    
    # Create a list of labels for the legend
    labels = [f'{name} (The species {i+1})' for i, name in enumerate(species_names)]
    
    # Plot the data points with different colors and labels
    for num, rate in zip(num_species_range, decomposition_rates):
        plt.plot(num, rate, 'o', label=labels[num-1])
    
    plt.xlabel('Number of Fungal Species')
    plt.ylabel('Decomposition Rate')
    plt.title(f'Impact of Fungal Diversity on Decomposition Rate ({env_type.capitalize()})')
    
    # Add a legend to the figure
    plt.legend()
    plt.savefig('diversity.png')
    plt.show()

We thereafter define the fungal species and initial conditions:

# Create fungal species instances
trichoderma = FungalSpecies("Trichoderma harzianum", 0.1, 0.6, 0.02, 100, [1, -0.2, -0.1])
aspergillus = FungalSpecies("Aspergillus niger", 0.05, 0.4, 0.01, 80, [-0.3, 1, -0.4])
penicillium = FungalSpecies("Penicillium chrysogenum", 0.08, 0.8, 0.015, 120, [-0.1, -0.2, 1])

The three fungal species are defined with their respective attributes:

  • Trihoderma harzianum: fast growth, moderate moisture tolerance, high decomposition rate, high carrying capacity, and interactions with other species.
  • Aspergillus niger: slower growth, lower moisture tolerance, lower decomposition rate, lower carrying capacity, and interactions with other species.
  • Penicillium chrysogenum: moderate growth, high moisture tolerance, moderate decomposition rate, high carrying capacity, and interactions with other species.

Next, we create the litter decomposition model and simulate the system:

# Create the litter decomposition model
model = LitterDecompositionModel([trichoderma, aspergillus, penicillium], [0.1, 0.1, 0.1], 100)

# Simulate the model for different environments
num_days = 365
num_time_steps = 365
time_steps = np.linspace(0, num_days, num_time_steps)
simulate_environments(model, time_steps, num_days)

The output of the simulation is a set of plots showing the growth of each fungal species in different environments.

Fungal Biomass in Different Environments

From the four plots, we observe how the growth of each fungal species varies in different environments. The fast-growing species, Trichoderma harzianum, dominates in the semi-arid environment, while the slower-growing species, Aspergillus niger, thrives in the arid environment. Penicillium chrysogenum shows a more stable growth pattern across temperate and tropical environments.

We then create a sensitivity analysis to investigate the impact of growth rates on decomposition rates:

parameter_name = ['growth_rate', 'decomposition_rate']
range_ = np.linspace(0.05, 0.15, 100)
for parameter in parameter_name:
    sensitivity_analysis(model, parameter, range_, time_steps, generate_moisture_series('temperate', num_days))

The output of the sensitivity analysis is a set of plots showing the relationship between growth rates and decomposition rates for each species.

Sensitivity Analysis of Growth Rate

The sensitivity analysis reveals that growth rates have a significant impact on decomposition rates for some special numbers of growth rates and decomposition rates.

Finally, we simulate the impact of fungal diversity on decomposition rates in a temperate environment:

# Simulate the impact of fungal diversity on decomposition rate
env_type = 'temperate'
simulate_diversity(model, time_steps, num_days, env_type)

The output of the diversity simulation is a plot showing how the decomposition rate changes with the number of fungal species in a temperate environment.

Impact of Fungal Diversity on Decomposition Rate

The plot illustrates how the decomposition rate varies with the number of fungal species present in the system. The results suggest that higher fungal diversity generally leads to increased decomposition rates, highlighting the importance of biodiversity in litter decomposition.

Discussion

The model we developed provides valuable insights into the dynamics of fungal decomposition in different environments. By incorporating fungal growth, decomposition, and environmental factors, we were able to analyze the interactions between fungal species and predict their relative advantages and disadvantages in various ecosystems.

Our analysis revealed that fast-growing fungi dominate in stable environments, while moisture-tolerant species are more resilient to fluctuations. Fungal diversity was shown to enhance decomposition efficiency, especially in variable environments, emphasizing the importance of biodiversity in ecosystem functioning.

Conclusion

Our mathematical model provides valuable insights into fungal decomposition dynamics:

  • It highlights the importance of considering multiple fungal traits in ecosystem models.
  • It demonstrates how environmental variability impacts fungal community composition and decomposition rates.
  • It underscores the value of fungal biodiversity in maintaining resilient decomposition processes.

These findings have important implications for ecosystem management and carbon cycle modeling in the face of climate change.

Comments

You May Also Enjoy

Climbing Invisible Mountains: Heuristics on Fitness Landscapes

4 minute read

Published:

In the last post, we embraced “good enough.” This sequel goes deeper into the math of how heuristics move across rugged fitness landscapes—the invisible mountains where every coordinate encodes a solution and its height is the objective value. If you want a warm-up, skim the intro first:

Cooperative Energy Trading in Smart Microgrids: Integrated Differential Equations, Game Theory, and Multi-Agent Optimization for Renewable Energy Markets

35 minute read

Published:

As renewable energy sources become increasingly integrated into power systems, a fundamental challenge emerges: how can communities design fair and efficient energy trading mechanisms that maximize renewable energy usage while ensuring all participants have incentive to participate, maintain grid stability, and operate profitably? This question becomes particularly critical as residential solar installations proliferate, creating distributed networks of energy producers and consumers who must make real-time decisions about energy trading. Unlike traditional centralized power grids, smart microgrids require sophisticated mathematical frameworks that integrate physical energy dynamics, strategic agent behavior, pricing mechanisms, and system stability—all operating under uncertainty from weather patterns, equipment failures, and market fluctuations.

Solar Tracker Cost-Effectiveness: Integrated Economic Modeling and Solar Geometry for Optimal PV System Design

42 minute read

Published:

As solar photovoltaic (PV) systems become increasingly prevalent worldwide, a fundamental question emerges: are solar panel trackers—electromechanical devices that rotate panels to follow the sun—cost-effective investments, or are they merely expensive “toys” that fail to justify their additional cost? This question becomes particularly critical as installation costs, electricity prices, and system configurations vary dramatically across geographic locations, installation types (residential, commercial, industrial), and market conditions. A comprehensive mathematical framework is needed to evaluate whether the energy gains from tracking justify the additional capital costs, operational expenses, and maintenance requirements.