Vol 17, No 1 (2026)

Introduction. Biokinetic models are used to elucidate the characteristics of chemical and physical processes occurring in living organisms. In various fields of biology, including ecology, formal methods for constructing biokinetic mathematical models are well-developed. The most common models comprise systems of differential equations. Two types of problems can be formulated for such systems: direct and inverse problems. The direct problem involves finding a solution to the system given specific parameter values, initial conditions, or boundary conditions. The inverse problem most often consists of model parameter identification: finding numerical parameter values for which the system's solution best fits the available experimental (observed) data.

Typically, realistic and practically relevant systems of equations (mathematical models) lack analytical solutions. While numerical methods for solving direct problems are well-developed, universal and effective methods for solving inverse problems still seemingly do not exist. The aim of this lecture is, through concrete examples, to introduce students to some frequently used methods for the parametric identification of biokinetic mathematical models.

Problem Statement. A typical biokinetic problem is the parametric identification for a system of several interacting components whose concentrations change over time. For example, a system of microorganisms consuming a substrate can be considered. As a result of the interaction between microorganisms and the substrate, microbial concentration will increase over time, while substrate concentration will decrease. If component concentrations at various time points are known, to determining unknown parameters one can minimize a function expressing the sum of squared deviations between observed and predicted values. The primary requirement for the sought model parameters is that their variation must have a noticeable effect on this function, called the residual function. These problems are known as well-conditioned and will be considered in this lecture.

There are several types of residual functions. For example, they may incorporate a scaling parameter, also called a weighting factor or weight. Weights determine the importance of including certain deviations between experimental data and calculated values of variables into the corresponding measure of disagreement. The larger the weights, the more important the corresponding deviation (i.e. data point) is considered, which should influence the results of the identification problem. Weights are often set as the inverse of the variance of observed values. Another type of scaling parameter is the concentration of a characteristic component. In this case, the sum of squared relative deviations, rather than simple deviations, is minimized. Besides the least squares criterion, other distance measures between observed data and predicted values can be used, such as the criterion of absolute deviations or the Chebyshev (minimax) criterion.

Methods for Solving Parameter Identification Problems. The main issue of parameter identification lies in the fact that in real-world problems, the residual function can have several local minima. Unfortunately, there are currently no universal and effective methods for finding the global minimum of a function when it has multiple local minima. Consequently, several types of parametric identification problems are distinguished, for which specific computational methods are applied. Two main types are based on whether analytical expressions exist for functions describing the time-dependent properties of the system. If such expressions exist, linearization of input data is applied in the simplest cases, while nonlinear approaches are used in more complex ones.

1. Linearization is the transformation of nonlinear formulas into a linear expression through specific transformations. Models that can be transformed into linear ones are called intrinsically linear models, as opposed to intrinsically nonlinear models. As a first example, an intrinsically linear function representing the simplest mathematical model expressed by a first-order linear differential equation is considered. Logarithmic transformation allows rewriting this equation in linear form. Special weights must be used to account for the specific transformation applied to the initial nonlinear equation. As a second example, a system of equations describing microbial biomass growth during the consumption of 4-chloroaniline is considered. The method of approximate linearizing transformation is applied for parameter identification in this model.

2. Substitution method. Using the aforementioned 4-chloroaniline consumption model as an example, it is shown how the substitution method allows obtaining a fully satisfactory description of the experimental results by the model within a certain concentration range. In this case, just two points from the original data were judiciously selected for calculations. Furthermore, it is demonstrated how the substitution method can be implemented to describe 4-chloroaniline consumption using more realistically justified functions, specifically the Monod equation. In this case, the function describing 4-chloroaniline concentration dependence on time becomes implicit, and three points from the original data are required to identify the model parameters. Notably, obtained parameters result in satisfactory model predictions across the entire range of 4-chloroaniline concentrations. The substitution method can also be implemented in cases when differential equations cannot be solved analytically. In such cases, function extrema, where derivatives are zero by definition, can be used. Substituting experimental data at extremum points into the system equations allows determining at least some model parameters or their ratios. This approach is illustrated using a model of biomass growth in a continuous bioreactor (chemostat), accounting for the lag due to the ribosome synthesis.

3. Elimination of the independent variable. Most often for biokinetic models, the eliminated independent variable is time. It is demonstrated how eliminating time by dividing one differential equation of the system by another allows computing one parameter for the aforementioned model of biomass growth in a continuous bioreactor. Another fairly simple and limited-use method for solving inverse problems is model simplification. In some cases, a complex model can be well approximated by a simple model that allows analytical integration in a form to which a straightforward linearizing transformation can be applied. It is shown how this approach can be applied to the aforementioned model of biomass growth in a continuous bioreactor.

4. Representing the model equations as a Taylor series. This approach is used when dependent variables change almost linearly or exponentially over a certain time interval. It is demonstrated how this approach is applied to the aforementioned model of biomass growth in a continuous bioreactor, yielding two dependencies between model coefficients.

5. Linearization can also be applied to complex systems of nonlinear differential equations. There are no universal algorithms for applying this method; a creative approach is necessary. It is shown how linearization is applied to the aforementioned model of biomass growth in a continuous bioreactor. All the variables are divided by their maximum experimental values, after which new notations for variables and parameters are introduced. The concentration equation is then rewritten to enable partial integration with the new variables. Then linear regression is applied to the left-hand side (normalized concentration) to identify new introduced parameters, while numerical integration using observed concentrations is applied for the right-hand side. As a result, another model parameter is found.

6. Using the two found parameters and three parameter ratios, the remaining parameter of the model of biomass growth in a continuous bioreactor is determined by means of one-dimensional nonlinear minimization of the residual function. The residual minimum is searched using a MATLAB program. The sum of absolute deviations between the observed and predicted values was used as the residual function, with the corresponding average concentration values as weights. It is shown that the obtained parameter is determined with some ambiguity related to the ill-posed nature of the problem. Comparing the results of identifying different parameters shows that parameters can be determined with significant error, reaching several tens of percent. Solving optimization problems in the MATLAB environment is described in detail including various minimization search functions, as well as settings regulating search accuracy, output of results, and other operation parameters.

7. Finally, it is shown how all parameters in the model of biomass growth in a continuous bioreactor could be identified simultaneously from the experimental data. This inverse problem is solved using multidimensional optimization via the fminsearch function in the MATLAB. The search was conducted in the region of non-negative parameter values. Parameter values obtained earlier by other methods were used as initial guess. It is demonstrated that experimental data can be equally well fitted by different sets of model parameters, which vary widely with relative errors of tens and hundreds of percent. This indicates the inherent ill-posed nature of the formulated problem (i.e. of the particular model used to describe the experimental data).

Northern peatlands represent one of the largest terrestrial carbon reservoirs, playing a crucial yet potentially vulnerable role in the global carbon cycle under ongoing climate change. The stability of this vast carbon stock is intrinsically linked to the thermal regime of the peat soil, which controls key biogeochemical processes such as microbial decomposition, methane production, and plant productivity. However, peatlands are not thermally uniform; they are characterized by a pronounced microtopography, typically featuring a mosaic of elevated, drier features (ridges, hummocks) and water-saturated hollows. Understanding the spatiotemporal dynamics of temperature within this micro-landscape is therefore fundamental for accurate prediction of peatland response to warming. While the general influence of microrelief on temperature is recognized, there is a significant lack of detailed, high-frequency, and multi-year datasets that simultaneously capture the thermal behavior of all key microform elements in continental boreal peatlands, particularly in the vast and critically important region of Western Siberia. This study aims to fill this gap by providing a comprehensive, quantitative analysis of the soil and surface temperature regime in a typical ridge-hollow complex.

The research was conducted at the Mukhrino Bog, a large oligotrophic mire in the Middle Taiga zone of Western Siberia. The study site features a classic ridge-hollow complex (RHC) with well-defined shrub-Sphagnum-dominated ridges (1-3 m wide, up to 60 cm high) and water-saturated sedge-Sphagnum hollows. To investigate the temperature regime, a dedicated monitoring system was deployed in July 2020. The setup included an array of automatic temperature sensors (DS18B20) installed along a 27-meter transect crossing a sequence of a northern hollow, a ridge, and a southern hollow. A total of 11 soil temperature profilers were placed on characteristic microfeatures (hummocks, depressions, slopes), each measuring temperature at depths of 0, 2, 5, 10, 15, 20, 40, and 60 cm. One additional deep probe monitored temperatures down to 320 cm at a ridge location. Air temperature was measured at 2 m and 15 cm above the surface. Data were logged hourly from July 2020 to November 2022. Complementary field measurements included precise leveling of the microrelief and, in April 2023, sampling of frozen peat monoliths from both a ridge and a hollow for subsequent laboratory determination of natural moisture content and absolute dry peat density.

The three-year monitoring period captured a significant range of meteorological conditions, including an extremely cold and low snow accumulation winter in 2020/2021 and a milder with high snow accumulation winter in 2021/2022, as well as contrasting summer conditions. The results reveal a persistent and clear spatial pattern: the saturated hollows were consistently warmer than the elevated ridges throughout the annual cycle. However, the underlying physical drivers of this thermal contrast were seasonally distinct.

During the warm season (April-October), the primary mechanism is the difference in the thermophysical properties of the peat. Water-saturated hollow peat has high volumetric heat capacity and thermal conductivity. This leads to efficient absorption, deeper penetration, and slower release of heat. Consequently, diurnal temperature amplitudes are strongly dampened with depth in the hollows. In contrast, the aerated, drier peat of the ridges has a high pore air content, which acts as an effective insulator. This results in extreme surface heating (up to +34.4°C) and cooling, but a very sharp attenuation of these fluctuations with depth. The ridge peat essentially functions as a "reverse thermos," inhibiting heat transfer into the deeper layers. Data from July 2020, the warmest month, quantitatively illustrate this: the mean monthly surface temperature was +21.2°C in the northern hollow versus +19.6°C on the ridge. More importantly, the temperature difference increased with depth, reaching +2.2°C at 60 cm. The southern hollow was consistently 1-2°C warmer than the northern hollow, likely due to higher water saturation.

In winter (November-March), the dominant factor shifts to the snow cover distribution. The microtopography dictates snow accumulation: deeper snowpacks form in the hollows, while wind exposure keeps the ridges relatively snow-free. Snow is an excellent insulator. Therefore, the thick snow layer over the hollows effectively decouples the soil from the extreme cold air temperatures, maintaining surface temperatures close to 0°C. Over the ridges, the thin snow cover provides minimal insulation, leading to intense soil cooling. In December 2020, the mean surface temperature on the ridge was -9.8°C, while in the hollow it was only -0.9°C—a difference of 9°C. This snow-mediated effect directly controls the depth of seasonal frost. During the harsh winter of 2020/2021, the frost depth reached >60 cm on the ridge but only about 30 cm in the hollow. In the following, milder winter, maximum frost depths were ~55 cm and ~12 cm, respectively. The latent heat released during freezing of the water-saturated hollow peat further moderates cooling.

High-resolution data from a 13-day period in June 2022 further elucidated the diurnal dynamics. Under cloudy, rainy periods, temperatures were uniform across the microrelief. With the onset of clear, anticyclonic conditions, strong diurnal contrasts emerged. During the day, the northern hollow heated most strongly (up to +33.5°C), while shaded areas on the ridge remained cooler. At night, the hollows (especially the southern one) cooled more slowly than the ridge, maintaining a higher temperature. These near-surface patterns persisted to a depth 10 cm. At 20 cm depth, diurnal cycles were almost absent in the hollow but remained pronounced (5-7°C amplitude) on the ridge. Below 40 cm, diurnal variations vanished everywhere, revealing the persistent background spatial pattern: warmer hollows and cooler ridges.

Analysis of peat properties confirmed the foundational differences: the hollow peat maintained a very high natural moisture content (94-98%), while the ridge peat showed lower and more variable moisture (85-95%) and slightly higher dry density in the surface layer.

In conclusion, this study demonstrates that the thermal regime of a boreal peatland is governed by a dynamic interplay between microtopography, moisture content, and snow cover, with seasonally switching dominant mechanisms. The water-saturated hollows act as thermally buffered, energy-accumulating elements, while the aerated ridges experience thermal extremes and function as insulators. The spatial pattern of temperature—warmer hollows, cooler ridges—is a robust feature sustained year-round. The quantified relationships and the extensive dataset presented here are essential for improving process-based models of heat and water transfer in peatlands. This, in turn, enhances our ability to forecast the fate of the massive carbon stored in these ecosystems under changing climatic conditions, particularly for the extensive and vulnerable peatlands of Western Siberia. The observed thermal heterogeneity underscores the necessity of representing microtopographic diversity in landscape-scale models to prevent significant biases in forecasting carbon cycle feedbacks.

Information of the publication of the book "Classification of Vegetation of Russia". Eds. N.B. Ermakov, O.V. Morozova, P.V. Krestov, Yu.V. Plugatar. Volume I. Vegetation of polar deserts, tundra, alpine belt, rocks, screes, snow-covered, aquatic and near-water habitats, treeless and sparsely forested swamps / Eds. O.V. Morozova, O.V. Lavrinenko, Yu.A. Semenishchenkov.