Spatiotemporal Modeling of Coniferous Forests Dynamics along the Southern Edge of Their Range in the Central Russian Plain

Abstract

Forests with predominance of Norway spruce (Picea abies (L.) H. Karst.) and Scots pine (Pinus sylvestris L.) within the hemiboreal zone are considered as secondary communities formed under long-term human activity (logging, plowing, fires and silviculture). This study raises the question—how stable is current state of coniferous forests on the southern border of their natural distribution in the center of Eastern Europe using the example of the Moscow region (MR)? The object of the study are spruce and pine forests in different periods of Soviet and post-Soviet history within the Moscow Region (MR). The current proportion of spruce forests is 21.7%, and the proportion of pine forests is 18.5% from total forest area according to our estimates. The direction and rate of forest succession were analyzed based on current composition of populations of the main forest-forming species (spruce, pine, birch, aspen, oak, linden, and ash) based on ground-based research materials collected in 2006–2019. This allowed to develop the dynamic model (DM) of forest communities with the participation of Norway spruce and Scots pine for several decades. Assessment of the spatial distribution of coniferous communities is based on field data and spatial modeling using remote sensing data—Landsat 8 mosaic for 2020. In parallel, a retrospective model (RM) of the spatial-temporal organization of spruce and pine forests for a 30-year period was developed using two Landsat 5 mosaics. For this, nine different algorithms were tested and the best one for this task was found—random forest. Geobotanical relevés were used as a training sample combined with the 2006–2012 mosaic; the obtained spectral signatures were used for modeling based on the 1984–1990 mosaic. Thus, two multi-temporal spatial models of coniferous formations have been developed. Detailed analysis of the structure of spruce and pine forests based on field data made it possible to track trends of successional dynamics for the first time, considering the origin of communities and the ecological conditions of habitats. As a result, ideas about the viability of spruce and pine cenopopulations in different types of communities were formulated, which made possible to develop a dynamic model (DM) of changes in forest communities for future. Comparison of the areas and nature of changes in the spatial structure of coniferous formations made possible to develop the RM. Comparison of two different-time models of succession dynamics (DM and RM) makes possible to correct the main trends in the transformation of coniferous forests of natural and artificial origin under the existing regime of forestry. A set of features was identified that indicates risk factors for coniferous forests in the region. A further decrease of the spruce and pine plantations and increase of the spruce-small-leaved and deciduous formations are expected in the study area. The proportion of pine-spruce forests does not exceed 3% of the area and can be considered as the most vulnerable type of forest.

1. Introduction

The centuries-old history of anthropogenic activity has led to highly fragmented and relatively young forest cover in the central part of eastern Europe, since almost the entire territory was impacted by various types of agriculture, forest felling, and silviculture. Taking into account the human economic activity, there is a reason to interpret the factors of the land cover change and to consider two directions of distinguishing territorial units: (a) according to “stable” natural characteristics [1,2], and (b) as a result of the interaction of natural and socio-economic relations [3,4,5]. It is believed that the changes of the vegetation cover are subordinate to the anthropogenic factor in most regions of Europe at local level (the lower level of the spatial scale). It was found that, for background intact areas, the features of the structure and functioning of vegetation are largely determined by natural factors, including regional and local climate conditions [6,7,8], the composition of tree dominants, and the set of local flora [9,10]. Assessment of the joint impact of the factor’s combination for active land use areas (including the previous types of anthropogenic impact) is necessary for understanding the modern dynamics of forest vegetation at a regional level.

The current composition of temperate forests in Europe, including the center of the Russian Plain, is largely determined by forestry practices and modern European coniferous forests are the result of long-term human activity [11]. Fires in coniferous forests occurred not only from atmospheric electrical discharges, but also due to human activities. Pine is highly resistant to fires. Its thick bark is less damaged by ground fire and its deep taproot makes it more resilient compared to other conifers. Therefore, light ground fires sometimes contribute to pine forests restoration [12]. In this regard, the pine forests in moraine landscapes of the center of the Russian Plain are often either plantations or the result of natural disturbances. Pine forests are sustainable in extreme habitats predominantly [13,14]. The following questions remain unanswered: what the nature of coniferous forests dynamics is, and what are the prospects of their development within the transition zone from coniferous–broad-leaved to broad-leaved forests.

The forest area of the central regions in Russia was about 26% in the early 19th century. The forest cover proportion of the Moscow governorate has almost doubled over the past 100 years, largely due to silviculture in clearings and burned-out areas [15,16]. Additionally a massive process of natural reforestation is actively taking place as abandoned agriculture lands turns into overgrown lands in the central part of the Russian Plain [17]. Forest cover increased by 4.7% in most of the eastern European countries over the period of 1985–2012 according to [18]. However, forest cover decreased in the most populated administrative regions of Russia (Moscow and Leningrad). Forest loss area amounts to more than 26 thousand hectares in the MR during the period from 1985 to 2012 due to annual violations of forest cover. Forest loss during the period 2007–2012 was 25% larger than during 1985–2006.

The degree of various factors that impacted forest cover in the MR changed in different periods of post-Soviet history (1990–2000 and 2000–2020) depending on the socio-economic priorities of land use. Development of cottage construction and general infrastructure has come out on top at the present time. Agriculture has lost its relevance. The complicated history of the forest cover of the Moscow Region (MR) and the current condition of Moscow green belt leaves many questions regarding to the assessment of forest condition and dynamics as well as the possibility of relevés forecasting the ecological and recreational potential of forest stands of natural and artificial origin. The development of dynamic models is needed as field and remote information on the main patterns of forest cover development. The vegetation cover of the MR is researched well enough [19,20,21,22,23,24,25,26,27], including the map of vegetation cover based on forest inventory data and field studies [28]. The active use of mathematical methods and data analysis allows for a fresh look at the issues of forest’s biodiversity assessment and vegetation mapping.

One of the most popular interdisciplinary topics between studies of vegetation cover and remote sensing is the land use and land cover change (LULCC) which offers a set of methods for spatial and temporal analysis of vegetation cover dynamics [29]. Quantification of vegetation cover change is considered as still underactive stream; therefore, some hardened beliefs may become a subject of disproof. Contrary to the prevailing view that forest area has declined globally [30]—tree cover has increased by 2.24 million km² (+7.1% relative to the 1982 level), which is shown based on multiple satellite sensors from 1982 to 2016 [31]. Studies of LULCC have crucial importance in national and regional context to produce summary of vegetation dynamics for the period of most rapid changes of vegetation cover. Particularly, LULCC is used to evaluate and to benchmark results of ambitious national projects for reforestation [32] and flooding and soil erosion mitigation programs [33]. Studies of LULCC in Eastern Europe and especially in post-soviet countries is of special interest because of significant and multidirectional changes in trends of land use [34].

One of the significant phenomena of this period is the launch and maintenance of satellites constellation which continuously surveys the Earth’s surface in multispectral bands. Data are continuously updated as well as mapping methods are being updated [35,36]. Remote sensing data are widely used in best practices to assess the structure and properties of vegetation, ecosystems management, and forestry optimization [37,38,39,40,41]. Landsat 5, 7, and 8 collections of imagery provide an opportunity to trace dynamic changes in vegetation cover in the periods 1982–2012, 1999–2003, and 2013–ongoing, respectively. Spatial resolution of 30 m can be considered optimal for the wide range of medium-scale regional studies. LULCC surveys are very effective when they implement Landsat time-series, spectral-temporal modeling, and modeling based on spectral signatures of land/vegetation cover types [42,43,44]. We use this approach in remote sensing and modeling section of methodology of the current study. Coniferous forest formations within the study area were not previously the subject of quantitative analysis of successional dynamics based on remote sensing time-series according to our literature search.

The development of retrospective models (RM) of forest cover is important for understanding the mechanisms of its recovery dynamics. To solve this problem, information is needed on the main regularities of the successional development of forest ecosystems in a particular region under the influence of multidirectional anthropogenic factors.

The aim of current study is to perform vegetation mapping and to identify coniferous forests dynamics along the southern edge of their range in the central Russian Plain in the zone of large metropolis influence (case study in the MR). This study is based on both field and remote sensing data. The research objectives are as follows: (1) the assessment of the current composition and distribution of coniferous forests in 2020; (2) the development of a dynamic model (DM) of coniferous forest development for several decades; and (3) the development of a RM of spruce and pine forests dynamics for the time slices of 1990 and 2010. The results of current study will contribute to the development of plans for sustainable management and conservation of forest biodiversity under different management scenarios.

2. Materials and Methods

2.1. Study Area

The MR is located in the Central part of East European (Russian) Plain—35°10′–40°15′E, 54°12′–56°55′N, its area is 4.58 million ha (taking into account the territory of “New Moscow” with an area of 0.15 million ha) and its population is more than 20 million people (about 8 million people live in the MR and more than 12 million people live in Moscow) (Figure 1).

According to the map of climatic regions the MR is assigned to moderately continental region [46]. Average annual air temperature is 2.7–3.8 °C, precipitations are 560–640 mm [47]. Relief is gently hilly, heights vary from 90 to 320 m, on average—174 m above sea level and the average slope is 2.06° (0–30.9°).

There are several important landscape and botanical–geographical boundaries in the MR due to climatic gradients and due to diverse history of glaciation during the Quaternary period [48].

In accordance with the geobotanical zoning scheme, the main part of the MR is located within coniferous–broad-leaved forest zone, in the south of the MR there is a border with broad-leaved forests zone [45,49,50] (Figure 1). Zone boundary is associated with the change of geomorphological conditions caused by the boundary of Moscow glaciation: poorly drained loamy plateaus in the northwest are replaced towards the southeast by well drained eroded interfluves.

Primary forests have not survived in the MR by now. Conditionally indigenous communities are widespread in the MR. Under “Conditionally indigenous” we suppose the forests that are close to indigenous analogs in the composition of tree and subordinate layers, but significantly differ from them in age structure. Consequently, the forest cover of the MR is represented by the succession mosaic of forests of different composition, age, and origin. The proportion of coniferous forests in forest cover area is 45% currently [51], of which the spruce forests occupy 25%, the pine forests occupy 20%. The proportion of plantations is about 22% of the Stare Forest Fund of the MR [51]. At the same time plantations do not include forests over 60–80 years according to the rules of the State Forest Inventory. They are renamed to “natural forest stands” after reaching 60 years age. Because coniferous stands are mainly silvicultural, at the same time oak, linden, fir, and Siberian larch are relatively sparse species, it can be concluded that the predominant part of “mature” coniferous forests (60–80%) in the MR are plantations. These data are consistent with the opinion of other experts [52,53].

2.2. The Main Trends of Forest Management and History of Coniferous Forests Formation

The forest cover of the MR has varied greatly over the past centuries in accordance with social and historical conditions in Russia. Coniferous, coniferous–broad-leaved, and broad-leaved forests dominated here in the pre-agricultural period (up to the 10th century) [20,54]. Huge areas of forests began to be reduced for plowing in the 10–13th centuries. The area of arable land approached its maximum possible level by the beginning of the 16th century [55].

In the 18th and 19th centuries the dynamics of forest cover was primarily due to the peculiarities of agriculture. In the 18th century, the three-field system of agriculture already existed in the Central Russia, and the spatial structure of land use was established, including large fields that exist up to the present time [56]. Agricultural practices were mainly extensive and soil nutrient exhausting. They plowed up new rested lands, which had time to be overgrown with forest. Thus, the land rotation was combined with three-field agriculture in the MR at least until the second half of the 19th century [57,58]. Forest and non-forest areas often changed their spatial location. The composition of forests in watersheds was represented mainly by early successive small-leaved species (birch and aspen) aged 30–35 years. The total forest cover area of the MR at the end of the 18th century was 45.5%, and by the middle of the 19th century it had decreased to 38.8% [15].

In the late 19th—early 20th centuries, after the abolition of serfdom, the economic situation of peasants and landowners in the center of European Russia generally worsened [59,60]. Landowners often felled their forests in large quantities for sale if their economy was not generating the required income. The population of peasants increased during this period, but the methods of farming remained old, the yield was extremely low, so the peasants were forced to plow up almost all of their land, even part of the hayfields [59]. As a result, by 1914, the forest cover of the MR decreased and amounted to 26.3%.

The forest cover increased to 32.3% in the postwar years (1956) [61] due to planned management of forestry, mainly by planting coniferous plantations. The forest cover reached 42.7% by 2020 due to the creation of forest plantations and spontaneous overgrowing of abandoned lands [62].

The practice of silviculture started in MR in the first half of 19th century [63]. Forestry technology was developing slowly and the area of plantations created before 1917 was insignificant—14.6 thousand hectares (1.7% of total forest area) [15]. Silviculture became widespread only in the second half of the 20th century due to the need to eliminate the timber shortage which was formed during the war years. Pine silviculture was widespread in the postwar years, in the late 1940s. There was a transition to spruce plantations since the 1960s because the spruce plantations are less damaged by ungulates, diseases and pests [64,65].

Forests of MR are periodically exposed to spontaneous natural influences: fires, invasions of insect pests, unfavorable weather conditions, and forest diseases [66]. The degree of natural impacts increases under inappropriate environmental management practices. Thus, uncleared windfalls triggered an outbreak of reproduction of the bark beetle (Ips typographus) in the 2000s, which caused the mass destruction of coniferous forests. In the 2000–2006 the bark beetle continued to be the main cause of forest stands dieback. The forest pathological situation of the MR forests was relatively sustainable in the 2007–2009 (the area of the dead-wood was 666–1061 ha annually). The area of dead-wood increased sharply to 21.4 thousand hectares in 2010 due to abnormally dry and hot summers and persistent ground fires (20 times compared to 2009). The bark beetle outbreaks in spruce forests became especially typical in recent years [67,68,69] (Figure 2). The territories after cuttings were grown with plantations.

2.3. Methods

2.3.1. Field Data and Classification

The main source of data is the database of geobotanical relevés collected for the 2006–2019 within the MR. The total number of relevés is 1608, of which 906 are in coniferous forests (Table 1). Field data were collected using standard methodology [70] in forest sample plots of 400–625 m² area. The composition and vertical structure of communities were assessed by layers to study the organization of different types of coniferous forests: A—tree layer; B—undergrowth (shrub layer and undergrowth 1–10 m high); C—herb-shrub layer (below 1 m); D—moss-lichen layer.

The classification of communities was done in accordance with ecological-phytocoenotic approach [71]. Syntaxons in the ranks of formation and association group are distinguished based on dominant tree species and representation of ecological and morphological groups in an underlayer. We selected relevés with dominance of spruce or pine to study the composition and the structure of formation with the participation of coniferous stands. We selected relevés where spruce accounted for 75% and more as the spruce forests formation. Respectively, we selected relevés where spruce and aspen/birch stands were almost equal as the spruce–small-leaved formation. Similarly, we selected relevés where pine accounted for 75% and more as the pine formation, and relevés where the proportions of spruce and pine were almost equal as the pine–spruce formation.

Each formation is represented by communities with different combinations of common dominating species in an underlayer. The identification of syntaxa at the level of association group was performed basing on the prevailing ecological and morphological groups of plants in subordinate layers. Each spruce and spruce–small-leaved forest formation is represented by four or six association groups [72].

• dwarf shrubs–small herb–green moss (#1—DshShG);
• small herb (#2—Sh);
• small herb–broad herb (#3—ShBh);
• broad herb (#4—Bh);
• meadow herb (#5—Mh);
• dwarf shrubs–herbal-sphagnum (#6—DshHSh).

The detailed analysis of the composition of spruce and pine formations makes it possible to track the main trends of the natural dynamics of communities, considering their origin and ecological conditions of habitats.

A detailed analysis of the composition of all layers of communities was carried out in terms of the activity (A) of the main tree species within the framework of the identification of association groups [73].

A = √(F × D)

(1)

where F—is the relative abundance of tree species at all relevés; and D—is the average value of the abundance of tree species (%) for relevés where this species was recorded. This makes it possible to assess the prospects of the renewal of pine and spruce communities based on the cenopopulation dynamics of the main species of coniferous forest forming species (Picea abies, and Pinus sylvestris), small-leaved species (Betula pendula, B. pubescens, and Populus tremula) and broad-leaved species (Quercus robur, Tilia cordata, Acer platanoides, and Ulmus spp.).

Indirect ordination methods were applied to interpret the ecology of groups—the nonmetric multidimensional scaling (NMDS ordination) in the R software environment [74], using the square root transformation, the Wisconsin double standardization, and the Bray–Curtis distance [75]. Values of the Ellenberg’s ecological factors (given in scores for each species according to Ellenberg) were calculated for axes during ordination and interpretation [76,77]. The indicators such as light (L), nutrients (N), soil pH (R), and moisture (M) were calculated using software Juice 7.0 [78] for each relevé (weighted by species cover). The factors values were used to position ecological factors vectors in an ordination space and to assess the differences between community groups. Note that some uncertainty is typical for all ecological scales. However, due to the lack of field measurements of environmental variables, we used literature data to interpret ecological differentiation of communities. The differentiation of groups was analyzed basing on the composition of subordinate layers B, C, and D.

2.3.2. Use of Remote Sensing Data

An analysis of the dynamics and the current distribution of coniferous forests was carried out basing on spatial and temporal modeling using the mosaics of historical Landsat 5 images (1990 and 2010) and spectral indices calculated on their basis. Landsat 5 is the most suitable sensor because it has the longest period of observations. It starts in 1984 and finishes in 2012. We collected data in the beginning of this period and at very end. The following workflow was applied:

• Landsat 5 mosaicking for two periods: 2006–2012 (hereinafter 2010) and 1984–1990 (hereinafter 1990);
• Manual spatial rarefication of relevés which subject of forest loss and degradation;
• Testing of different models for best prediction of forest formations for the 2010;
• Creation of spectral signatures for the 2010 (training the classifier);
• Using the classifier on the 2010 mosaic;
• Using the classifier on the 1990 mosaic.

There are several limitations of remote sensing-based modeling and mapping in the Moscow region. The main issue is that there are only 13 clear days without cloud cover during the period May–September (

Table A6. Distribution of clear and cloudy days within the year.

Month	1	2	3	4	5	6	7	8	9	10	11	12	Total
Total cloud cover
Clear days	1	2	4	3	4	2	2	3	2	2	1	1	27
Cloudy days	8	9	13	15	19	20	21	19	16	11	7	7	165
Overcast days	22	17	14	12	8	8	8	9	12	18	22	23	173
Lower cloud cover
Clear days	5	7	10	9	8	5	6	8	8	5	3	3	77
Cloudy days	11	12	13	16	20	22	22	19	16	14	9	11	185
Overcast days	15	9	8	5	3	3	3	4	6	12	18	17	103

). It means that Landsat mosaic without clouds and shadows can be collected during the 4–6 year period. However, we assume that forest succession changes in boreal forest are enough slow to disregard them in this period. Rapid forest disturbances such as cuts, fires, etc., are accounted during manual relevés rarefication.

Mosaics were created using the script for Landsat 5 the Google Earth Engine environment. The script includes preliminary preparation such as date filtering of the images and atmospheric correction for Landsat 5 Surface Reflectance Tier 1 image collection. Then the script performs cloud masking based on different spectral brightness of (i) blue and cirrus bands, (ii) visible bands, and (iii) infrared bands and separates clouds from snow based on green and swirl bands. Finally the script calculates median value between unmasked pixels. [79,80,81]. The script is available in Google Earth Engine repository [82]. Mosaics were collected for the two periods:

• 1984–1990, days of year: 120–260 (May 1–September 15), total 224 images (further referred as the “1990”), and
• 2006–2012, days of year: also 120–260, total of 175 images (further referred as the “2010”).

All spectral bands were used, except for thermal and panchromatic, the spectral indices NDVI, EVI, MSAVI, SAVI, NDMI, and NBR were calculated [83]. Six forest formations were modeled using the training sample: 4 formations of coniferous forests, the deciduous forests and non-forest areas. The deciduous forests include 2 categories—small-leaved and broad-leaved forests [84].

The training samples (points of geobotanical relevés) are combined with spectral reflectance and indices for the 2010 mosaic, and then the classifier was trained. The classifier then was used on 2010 mosaic and 1990 mosaic. The quantity of relevés used in the training sample was reduced during spatial rarefication and amounted to 1542. We performed preliminary visual analysis of forest disturbances and spatially rarefied forest relevés which are in later disturbed areas, and we also consider that slow successional changes of forest formations are insignificant during period 2012–2019.

2.3.3. Specific Classification Algorithms

Seven classification algorithms have been tested: LibSVM, boost, decision tree, normal bayes, random forest, and KNN, shark random forest using OrfeoToolbox software [85].

LibSVM classifier. Support vector machine works by drawing a line between different clusters of points to be grouped into classes. On one side of the line there will be points belonging to one class, on the other side—to another class. The classifier seeks to increase the distance between the lines being drawn and the points on different sides to increase its “confidence” in the class definition. When all points are plotted, the side to which they are projected is the class to which these points belong. The best or optimal hyperplane separating the two classes is the line with the largest difference. Only these points are of importance in defining a hyperplane and in constructing a classifier. These points are called support vectors. To determine the values of the coefficients that maximize the difference, special optimization algorithms are used [80,86].

Boost classifier. Boosting is a family of ensemble algorithms, the essence of which is to create a strong classifier based on several weak ones. To do this, first one model is created, then another model, which tries to correct the errors in the first one. Models are added until the training data is perfectly predicted or until the maximum number of models is exceeded. AdaBoost is used in conjunction with short decision trees. After the first tree is created, its effectiveness is tested on each training object to understand how much attention the next tree should pay to all objects. Data that are difficult to predict are given more weight, and those that are easy to predict are given less weight. The models are created sequentially one after the other, and each of them updates the weights for the next tree. Once all the trees have been built, predictions are made for the new data, and the performance of each tree is determined by how accurate it was on the training data. The method is sensitive to the quality of the data and to the presence of anomalies in them [86].

Decision tree and random forest. This classifier splits data into smaller and smaller subsets based on different criteria, that is, each subset has its own sorting category. With each division, the number of objects of a certain criterion decrease. A decision tree can be represented as a binary tree. Each node represents an input variable and a split point for that variable. Leaf nodes are an output variable that is used for prediction. Predictions are made by traversing the tree to a leaf node and printing the class value at that node. Trees learn quickly and make predictions. In addition, they are accurate for a wide range of tasks and do not require any special data preparation. The classification will come to an end when the network reaches a subset with only one object. [86,87].

Random forests classifier. Random forest is a very popular and efficient machine learning algorithm. This is a kind of ensemble algorithm called bagging. Decision trees are most commonly used to evaluate all statistical models. The training data is split into multiple samples, for each of which a model is created. When a prediction needs to be made, each model does it, and then the predictions are averaged to give a better estimate of the output. In the random forest algorithm, decision trees are built for all samples from the training data. When constructing trees, random features are selected to create each node. Taken together, the resulting models are not fully accurate, but when they are combined, the prediction quality improves significantly. If an algorithm with high variance, such as decision trees, performs well, then this result can often be improved by applying bagging [86,88].

Normal (Naive) Bayes classifier. Such a classifier calculates the probability of an object belonging to a certain class. This probability is calculated from the chance that an event will occur, based on events that have already occurred. Each parameter of the classified object is considered independent of other parameters. The model consists of two types of probabilities, which are calculated using the training data: (i) the probability of each class, and (ii) the conditional probability for each class at every value of x. After calculating the probabilistic model, it can be used to make predictions with new data using Bayes’ theorem. If you have real data, then, assuming a normal distribution, it is not too difficult to calculate these probabilities. Naive Bayes is so called because the algorithm assumes that each input variable is independent. This is a strong assumption and does not match the actual data. Nevertheless, this algorithm is effective for several complex tasks [89].

K nearest neighbor (KNN) is nonparametric algorithm which includes a set of techniques for estimating a regression curve without making strong assumptions about the shape of the true regression function. k-NN is a type of classification where the function is only approximated locally and all computation is deferred until function evaluation [90].

2.3.4. Quality Assessment and Additional Sources of Data

Classification quality is assessed using the confusion matrix and the overall accuracy [91]. The following strategy was used to select test samples. We selected test sample from the initial database of relevés and we did not use them for modeling only for testing. Since the amount of relevés is limited and non-equal between different formations we tried to make the proportion of test sample about 10–15% from total quantity of relevés in each formation. As a result, we obtained test sample proportions 8–19% and amounts 10–100 relevés (Table 1).

The spatial distribution of forest formations in 2020 is taken from our previous study [84]. It uses the data on forest cover for the 2020 [84]. It this study we used multi-seasonal mosaics of Landsat 8 (March, May, June, July, and September), spectral indices for each season, DEM with morphometric variables and HH/HV polarization of Palsar radar image. Field data were same relevés used in the current study. The modeling method was maximum entropy with software Maxent and SDM toolbox. SDM toolbox allows for multiple models calibration, comparison, and best model selection [92], spatial jackknifing [93], and testing feature class combinations and regularization multipliers. These features allow to avoid overfitting together with best model selection.

Additionally, for benchmarking our modeling results we use the Global Forest Watch (GFW) data on forest cover using the layer “percent forest cover” for two periods: 2010 and 2020 by calculation of forest cover in 2000 minus forest loss (cumulative for 2010 and 2020) and plus forest gain (cumulative for 2012). The forest cover mask is produced according to recommendations of Global forest watch team: >30% tree canopy [17]. Water mask is taken from Landsat quality band. Settlements mask is taken from vectorized topographic maps.

3. Results

3.1. Characteristics of Present Spruce and Pine Forests

Four forest formations were identified through relevés classification: Spruce, Spruce–small-leaved, Pine–spruce, and Pine forest formations. The following numbers of association groups were identified within each formation: four in the Spruce, Spruce–small-leaved, and Pine-spruce formations; and six association groups within the Pine formation (Table 1).

Table 1. Distribution of points of relevés by formations and association groups.

Formation	Post Field Classification		Modeling
	Association Group	Number of Points	Number of Points	Train Sample	Test Sample	Proportion %
Spruce	DshShG	37	349	299	50
	Sh	40				14
	ShBh	148
	Bh	148
Spruce–small-leaved	DshShG	32	219	179	40
	Sh	22				18
	ShBh	81
	Bh	102
Pine–spruce	DshShG	32	122	107	10
	Sh	16				8
	ShBh	44
	Bh	42
Pine	DshShG	46	216	176	40
	Sh	23
	ShBh	35				19
	Bh	64
	Mh	15
	DshHSh	45
Deciduous forests			571	471	100	18
Non-forest			65	60	5	8

Table 1. Distribution of points of relevés by formations and association groups.

Formation	Post Field Classification		Modeling
	Association Group	Number of Points	Number of Points	Train Sample	Test Sample	Proportion %
Spruce	DshShG	37	349	299	50
	Sh	40				14
	ShBh	148
	Bh	148
Spruce–small-leaved	DshShG	32	219	179	40
	Sh	22				18
	ShBh	81
	Bh	102
Pine–spruce	DshShG	32	122	107	10
	Sh	16				8
	ShBh	44
	Bh	42
Pine	DshShG	46	216	176	40
	Sh	23
	ShBh	35				19
	Bh	64
	Mh	15
	DshHSh	45
Deciduous forests			571	471	100	18
Non-forest			65	60	5	8

3.1.1. Spruce and Spruce–Small-Leaved Formations

The formation composition of spruce forests is complex (combinations of spruce with birch, aspen, pine, and broad-leaved species) and largely characterizes the composition of undisturbed forests within coniferous–broad-leaved zone (Figure 3A). The area of plantations is high (mainly monodominant spruce stands). The vegetation of subordinate layers characterizes full aspect of the transition of representatives of a wide range of species—from boreal to nemoral composition. Spruce–small-leaved forests are considered here as the successional stage of spruce forests. Their reforestation occurs according to two scenarios with the participation of small-leaved species—birch (Betula pubescens, B. pendula) and aspen (Populus tremula). According to the first scenario, the spruce population is formed during the initial overgrowing of felling or arable land with small-leaved trees in the first 20–30 years or with late regeneration of spruce, gradually entering main layer. According to the second scenario, a mixed composition is formed as the result of birch active renewal in spruce plantations with inadequate quality of care (Figure 3B).

The specified classes are clearly distinct according to the diagnostic criteria of subordinate layer vegetation–indicator species and ecological-coenotic groups of species.

The indicator species of the DshShG group of communities (#1) with high IndVal (IV) values include green mosses (Pleurozium schreberi and Hylocomium splendens), small boreal shrubs such as blueberry and cowberry, and typical boreal small herb species. The ground layer of spruce Sh forests (#2) has only two indicator boreal species: Oxalis acetosella and Mycelis muralis. Spruce ShBh forests (#3) have Corylus avellana in the shrub layer and the nemoral species Ajuga reptans in the herb layer. The large number of diagnostic nemoral species in spruce forests (#4) shows that they are mature spruce communities (

Table A2. Indicator values (IV) for species in community groups of spruce and spruce–small-leaved forests.

Community Group
#1—DshShG		#2—Sh		#3—ShBh		#4—Bh
Species	IV	Species	IV	Species	IV	Species	IV
Pleurozium schreberi	69	Oxalis acetosella	34	Corylus avellana (B2)	28	Lamium galeobdolon	45
Vaccinium myrtillus	60	Mycelis muralis	25	Oxalis acetosella	27	Aegopodium podagraria	43
Hylocomium splendens	52			Ajuga reptans	26	Carex pilosa	43
Frangula alnus (B2)	44					Ranunculus cassubicus	41
Vaccinium vitis-idaea	39					Corylus avellana (B2)	37
Calamagrostis arundinacea	37					Pulmonaria obscura	32
Orthilia secunda	35					Dryopteris filix-mas	31
Luzula pilosa	32					Stellaria holostea	31
Maianthemum bifolium	30					Asarum europaeum	30
Trientalis europaea	30

¹ We show here species with IndVal > 25 and significance value < 0.05.

).

The composition of trees different age range provides the information on the viability of spruce cenopopulations in the different types of communities. Analysis of the participation of other tree species (pine, birch, aspen, oak, linden, and ash) in different vegetation layers, in addition to spruce, indicates the direction of the successional change of communities. In spruce forests, the proportion of spruce in main tree layer (A), undergrowth (layer B), and herb layer (C) prevails in comparison with small-leaved and broad-leaved tree species in all groups of communities (Figure 4A). In spruce–small-leaved forests, both spruce and small-leaved tree species participate in similar proportion in tree layer (A), but in the undergrowth (layer B), spruce is replaced to small-leaved tree species in groups #1–3, except the broad-leaved group of communities (#4). This fact indicates the possibility of developing not spruce, but spruce–broad-leaved or broad-leaved forests in future. In herb-shrub layer (C), the participation of spruce, small-leaved and broad-leaved tree species practically does not differ, which indicates equal starting conditions for the spread of ovules, their germination, and the survival of young plants of different species (Figure 4b).

3.1.2. Pine and Pine-Spruce Formations

The ecological range of pine is wide. Pine forests can be found on soils of different texture and moisture regime. Typology diversity and accordingly, the number of community groups in comparison with spruce forests increases: 1—dwarf shrubs-small herb-green moss (DshShG); 2—small herb (Sh); 3—small herb-broad herb (ShBh); 4—broad herb (Bh); 5—meadow herb (MH); and 6—dwarf shrubs–herbal-sphagnum (DShHSh) (Figure 5).

The composition of indicator species of groups 1–4 of pine forests is quite like groups 1–4 of spruce forests (Table A2 and

Table A3. Indicator values (IV) for species in community groups of pine and pine-spruce forests.

Community Group
#1—DshShG		#2—Sh		#3—ShBh		#4—Bh		#5—Mh		#6—DshHSh
Species	IV	Species	IV	Species	IV	Species	IV	Species	IV	Species	IV
Pleurozium schreberi	76	Mycelis muralis	43	Corylus avellana B	38	Athyrium filix-femina	48	Trifolium medium	67	Eriophorum vaginatum	71
Vaccinium myrtillus	59	Oxalis acetosella	42	Dryopteris carthusiana	36	Ranunculus cassubicus	36	Calamagrostis arundinacea	61	Sphagnum magellanicum	71
Hylocomium splendens	52	Circaea alpina	31	Paris quadrifolia	32	Lamium galeobdolon	32	Agrimonia eupatoria	58	Ledum palustre	57
Dicranum polysetum	50	Sorbus aucuparia B	30	Viburnum opulus	31			Knautia arvensis	58	Vaccinium uliginosum	57
Picea abies C	39							Leucanthemum vulgare	55	Sphagnum angustifolium	57
Melampyrum pratense	32							Veronica officinalis	51	Vaccinium vitis-idaea	55
								Clinopodium vulgare	46	Carex globularis	54
								Carex pallescens	44	Oxycoccus palustris	43
								Vicia cracca	42	Polytrichum strictum	43
								Campanula persicifolia	41	Aulacomnium palustre	41
								Fragaria vesca	40	Betula pubescens B	38
								Pinus sylvestris C	37
								Lathyrus vernus	35
								Melica nutans	34
								Antennaria dioica	33
								Astragalus glycyphyllos	33
								Viola hirta	33
								Chamaenerion angustifolium	32

¹ We show here species with IndVal > 30 and significance value < 0.05.

). The ratio of ecological groups of plants in ground layers in groups 1–4 of pine formation (DshShG, Sh, ShBh, and Bh) is also similar to groups 1–4 of spruce formation, where the proportion of boreal species gradually decreases in groups 1 to 4 from 80 to 20%, and the proportion of nemoral and nitrophilic-wet herbal species, on the contrary, increases in the same pattern. However, in pine forests of the Mh (#5) and dwarf shrubs–herbal-sphagnum group (#6), there is a set of species that indicate unique habitat conditions. In the Mh group of communities (#5), species of the meadow and the forest edge group prevail (50% of the total composition). In the DshHSh group (#6), which is widespread in humid habitats, the species of oligotrophic-bog group (60%) naturally dominate; the herb-bog group is 6%.

Pine-spruce forests are generally plantations with a two-layer stand, where the upper canopy is formed by pine, and the second layer—by spruce. Presumably, we consider them as a transitional successional stage to spruce communities. Considering the participation of pine in the different layers of these communities, there is no pine in undergrowth (layer B) and in herb-shrub layer (C), which confirms our assumption (Figure 6A). The presence of the significant proportion of broad-leaved species (oak, linden, and ash) in groups (#3 and #4) with the nemoral composition of species of subordinate layers indicates the possibility of demutation of broad-leaved types of communities.

In pine forests, the absence of pine regeneration in communities (#1–4) of automorphic habitats indicates the secondary origin of pine forests after fires and onsite cutting, as well as their artificial origin (Figure 6B).

In one case the reforestation is accompanied by the active demutation of spruce forests; in another case the reforestation takes place in habitats with nutrient-rich soils with broad-leaved species, which push out pine and pine-spruce communities during the next few decades. Pine forests of the Mh group (#5) with species of more southern flora have viable juvenile seedlings in the ground layer and thus can survive on steep river valley slopes. It is more obvious that pine is represented in the wide age range of communities of the DshHSh group (#6), which indicates the stability of the pine cenopopulation under hydromorphic conditions and without strong competition from other species.

The direction of successional dynamics is obviously associated with certain environmental conditions. The analysis of the division of six groups of communities of pine and spruce forests according to the species composition of ground layers in ordination space was performed, which confirmed their differences in the ecological conditions of habitats (Figure 7). The results of correlation with the NMDS ordination axes in an ecological space showed that changes in light, soil pH, and nutrients abundance are associated with the first axis of variation, and soil moisture is associated with the second axis. The highest correlation is observed with soil pH (R² = 0.67) and soil nutrients (R² = 0.74) (

Table A1. Direction of ecological factors vectors and their correlation with points coordinates in NMDS ordination.

Eecological Factors	NMDS1	NMDS2	R2	pr (>r)
L	0.95145	0.30779	0.3978	0.01
R	−0.91634	−0.40041	0.6658	0.01
M	0.48422	−0.87495	0.0288	0.01
N	−0.98414	−0.17742	0.7424	0.01

¹ Ecological factors: L—light, R—soil pH, M—moisture, N—nutrients.

, Appendix A).

The main ecological factors of differentiation of coniferous forests are nutrients abundance and soil pH. However, the role of moisture and insolation in pine forests is significantly higher than that in spruce forests. Pine forests are found both in bogs and on dry sandy soils. Tree canopy depending on the participation of other tree species and environmental conditions, forms different insolation conditions for subordinate layers.

3.2. Dynamic Model of Coniferous Forests

The DM of the identified types of formations is prepared basing on the analysis of age range of cenopopulations of main tree species of pine and spruce communities (Figure 8). Coniferous forests aged 60–100 years are considered, which, depending on ecotopic habitat in the conditions of the existing forest management regime, have a vector of development towards either coniferous or broad-leaved forests. If the duration of qualitative stage of change in composition of forest communities is considered on average 20 years, then the composition of spruce and spruce–small-leaved communities, according to our estimates, approaches the coniferous–broad-leaved forests of the zonal type in 40–60 years, spruce-pine forests—in 50–70 years, pine—for 60–80 years, depending on the correspondence of plantations to optimal habitat conditions. In accordance with the pattern of communities’ distribution in ordination space, there is a general trend in the development of communities towards an increase in soil richness and an increase in soil pH. In some cases of spruce plantations, there is a retrograde movement towards mixed spruce–small-leaved communities, as mentioned above, due to poor quality of tending. Beyond the framework of the DM the meadow non-forest areas, the small-leaved secondary forests, the herb pine forests, and the sphagnum pine forests. The changes of composition and distribution over the area of such pine forests are assumed to be quasi-stationary.

The principles of dynamics of coniferous communities are used in interpretation of the RM of forest cover.

3.3. The Results of Retrospective Modeling 30 Years Ago

In accordance with the proposed methodology, several models of forest formations for the period of 2010 were tested (

Table 2. Kappa and overall accuracy of different models of forest formations.

Classifier	LibSVM	Boost	Decision Tree	Normal Bayes
Parameters	Kernel type: Gaussian radial basis function Model type: C-support vector classification	Gentle AdaBoost	Max depth of tree 100 Min number of samples in each node 100	n/a
Kappa	0.40	0.04	0.39	0.34
Overall accuracy	0.58	0.21	0.55	0.49
Classifier	Random Forests	KNN	Shark Random Forest
Parameters	Max depth of tree 25 Min number of samples 10 Max number of trees 100	Number of neighbors 32	Max number of trees 50 Min size of the node 25
Kappa	0.47	0.41	0.42
Overall accuracy	0.62	0.58	0.59

) and compared. The best result is shown for random forest model. Overall accuracy for an independent test sample (20% of the selected points) was 0.63 (Kappa = 0.47) which is moderate [94]. This is consistent with our study of modeling the MR forest formations [84]. The confusion matrix is shown in

Table 3. Confusion matrix test sample (random forest).

	Test Sample
	Formations (Sample Plots)
Formations (Modeling)	Spruce	Spruce–Small-Leaved	Pine-Spruce	Pine	Deciduous	Non-Forest	User Accuracy
Spruce	35	3	4	3	5	0	0.70
Spruce–small-leaved	13	10	1	1	15	0	0.25
Pine- spruce	3	0	0	2	0	0	n/a
Pine	8	0	2	25	4	1	0.63
Deciduous	11	5	0	3	79	2	0.79
Non-forest	1	0	0	1	2	1	0.20
Producer accuracy	0.30	0.56	n/a	0.71	0.75	0.25
Kappa							0.47
Overall accuracy							0.63

. The deviation between Kappa and Overall accuracy is caused by two factors which influence Kappa sensitivity: (i) different size of samples between different formations, and (ii) non-equal level of agreement for different formations. Both factors are especially typical for pine-spruce formation. This formation has the smallest amount of relevés and, consequently, small percentage of test sample.

The 2010 mosaic was used to obtain spectral signatures and to train classifier. Then the classifier was used for 2010 mosaic and 1990 mosaic, maps for the corresponding periods were made (Figure 9).

There is a clear decline in the area of spruce forests, especially after the 2010. The area of spruce–small-leaved forests grows symmetrically. Pine-spruce forests have weakly growth dynamics. The area of pine forests increases slightly by the 2010 and then falls below the level of 1990. The area of deciduous forests loses a few percent by the 2010 and then compensates for losses by the 2020 (

Table 4. Composition of coniferous formations in the periods 1990, 2010, and 2020, %.

Type of Forest	1990	2010	2020
Forest cover (our data)	56.85	55.87	48.57
Forest cover (Global Forest Watch)	n/a	52.75 1	48.57 2
Area of formations (% from forest cover according to our data)
Spruce	19.61	14.21	7.04
Spruce–small-leaved	6.37	7.33	18.09
Pine–spruce	1.90	2.65	2.51
Pine	17.30	22.13	15.96
Deciduous forest	54.83	51.95	56.40

¹ Our results independent from Global forest watch (GFW). ² GFW forest mask used in study for Forest formations in 2020.

).

Below is the analysis of cross-tabulation for four coniferous forest formations between the periods 1990 and 2020 (

Table A4. Changes in the composition of coniferous formations in the period 1990-2020, loss of the area of forest formations relative to the 1990 level, %.

		2020						Loss of Forest Formation from Level 1990
	Formations	Spruce	Spruce–Small-Leaved	Pine–Spruce	Pine	Deciduous	Non-Forest
1990	Spruce	13.1	16.3	6.1	28.0	14.7	21.7	86.9
	Spruce–small-leaved	12.5	38.5	2.1	9.3	30.5	7.1	61.5
	Pine-spruce	10.5	7.5	10.9	29.0	12.6	29.4	89.1
	Pine	4.4	6.9	2.1	26.7	23.3	36.7	73.3
	Deciduous	2.9	14.4	0.3	3.5	58.4	20.5	41.6
	Non-forest	0.2	1.0	0.1	0.8	12.3	85.7	14.3

and

Table A5. Changes in the composition of coniferous formations in the period 1990-2020, gain of the area of forest formations relative to the 2020 level, %.

		2020
	Formations	Spruce	Spruce–Small-Leaved	Pine–Spruce	Pine	Deciduous	Non-Forest
1990	Spruce	42.8	20.7	56.1	40.3	6.0	5.9
	Spruce–small-leaved	13.3	15.9	6.1	4.4	4.0	0.6
	Pine-spruce	3.3	0.9	9.7	4.0	0.5	0.8
	Pine	12.5	7.7	17.2	33.8	8.4	8.8
	Deciduous	26.2	51.1	7.8	14.2	66.4	15.5
	Non-forest	1.9	3.7	3.2	3.3	14.7	68.4
Gain of forest formation from level 2020		57.2%	84.1%	90.3%	66.2%	33.6%	31.6%

). The original spruce forests lost about 87% of the area, mainly due to the transition to pine (28%) and deforestation (22%). The growth of spruce forests amounted to 57%—mainly due to transition from deciduous (26%). Spruce–small-leaved forests lost 62% of the area mainly due to transition to deciduous (30%). The growth of spruce–small-leaved forests was 84%, mainly due to transition from deciduous forests (51%). Pine-spruce forests lost 89%, equally due to the transition to pine forests or deforestation (29% each). The growth of pine-spruce forests was 90%, mainly due to transition to spruce forests. Pine forests lost 73% of the area, mainly due to deforestation (37%) and transition to deciduous forests (23%). The growth of pine forests was 66%, mainly due to transition from spruce forests. Deciduous forests lost 42% of their area, mainly due to deforestation (20%), and gained 34% mainly due to overgrowth (15%). Basing on the results of cross tabulation the diagram of the successional dynamics of the identified types of formations is prepared (RM) (Figure 10).

4. Discussion

We generalized the accumulated experience of numerous studies carried out in the MR over the past years with respect to the typological composition of coniferous communities [23,28,95,96,97]. According to our estimates [84], the area of spruce forests amounts 21.7% (including spruce–small-leaved forests) and the area of pine forests amounts 18.5% (including spruce-pine forests). This proportion is close to official statistics [98].

The current state of MR forest cover reflects the results of the post-Soviet socio-economic extensive forest management model on the one hand, and reorientation towards the conservation of biodiversity and ecosystem services of the forest cover on the other hand. Due to the multidirectional action of socio-economic factors in the MR (deforestation, on the one hand, and silviculture and fire resistance, on the other), the natural cycles of autogenous successions are disrupted. Most of the MR coniferous forests are plantations, with the most active silviculture activities after WWII. In the so-called “complex” pine and spruce forests with the participation of linden and oak of small herb–broad herb layer, the percentage of plantations (according to geobotanical relevés) is about 80%.

The main trend of forest management in the current period of time differs from previous eras by reorientation towards the conservation of forest biodiversity and ecosystem services [99]. In accordance with the Forest Code of the Russian Federation, the MR forests belong to protection group with the restriction of commercial felling and thinning. The proportion of specially protected natural areas (SPNA) according to our data is 6.83% [100] and is constantly increasing. Nevertheless, the condition of forest plantations is considered unsatisfactory due to insufficient sanitary care.

In scientific literature, when modeling the successional dynamics of forests, the possibility of using “Markov chains” has long been discussed [101,102,103,104]. We did not use this approach, since the main disadvantage of dynamic states of communities of models of “Markov chains” is the deliberate omission of any causal mechanisms and the independence of the future state from the prehistory of formation process [105].

Within the framework of identified types of communities, based on the analysis of the age range of cenopopulations of the main forest-forming tree species of pine and spruce communities, the origin of forests, the scheme of succession dynamics was drawn up. In accordance with the pattern of distribution of communities in ordination space, there is a general trend in the development of communities towards an increase in soil richness and a decrease in acidity.

If the duration of qualitative stage of the transformation of forest communities is taken on average 20 years, then the composition of spruce and spruce–small-leaved communities approach to the stage of coniferous–broad-leaved forests of the zonal type in 40–60 years, spruce-pine forests—in 50–70 years, pine—in 60–80 years, depending on the conformity of the forest plantations to optimal habitat conditions. This rate of transition of successional stages is close to the pattern of forest succession (overgrowing of fallows) in the series of spruce and mixed (with oak, linden and spruce undergrowth) communities in the south of the MR [106]. Based on the DM, it is expected the redistribution of pine and spruce and mixed forests as the result of restorative successional dynamics in coniferous plantations and on overgrown agricultural lands, as well as the reduction of pine plantations due to the replacement of pine forests with spruce and spruce–broad-leaved in automorphic habitat.

The combination of formation maps of 1990 and 2020 made it possible to reveal the RM of the organization of vegetation cover, indicating the directions of successional dynamics of forest formations in the study area.

Firstly, the aspects of successional dynamics of plantation monodominant forests are identified. For spruce forests, the noticeable decrease (from 19.6 to 7%) of area is observed, as well as partial replacement by pine forests. This may be due to aging and natural loss of monodominant spruce plantations, which, as mentioned above, played significant role in the second half of the 20th century. Thus, hurricanes at the end of last century triggered an outbreak of the bark beetle Ips typographus in the 2000s, which caused the mass defoliation of spruce forests. In 2000–2006, the harmful insects continued to be main cause of forest stands dieback. The forest pathological situation of the MR forests was relatively sustainable in the 2007–2009 (the area of the dead forests was 666–1061 ha annually). The area of dead forest stands increased sharply to 21.4 thousand hectares in the 2010 due to abnormally dry and hot summers and persistent ground fires (20 times compared to 2009) [66]. Outbreaks of the bark beetle in spruce forests have recently been especially typical [67,68,69]. Spruce–small-leaved forests show a strictly opposite trend. During the study period, their proportion increased from 6 to 18%, which is explained by the lack of care for spruce plantations in recent years and their overgrowth primarily with birch.

The modeling approaches for periods of the 1990/2010 and the 2020 differ in many terms—different remote sensing scanners (TM and OLI), different combinations of environmental variables, and different modeling algorithms. However they use the same field data and demonstrate relatively similar spatial pattern in the MR which is consistent with the official data [51] and earlier studies [28].

Landsat 5 (TM scanner) mosaics can be produced only for the 1982-2012 period. Landsat 7 (ETM+ which is almost similar to TM) performed well only until 2003 which makes it not suitable to use in long-term time-series unfortunately. Landsat 8 (OLI scanner) differs from Landsat 5 in terms of spectral bands, but it is the only mid-scale source of remote sensing data which can be used to compare current and historical land/forest cover. Considering the limitations mainly caused by the differences between spectral bands of Landsat 5 and 8 we may assume that our approach is suitable for comparing the land/forest cover of periods 1990/2010 and 2020. The overall accuracy of forest formations modeling is under influence of numerous factors. First, the fluctuating balance between deforestation (recreational impact, road, and construction infrastructure) and afforestation (tillage abandoning, silviculture). This factor significantly disrupts the natural spatial structure of forest formations. The MR is located in the transition zone between deciduous-coniferous forests in the center and North and broadleaf forests in the South, consequently mixed polydominant forests dominate, which are difficult to separate [107]. Increase of deciduous forests due to climate warming is also noted [108].

Second, environmental variables cause uncertainties such as nadir angle, radiometric correction and atmospheric transparency, each of them varies 5 to 7.5% [109] and one more important issue. The area of the Landsat pixel exceeds the sample plot which means the Landsat is too coarse in terms of statistics and sampling [110].

Third, the uneven spatial distribution of field data is critical for modeling and needs preliminary stratification and sample equalization. In our case we are constrained to use already collected field dataset with bias of spatial location of relevés [84].

According to our model, the restoration of pine, presumably plantation was noted for some spruce forests. For plantation pine forests, in turn, there is a partial replacement by deciduous species, mainly broad-leaved. Such picture as a whole is typical for the group of pine forests with participation of linden, small herb–broad herb oak (#3) and broad herb (#4) or mixed pine forests of the MR [111], and indicates “nemoralization” of middle-aged plantations during the transition to the stage of more “mature communities” [112]. In the composition of small herb–broad herb group of communities the largest proportion of pine plantations (66%) is noted.

Secondly, the equilibrium mutual succession transitions are observed for mixed forests. Mixed pine-spruce forests tend to transition to and from pine forests. The same can be said for a pair of spruce–small-leaved–deciduous forests.

Thirdly, there is a pattern of unilateral transition without mutual compensation. This is the case of deciduous forests, which successively develop into spruce and/or spruce–small-leaved forests, and spruce forests, which become more complex and turn into pine forests. The small area of pine forests has also been identified. These forests are successively transformed into deciduous (broad-leaved) forests. However, the area of such one-sided transformations is small. Probably, this is pattern of natural dynamics.

As a result, on the territory of the MR, there is the obvious decrease of the area of spruce and pine formations and the increase of the area of spruce–small-leaved and deciduous formations. The formation of pine-spruce forests does not exceed 3% of the area; it varies greatly in time and space and can be considered as a vulnerable formation.

5. Conclusions

The development of forest maps on the example of the study area—the Moscow region—made it possible to reveal the dualism of anthropogenic and natural factors: mosaic structure and significant disturbance of forest cover is determined at local level, and correspondence of its organization to general natural and geographical laws at regional level. Despite the current forest management regime in the MR (high rates of urban development in the central regions, recreational load, prohibition of final cutting and insufficient care of plantations in the last 20 years), and the artificial origin of most coniferous forests, their composition is close to the nominally primary of forest communities.

Two models of succession dynamics (DM and RM) are compared. In general, two models are similar in the main directions of the dynamics of forest communities. At the same time, features were identified which indicate risks and threats to sustainable development of coniferous forests under the existing forest management regime in the MR. Reduction of spruce and pine forests is main trend in transformation of forest cover in the MR. This conclusion is obtained based on retrospective dynamic model and supported by the insufficient care for coniferous plantations.

Constant multidirectional changes in the forest cover of dynamically developing MR (sometimes degradation) presuppose the need of permanent update of remote and field data.